Datareductie en stabiliteitsanalyse van IC3328

Faculteit Wetenschappen Vakgroep Fysica en Sterrenkunde

Datareductie en stabiliteitsanalyse van IC3328

door Jelle Dhaene

Promotor: Prof. Dr. S. De Rijcke

Masterproef ingediend tot het behalen van de graad van Master in de Fysica en Sterrenkunde

Academiejaar 2011-2012

Dankwoord Deze masterproef kon onmogelijk tot stand komen zonder de hulp van een aantal mensen, die ik bij deze wil bedanken. Eerst en vooral zou ik mijn promotor Prof. dr. Sven De Rijcke willen bedanken. Het onderwerp van deze masterproef werd door hem voorgesteld. Ook nam hij de tijd om mij persoonlijk te begeleiden en mij te introduceren tot het ESO-MIDAS programma en andere Fortran- en python-scripts die door hem beschikbaar werden gesteld. Verder hielp hij mij steeds verder als ik te pas en te onpas zijn bureau kwam binnenvallen om hem opnieuw met vragen te bestoken. Ook dr. Mina Koleva zou ik willen bedanken voor het fitten van een stellaire populatie aan het spectrum van IC3328 en het bepalen van de massa-lichkrachtverhouding van deze galaxie. Verder wil ik mijn ouders bedanken, die het in de eerste plaats voor mij mogelijk gemaakt hebben om mijn masterdiploma te halen. Zowel voor hun financiële als morele steun. En verder wil ik mijn moeder bedanken om het gros van de typefouten uit deze masterproef te halen. Tenslotte wil ik mijn medestudenten bedanken, die mij gezelschap hielden in het serverlokaal in S9 waar we samen aan onze masterproef werkten. Vooral de laatste (warme) dagen waren de lastigste en dan was het leuk om jullie gezelschap te hebben.

Inhoudsopgave 1 Inleiding

7

2 Data verzamelen en verwerken 2.1 De CCD camera . . . . . . . . . . . . . . . . 2.2 Telescopen en detectoren . . . . . . . . . . . . 2.2.1 VLT . . . . . . . . . . . . . . . . . . . 2.2.2 HST . . . . . . . . . . . . . . . . . . . 2.3 Data verkrijgen . . . . . . . . . . . . . . . . . 2.4 Basisprincipes van CCD data reductie . . . . 2.4.1 Verschillende effecten van een CCD . . 2.4.2 ESO-MIDAS . . . . . . . . . . . . . . 2.5 Datareductie FORS . . . . . . . . . . . . . . 2.5.1 Reductie van de foto’s . . . . . . . . . 2.5.2 Bepaling van het fotometrisch nulpunt 2.6 Datareductie WFC . . . . . . . . . . . . . . . 2.7 Gereduceerde data en opmerkingen . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

9 9 12 12 14 16 16 17 18 19 20 23 26 27

3 Analyse van de data 3.1 Model fitten aan galaxie . . . . . . . . . . . . . . . . . . 3.1.1 Schijnbare vorm van de galaxie . . . . . . . . . . 3.1.2 Code . . . . . . . . . . . . . . . . . . . . . . . . . 3.2 Resultaten van het model . . . . . . . . . . . . . . . . . 3.3 Fitten van een profiel aan de oppervlaktehelderheid . . . 3.3.1 Profielen voor de intensiteit . . . . . . . . . . . . 3.3.2 De point-spread functie . . . . . . . . . . . . . . 3.3.3 Bepaling standaardafwijking van de PSF . . . . . 3.3.4 Convolutie van de PSF en de intensiteit . . . . . 3.3.5 Het fitten van de data . . . . . . . . . . . . . . . 3.3.6 De totale intensiteit en de schijnbare magnitude 3.4 Het dichtheidsprofiel . . . . . . . . . . . . . . . . . . . . 3.4.1 De massa-lichtkrachtverhouding . . . . . . . . . . 3.4.2 Dichtheidsprofiel van de sferische componenten . 3.4.3 Oppervlaktedichtheidsprofiel van de schijf . . . . 3.4.4 De totale massa van de galaxie . . . . . . . . . . 3.5 De potentiaal . . . . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . .

29 29 29 30 34 39 39 40 42 44 46 52 54 54 55 58 58 59

5

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

. . . . . . . . . . . . .

6

INHOUDSOPGAVE 3.5.1 3.5.2

De bijdrage van de materie . . . . . . . . . . . . . . . . . . . . . . . 59 De bijdrage van de donkere materie . . . . . . . . . . . . . . . . . . 62

4 De distributiefunctie 4.1 Kinematische data . . . . . . . . . . . . . . . . . . . . . . . 4.2 De distributiefunctie en haar momenten . . . . . . . . . . . 4.2.1 De distributiefunctie . . . . . . . . . . . . . . . . . . 4.2.2 De momenten . . . . . . . . . . . . . . . . . . . . . . 4.3 Integraal van de beweging en theorema van Jeans . . . . . . 4.3.1 Constante van beweging . . . . . . . . . . . . . . . . 4.3.2 Integraal van de beweging . . . . . . . . . . . . . . . 4.3.3 Theorema van Jeans . . . . . . . . . . . . . . . . . . 4.4 Model voor de schijf . . . . . . . . . . . . . . . . . . . . . . 4.4.1 De distributiefunctie in functie van de integralen van 4.4.2 Vorm van de distributiefunctie . . . . . . . . . . . . 4.4.3 De fit aan de data . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . beweging . . . . . . . . . . . .

. . . . . . . . . . . .

. . . . . . . . . . . .

. . . . . . . . . . . .

67 67 71 71 71 72 72 73 73 73 74 74 75

5 Stabiliteitsanalyse 5.1 De eerste-orde Boltzmannvergelijking . . . . . 5.2 De respons van de distributiefunctie . . . . . 5.3 De respons van de dichtheid . . . . . . . . . . 5.4 Zelf-consistente storingen . . . . . . . . . . . 5.4.1 Eigenmodes van de schijf . . . . . . . 5.4.2 Familie van potentialen en dichtheden 5.4.3 Zoeken naar eigenmodes van de schijf

. . . . . . .

. . . . . . .

. . . . . . .

. . . . . . .

81 81 82 85 86 86 87 88

. . . . . . .

. . . . . . .

. . . . . . .

. . . . . . .

. . . . . . .

. . . . . . .

. . . . . . .

. . . . . . .

. . . . . . .

. . . . . . .

. . . . . . .

. . . . . . .

. . . . . . .

6 Besluit

89

A De epicykelbenadering

91

Hoofdstuk 1

Inleiding Zoals de titel al laat vermoeden gaat deze masterproef over de datareductie en stabiliteitsanalyse van IC3328. Dit is een dwerggalaxie gelegen in de Virgocluster op een afstand van ongeveer 16 Mpc en die geklasseerd is als een elliptische dwerggalaxie van het type dE1. Op het eerste lijkt dit een normale dwerggalaxie te zijn, maar na nader onderzoek blijkt deze toch een verborgen aspect te bevatten. Een spiraalstructuur werd in deze dwerg gevonden door Jerjen et al (2000) zoals te zien is in figuur 1.1. De spiraalstructuur is gelijkaardig aan die van M51, maar de prominente HII gebieden, die duiden op de aanwezigheid van gas, en het stof zijn afwezig in de dwerg. De spiraalstructuur zou verklaard kunnen worden door een fenomeen dat ’galaxy harassment’ genoemd wordt, zoals beschreven door Moore (1998). Deze theorie stelt dat de waargenomen elliptische dwerggalaxieën de overblijfselen zouden zijn van lage-massa spiraalgalaxieën. Door interacties met grote massieve galaxieën in een cluster zouden de dwerg spiraalgalaxieën omgevormd worden tot elliptische dwerggalaxieën. Uitgaande van deze theorie zou IC3328 dan interacties ondergaan hebben met massieve galaxieën en zo van een dwerg spiraalgalaxie naar een elliptische dwerg aan het evolueren zijn. Hier wordt de basis gelegd voor een zoektocht naar de verklaring van het aanwezige spiraalpatroon door gebruik te maken van de Lin-Shu hypothese, die stelt dat de spiraalarmen in schijfgalaxieën dichtheidsgolven zijn. In deze masterproef wordt er data gebruikt van de HST (Hubble Space Telescope) en de VLT (Very Large Telescope). De data gebruikt van de VLT is dezelfde als deze door Jerjen et al gebruikt in 2000, maar wordt uitgebreid met HST data, die een veel hogere resolutie heeft, zodat een beter beeld te vormen is van het centrum van de galaxie. De HST data zal vooral handig zijn in het beschrijven van de nucleus van de galaxie. In hoofdstuk 2 zal er begonnen worden met de basisbenodigdheden nodig om te observeren te bespreken. Eerst zal kort uitgelegd worden hoe een CCD-camera precies werkt en worden de gebruikte telescopen en instrumenten besproken. Hierna wordt overgegaan op de datareductie en wordt een eerste blik op de data geworpen. Zo werd een astero¨ıde gespot die door het gezichtsveld van de telescoop bewoog gedurende het moment van de waarnemingen. Eénmaal de data gereduceerd is, wordt in hoofdstuk 3 overgegaan tot de effectieve analyse ervan. Er zal eerst een model, bestaande uit concentrische ellipsen, aan de galaxie gefit worden, om hieruit het oppervlaktehelderheidsprofiel te bepalen. De fit aan het oppervlaktehelderheidsprofiel gebeurt op zodanige wijze dat we de 3 componenten (nucleus, 7

8

HOOFDSTUK 1. INLEIDING

Figuur 1.1: Foto van de dwerggalaxie IC3328 in de R-band. Rechts is een bewerkte foto c Jerjen et al.) te zien waarop de spiraalstructuur duidelijk zichtbaar is. ( schijf en halo) van de galaxie uit elkaar zullen kunnen halen. Voor de fit worden Sérsicen exponentiële profielen gebruikt. Wanneer de intensiteit gekend is, kunnen we deze deprojecteren om zo de emmisiviteit van de galaxie te bepalen, die direct gerelateerd is, via de massa-lichtkracht verhouding, aan het dichtheidsprofiel van de galaxie. Hierna kan via de Poissonvergelijking de potentiaal uit de dichtheid bepaald worden. Hoofdstuk 4 heeft als doel een distributiefunctie te fitten aan de schijf van de galaxie met behulp van de tot nu toe verzamelde data. De dichtheid en de potentiaal volstaan echter niet om dit te doen. Hiervoor zal ook kinematische data gebruikt worden, verzameld door Simien en Prugniel (2002). Al deze data zal gebruikt worden om een distributiefunctie te reconstrueren. Dit is echter niet zo eenvoudig en dit kan enkel bereikt worden door tegenroterende sterren in de distributiefunctie in rekening te brengen. Op zich niet zo’n groot probleem, tot in hoofdstuk 5 overgegaan wordt op de stabiliteitsanalyse. De code die ontwikkeld werd door Dury et al (2008) is niet voorzien om rekening te houden met deze tegenroterende sterren. Dus voor verdere analyse zou een aangepaste versie van deze code nodig zijn. Niettegenstaande wordt toch de analyse besproken en een manier beschreven hoe op een relatief eenvoudige wijze op zoek kan gegaan worden naar de eigenmodes van de schijf.

Hoofdstuk 2

Data verzamelen en verwerken De basis van astronomische observaties begint bij de instrumenten die gebruikt worden om de sterren, galaxieën en talrijke andere objecten waaraan het universum rijk is waar te nemen. De eerste telescoop werd gemaakt in 1608 en sindsdien zijn een resem aan verschillende toestellen door de handen van astronomen gepasseerd en tegenwoordig zijn het gigantische toestellen geworden die in de hooggebergtes op aarde geplaatst worden of in een baan rond de aarde cirkelen. In het begin kon men enkel door een telescoop kijken zonder dat het beeld effectief opgeslagen kon worden, of men moest met de hand een tekening maken. Met de komst van de fotografische platen aan het einde van de 19de eeuw kon dit probleem verholpen worden en konden foto’s genomen worden van hetgeen men zag door de telescoop. In de jaren ’70 van de vorige eeuw werden dan de eerste elektronische detectoren gebruikt die als groot voordeel tegenover de fotografische platen hebben dat ze meer dan 1 keer gebruikt kunnen worden. Tegenwoordig worden CCD’s (Charged Coupled Devices) gebruikt om de beelden vast te leggen. In dit hoofdstuk zullen we eerst de basisprincipes van een CCD bespreken om daarna de telescopen die gebruikt werden voor het waarnemen van de data eens wat nader te bekijken en om tenslotte de datareductie te bespreken en een eerste blik op de gereduceerde foto’s te werpen.

2.1

De CCD camera

Basismateriaal van een CCD De basis voor elektronische apparaten zoals een CCD zijn halfgeleiders. De elektrische eigenschappen van halfgeleiders kunnen drastisch veranderen bij foto-excitatie door absorptie van UV, zichtbaar of infrarood licht. De elementaire halfgeleiders zijn germanium en silicium. Beide hebben 4 valentie-elektronen en beide zullen een kristal vormen met diamantstructuur (Fig. 2.1). In deze structuur zal ieder atoom 4 dichtste naburen hebben en zal het met elk van hen een valentie-elektron delen. Werking van een CCD Een CCD bestaat uit een matrix van identieke metaaloxide halfgeleider capaciteiten gevormd op een siliciumsubstraat. Ieder element van de matrix wordt een fotodetector genoemd en ieder inkomend foton raakt het silicium in zo een fotodetector en wordt ge9

10

HOOFDSTUK 2. DATA VERZAMELEN EN VERWERKEN

Figuur 2.1: De diamant roosterstrucuur, waaruit de germanium en siliciumkristallen ook c University of Illinois) bestaan. ( absorbeerd als ze de juiste energie bezitten. In een siliciumkristal is de energiesprong tussen de valantie- en conductieband 1.1 eV, dus enkel fotonen met een hogere energie kunnen elektron-gat paren creëren. Wanneer zo een paar gecreëerd is, zal het na ongeveer 100 µs terug recombineren. De elektronen moeten dus op hun plaats gehouden worden tot de CCD uitgelezen wordt. Iedere pixel heeft nog een onderliggende structuur zoals weergegeven op figuur 2.2. Verschillende spanningen kunnen aangelegd worden op de 3 elektroden die we poorten noemen en die te zien zijn op de figuur. Op deze manier wordt iedere pixel in 3 subpixels onderverdeeld. Ieder elektron, gecreëerd in 1 van de 3 subpixels, zal in dezelfde diepste potentiaalput terechtkomen. Op het einde van de belichting van de CCD ziet het systeem eruit zoals in de bovenste rij weergegeven op de figuur, zodat alle elektronen die gecreëerd werden in de pixel met behulp van een potentiaalput op hun plaats gehouden worden. Uitlezen van een CCD Nadat de CCD belicht is, moet hij ook nog uitgelezen worden. Het principe van het uitlezen van een CCD is schematisch weergegeven in figuur 2.3. Als alle elektronen opgevangen zijn, dan kunnen de potentialen van de 3 subpixels uit figuur 2.2 gewisseld worden zodat de elektronen een pixel kunnen opschuiven naar rechts zoals weergegeven in de onderste 2 figuren. Aan de rechterkant van de CCD, waarin de pixels in figuur 2.3 vervangen zijn door emmers, zit een rij pixels die niet belicht wordt, maar waarin de elektronen terechtkomen bij het uitlezen. Dus als alle elektronen 1 pixel naar rechts opschuiven, dan geraken er elektronen in de onbelichte rij pixels. De elektronen in deze rij kunnen dan met behulp van eenzelfde principe naar boven (of onder) verplaatst worden en zo kan voor iedere pixel in deze rij het aantal elektronen geteld worden. Nadat de rij uitgelezen is, kunnen opnieuw de spanningen van de matrix van pixels aangepast worden om de volgende rij elektronen in de uitleesrij te krijgen. Op basis van dit principe kan de CCD uitgelezen worden. De lading die uitgelezen wordt per pixel wordt gemeten als een spanning en wordt geconverteerd naar een digitaal nummer. De ruwe data die uit de CCD komt staat dan ook in ADU (Analog to Digital Unit) per pixel. Het terug converteren naar elektronen per pixel kan dan gedaan worden met behulp van de gain (meer hierover in het stuk over datareductie). Om de uitleestijden van een CCD te verkorten kan deze opgedeeld worden in verschillende delen die simultaan uitgelezen kunnen worden.

2.1. DE CCD CAMERA

11

Figuur 2.2: Schematische weergave van de onderverdeling van een pixel van een CCD camera en de potentiaalputten die de gecreëerde elektronen opvangen. De pakketjes elektronen kunnen van pixel naar pixel verplaatst worden bij het uitlezen van de camera c Observational techniques in door de spanningen over de 3 subpixels te veranderen. ( astronomy)

Figuur 2.3: Schematische weergave van het uitlezen van een CCD waarbij iedere emmer c Nikon, een pixel voorstelt die elektronen (voorgesteld door druppels water) opvangt. ( www.microscopyu.com)

12

2.2


Telescopen en detectoren

Voor deze masterproef werden datasets van 2 verschillende telescopen gebruikt. De eerste dataset werd bekomen met het FORS1 (FOcal Reducer and low dispersion Spectrograph) instrument van de VLT (Very Large Telescope array) die zich in Cerro Paranal, Chili bevindt (Fig 2.4(a)). Voor de tweede dataset werd gebruik gemaakt van WFC (Wide Field Channel) detector van het ACS (Advanced Camera for Surveys) instrument van de HST (Hubble Space Telescope) die zich in een baan rond de aarde bevindt ongeveer 600km boven het aardoppervlak (Fig 2.4(b)).

(a) VLT

(b) HST

c ESO) en (b) De Hubble ruimtetelescoop ( c NASA) Figuur 2.4: (a) De VLT in Chili (

2.2.1

VLT

Het observatorium op Paranal wordt gerund door de ESO (the European Southern Observatory). De ESO is een organisatie waar 14 Europese landen (waaronder België) deel van uitmaken. De organisatie beheert 2 sterrenwachten in Chili, één op La Silla en één op Cerro Paranal en momenteel wordt een derde observatorium op het Chajnantor plateau in de Atacama woestijn afgewerkt waar ALMA (Atacama Large Millimeter Array) volledig operationeel zou moeten zijn tegen het einde van 2012. De huidige opstelling op Paranal ziet eruit zoals in figuur 2.5 weergegeven. De VLT bestaat uit 4 grote telescopen, de Unit Telescopes (UT), genaamd Antu, Kueyen, Melipal en Yepun. Deze namen zijn in het Mapudungun, de taal die gesproken wordt door de Mapuche in centraal Chili en ze betekenen respectievelijk zon, maan, zuiderkruis en Venus. Deze telescopen zijn gemonteerd in alt-azimuth montage en hebben een primaire spiegel van 8.2m diameter. Verder zijn er nog 4 kleinere verplaatsbare Auxiliary Telescopes (AT) met primaire spiegels van 1.8m diameter. We zien ook nog de VST (VLT Survey Telescope) en VISTA (the Visible and Infrared Survey Telescope for Astronomy). Deze 2 telescopen worden gebruikt om de hemel af te scannen naar interessante objecten die dan later met de VLT kunnen bestudeerd worden. De Unit telescopen zijn ondergebracht in thermisch gecontroleerde gebouwen die synchroon meedraaien met de telescopen, dit om de condities van de observaties zo goed mogelijk te houden zoals bijvoorbeeld de luchtturbulentie in de telescoopbuizen. De Eerste

2.2. TELESCOPEN EN DETECTOREN

13

van de 4 UT’s werd in gebruik genomen voor wetenschappelijke doeleinden in april 1999. Vandaag de dag zijn alle 4 de UT’s en de AT’s in gebruik. De 4 UT’s kunnen gecombineerd worden in de VLTI (the Very Large Telescope Interferometer) maar ze worden meestal individueel gebruikt. Slechts een beperkt aantal nachten per jaar werken ze samen om de VLTI te vormen. De 4 kleinere AT’s zijn wel altijd beschikbaar om de VLTI iedere nacht te gebruiken.

Figuur 2.5: Een overzicht van de opstelling op de site op Cerro Paranal, VISTA bevindt c ESO) zich ongeveer een kilometer noordwaarts van de andere telescopen. (

FORS1 De door ons gebruikte data werden verzameld met het FORS1 instrument (Fig 2.6(a)) en was tijdens het verzamelen van de data gemonteerd in het Cassegrain brandpunt van UT1. Deze detector heeft een gezichtsveld ter grootte van 6.8 × 6.8 boogminuten in standaard resolutie en 4.25 × 4.25 boogminuten in hoge resolutie in een golflengtegebied tussen de 300 en de 1100 nm. 1 pixel van de CCD is ongeveer 15 µm groot en correspondeert met ongeveer 0.2 boogseconden aan de hemel in normale mode en 0.125 boogseconden in hogeresolutie mode. De FORS1 detector werd in 1998 op UT1 ge¨ınstalleerd en verplaatst naar UT2 in 2004 en sinds 2009 is het vervangen op UT2 door X-Shooter waarbij sommige observatiefuncties verplaatst werden naar FORS2 die nog steeds actief is op UT1. Voor de door ons gebruikte observaties werd de FORS1 gebruikt in imaging mode met de Rband filter met een effectieve golflengte van 655 nm en een FWHM (Full width at half maximum) van 165.0 nm, zoals weergegeven in figuur 2.6(b) .

14


(a) FORS1

(b) R-band Bessel filter

c ESO) en (b) de Figuur 2.6: (a) Het FORS1 instrument ge¨ınstalleerd op een UT ( c Spanish Virtual Observatory). gebruikte R-band filter (

2.2.2

HST

De Hubble ruimtetelescoop werd gebouwd door NASA (National Aeronautics and Space Administration) met bijdragen van ESA (European Space Agency). Het is samen met CGRO (Compton Gamma Ray Observatory), CXO (Chandra X-ray Observatory) en SST (Spitzer Space Telescope) één van de 4 grote ruimteobservatoria van NASA. Buiten de CGRO zijn ze nog allemaal in gebruik, waarvan de HST de eerste was die met observaties begon nadat hij gelanceerd werd op 24 april 1990. De 11 ton zware telescoop met een primaire spiegel van 2.4m werd gelanceerd met de space shuttle Discovery en in een baan om de aarde gebracht op een hoogte van 569 km met een orbitale snelheid van 7.56 km/s. Met deze snelheid duurt het ongeveer 96-97 minuten om 1 baan rond de aarde af te leggen. Het bereik van de telescoop reikt van het ultraviolet over het zichtbaar licht tot aan het nabij infrarood. Het beschikt voor deze waarnemingen over 4 hoofdinstrumenten die sinds de lancering in 1990 al door 4 service missies vervangen werden om nog beter waarnemeningen te kunnen uitvoeren. De door ons gebruikte observaties werden verricht met de ACS. De belangrijkste componenten van de HST zijn weergegeven in figuur 2.7 . ACS De ACS werd in maart 2002 door de space shuttle Columbia, die nog geen jaar later op 1 februari 2003 volledig vernietigd werd bij het terugkeren in de atmosfeer na een ruimtemissie, naar de HST gebracht om daar op 7 maart ge¨ınstalleerd te worden in de plaats van de FOC (Faint Object Camera), die tevens het laatste originele instrument was dat nog op de HST aanwezig was sinds zijn lancering. Het ACS bestaat uit 3 verschillende kanalen, het WFC (Wide Field Channel), het HRC (High-Resolution Channel) en het SBC (Solar Blind Channel). Het SBC is een systeem dat geoptimaliseerd is in het UV om lage achtergrond fotonen te tellen, maar dit is nu buiten gebruik sinds het faalde in 2004. Het HRC is momenteel ook buiten gebruik door een elektrische fout en werd zoals de naam het zegt gebruikt om hoge-resolutie foto’s te maken van heldere objecten. Het WFC werd gebruikt om de tot ons beschikbare data te meten en is tevens het meest gebruikte kanaal van ACS. De detector bestaat uit 2 CCD’s van 2048 × 4096 pixels, goed voor een totaal

2.2. TELESCOPEN EN DETECTOREN

15

c NASA). Figuur 2.7: De HST en zijn belangrijkste onderdelen ( van 16 megapixels. De grootte van 1 pixel is ongeveer 15 µm. Het gezichtsveld van het WFC is 202 × 202 boogseconden met een resolutie van 0.05 boogseconden/pixel. Verder kan de detector gebruikt worden over een golflengtegebied tussen 350-1100 nm. Voor onze data werd een filter gebruikt met een effectieve golflengte van 900.8 nm en een FWHM van 127.6 nm. De filter is weergegeven in figuur 2.8(b) .

(a) ACS

(b) Filter

c Wikipedia) en (b) de gebruikte Figuur 2.8: (a) De ACS toen hij nog op aarde stond ( c Spanish Virtual Observatory). filter van het WFC tijdens de observaties (

16

2.3


Data verkrijgen

In principe zijn er 2 manieren om data te verkrijgen. De eerste is dat je zelf een proposal indient en zo observatietijd krijgt voor een telescoop. Een tweede manier is de data gebruiken die iemand anders al verzameld heeft op een bepaalde telescoop. Omdat het moeilijk is om observatietijd te krijgen in de grote observatoria (5 tot 10 keer zoveel aanvragen als dat er observatietijd beschikbaar is) is het aangewezen om vroeger verzamelde data te gebruiken als die voldoet aan de eisen van je project en je enkel zelf een proposal indient als je een specifieke configuratie voor je waarnemeningen nodig hebt, of een specifiek instrument nodig hebt om je object(en) te observeren. Als je zelf een proposal indient en waarnemingen verricht is het vaak zo dat pas na ongeveer 1 jaar de dataset vrijgegeven wordt in de databases van de observatoria, dus dan heb je zelf ongeveer 1 jaar tijd om je datareductie en analyse te doen, voordat anderen toegang krijgen tot de data van je observaties. Voor deze masterproef werd gebruik gemaakt van eerder verzamelde data van IC3328 met de VLT en HST. Deze data is terug te vinden in het ESO archief (http://archive.eso.org /eso/eso archive main.html) en het HST archief (http://hla.stsci.edu/hlaview.html). Je kan in deze archieven zoeken naar de objecten die je wil bestuderen en verdere zoekcriteria opleggen zoals onder andere de filter die gebruikt wordt, de telescoop en het meetinstrument dat gebruikt wordt. In tabel 2.1 wordt wat basisinformatie gegeven over de datasets die hier gebruikt worden. Enkel de data van het te observeren object, IC3328, is opgenomen in de tabel. Verder zullen we ook nog bias frames, flatfields en observaties van standaardsterren nodig hebben om de kalibratie te kunnen doen. Meer info over de laatstgenoemde frames volgt verder. Tabel 2.1: De gebruikte observaties van IC3328(=VCC0856), genomen met FORS1(VLT) en ACS(HST). Voor een ruimtetelescoop is de Luchtmassa niet van toepassing (NVT). Object VCC0856 VCC0856 VCC0856 VCC0856 VCC0856 VCC0856

2.4

Instrument FORS1 FORS1 FORS1 ACS/WFC ACS/WFC ACS/WFC

Observatietijd 1999-07-13T23:27:11.952 1999-07-13T23:34:50.188 1999-07-13T23:42:26.038 2003-07-11T08:41:59 2003-07-11T08:45:47 2003-07-11T08:57:23

Belichtingstijd 399.983 399.983 399.981 90 560 560

Filter R BESS R BESS R BESS F850LP F850LP F850LP

Luchtmassa 1.369 1.396 1.424 NVT NVT NVT

Basisprincipes van CCD data reductie

Eerst beginnen we met het concept van datareductie te beschrijven voor astronomische foto’s, die tegenwoordig bijna allemaal met CCD camera’s getrokken worden. Als een telescoop een foto neemt van een bepaald object zullen er zich problemen voordoen als we de foto’s willen analyseren en de waarden voor de magnituden van onze objecten willen bepalen. De voornaamste problemen zijn hierbij de omzetting bij het uitlezen van de CCD naar een digitaal signaal, de elektronen die gecreëerd worden in het kristalrooster van de CCD zelf en het feit dat niet alle pixels van een CCD even gevoelig zijn voor inkomende fotonen. Deze moeilijkheden worden door de apparatuur zelf veroorzaakt. We zullen ook

2.4. BASISPRINCIPES VAN CCD DATA REDUCTIE

17

nog moeten kijken naar effecten die niet van de apparatuur afhangen zoals kosmische straling en de achtergrond op onze foto’s. Door een profiel aan deze achtergrond te fitten kan de foto hiervoor gecorrigeerd worden door de achtergrond ervan af te trekken. Om te corrigeren voor de effecten van onze apparatuur zelf zullen we gebruik maken van bepaalde frames. De bias frame, dark frame en een flatfield.

2.4.1

Verschillende effecten van een CCD

Bias frame De bias frame wordt uitgelezen met de sluiter van de camera gesloten en heeft een belichtingstijd van 0 seconden. Het is dus eigenlijk gewoon het uitlezen van de CCD camera. De output van de CCD wordt weergegeven in ADU (analog to digital unit), dit kan terug omgezet worden in elektronen per seconde met behulp van de gain die uitgedrukt wordt in elektronen per ADU. Deze gain kan voor ieder uitleescircuit van een CCD verschillen (In sectie 2.1 werd vermeld dat een CCD meerdere uitleescircuits kan hebben om de uitleestijd te beperken). De pixels waar geen licht op invalt vertonen namelijk een kleine variatie in helderheidswaarde. Om ervoor te zorgen dat er geen pixels zouden zijn die een negatieve waarde geven in ADU wordt er een kleine boost gegeven bij de convertie zodat als je een pixel uitleest, die niet standaard op 0 ADU staat, maar bijvoorbeeld op 200 ADU. De variatie van de pixels bij het uitlezen zal ervoor zorgen dat deze verdeeld zijn rond 200 ADU en dan heb je geen last meer van negatieve uitleeswaarden. De bias zal dus in iedere andere frame die getrokken wordt met de CCD aanwezig zijn en er zal dus altijd eerst voor de bias gecorrigeerd moeten worden. Het is zo dat een bias frame makkelijk te meten is, maar ook zo dat je met 1 bias frame niet voldoende hebt en je een statistisch gemiddelde moet nemen over verschillende frames. Het is dus best om een masterbias te maken door uit te middelen over een aantal frames en deze masterbias dan te gebruiken voor de verdere datareductie. Dark frame De tweede frame die normaal gebruikt wordt, is de dark frame. De dark current is een effect dat veroorzaakt wordt door elektronen gecreëerd in het kristalrooster van de silicium detector en dus ook gedetecteerd zal worden, want er kan geen onderscheid gemaakt worden door elektronen die ontstaan zijn door de absorptie van fotonen en thermische elektronen in het materiaal. Deze thermische elektronen zullen dan de dark current vormen. Net als de bias frames worden deze dark frames ook gemeten met een gesloten sluiter, maar met een van 0 verschillende belichtingstijd. De belichtingstijd van een dark frame moet normaal even lang zijn als deze van de foto’s die je ermee wil corrigeren, maar soms worden langere dark frames gemaakt door een lineaire toename te veronderstellen. De dark frames zijn dus een maat voor de thermische ruis van een detector. De dark current kan sterk onderdrukt worden door de temperatuur van de CCD te verlagen zodat er minder elektronen in het kristalrooster geproduceerd worden. Dark frames kunnen dus vermeden worden door de telescoop te koelen, wat bij de VLT en de HST inderdaad het geval is. Beide zijn ongeveer gekoeld tot -80◦ C ( ≈ 190K) en onze foto’s hebben relatief korte belichtingstijden, dus het is niet nodig om deze dark frames in rekening te brengen.

18


Flatfield Een flatfield wordt gebruikt om te corrigeren voor de pixel tot pixel variatie van de CCD. Dit is nodig omdat niet iedere pixel van de camera even gevoelig is voor het licht dat het opvangt. Voor een flatfield is het dus noodzakelijk dat de CCD door een uniforme bron belicht wordt. Deze flatfields worden daarom meestal tijdens de schemering genomen van een stuk aan de hemel waar nog zo goed als geen sterren zichtbaar zijn, zodat de camera uniform belicht is. Ook hier is het weer aangewezen om meerdere flatfields te nemen en deze dan terug te combineren tot één master flatfield en dit dan te gebruiken voor de verdere datareductie.

Kosmische effecten De hierboven beschreven effecten zijn van technische oorsprong, ze zijn intrinsiek aan de gebruikte apparatuur. Er zullen ook nog kosmische effecten optreden bij het maken van een foto. Een eerste effect is de kosmische straling. Deze hoogenergetische straling kan voor een heel heldere pixel zorgen als deze daarop invalt, wat zal te zien zijn op de foto als een ”witte”pixel. Nadat een camera een tijdje gebruikt wordt, zullen er pixels beschadigd zijn door deze hoogenergetische straling, wat we dode pixels noemen, met andere woorden dus een kapotte pixel waar we geen informatie meer kunnen uithalen. Deze beide effecten van kosmische straling kunnen verholpen worden door verschillende foto’s van het waar te nemen object te nemen en deze te combineren tot 1 foto. Het is hoogst onwaarschijnlijk dat een kosmisch deeltje met hoge energie 2 keer op net dezelfde plaats invalt op de foto. Dus door de mediaan van 3 of meer foto’s te nemen elimineer je deze hoogenergetische pixels en krijg je een foto zonder dat de kosmische straling nog zichtbaar is. Om te corrigeren voor de dode pixels kan eenzelfde principe gebruikt worden. Door in iedere foto de plaats van de galaxie op de CCD een klein beetje te veranderen zodat niet telkens hetzelfde punt van de hemel op dezelfde pixel terecht komt, kan je ook door een mediaan te nemen van de foto’s, de dode pixels elimineren. Een tweede effect waar we willen voor corrigeren is de achtergrond. De uiteindelijke doelstelling van de datareductie is dat we enkel het licht van onze galaxie zien en geen achtergrondlicht erbij hebben. Dit kan vrij makkelijk opgelost worden door een achtergrond te fitten door middel van eerste of tweede orde polynoom en deze dan van de foto af te trekken.

2.4.2

ESO-MIDAS

Om alle bovenstaande correcties door te voeren op de foto’s zullen we het programma ESO-MIDAS gebruiken. ESO-MIDAS is een programma van ESO ontwikkeld in de jaren ’80 om de bovenstaande datareductie uit te voeren. Het programma kan gerund worden in de terminal en is gebaseerd op commando’s die je moet ingeven. Het programma kan onder andere headers van fits files lezen en weergeven en de afbeeldingen weergeven in een display. Verdere informatie over wat het programma allemaal doet zal gegeven worden als we het gebruiken tijdens onze datareductie die hieronder beschreven wordt of kan gevonden worden op de ESO website (http://www.eso.org/sci/software/esomidas/).

2.5. DATAREDUCTIE FORS

2.5

19

Datareductie FORS

Eerst zullen we de datareductie van de FORS-data bespreken. Een ruwe onbewerkte foto van de geobserveerde galaxie is weergegeven in figuur 2.9. We zien dat de foto in 4 kwadranten verdeeld is. Dit komt doordat de CCD via 4 verschillende meetcircuits uitgelezen wordt. Voordeel hiervan is dat de CCD sneller wordt uitgelezen, anderzijds hebben de 4 verschillende meetcircuits wel een verschillende gain en verschillende bias waarmee er dus rekening moet gehouden worden. De gebruikte data is verzameld geweest op de nachten van 10 tot en met 14 juli 1999, en op 13 juli werd de galaxie waargenomen. De reden dat we data van deze 5 nachten gebruiken is omdat we alle standaardsterren die in deze nachten geobserveerd zijn later nodig hebben om de kalibratie van het fotometrisch nulpunt te doen. We kunnen met behulp van de biassen en de flatfields van deze 5 nachten een masterbias en master flatfield maken en deze gebruiken voor de verdere datareductie omdat beide over een korte periode constant blijven. De CCD wordt ook niet telkens met dezelfde gain uitgelezen, per uitleescircuit kan er een lage of hoge gain gebruikt worden, maar voor de 4 verschillende meetcircuits wordt wel telkens dezelfde gebruikt, dus ofwel gebruiken de 4 uitleescircuits hun lage gain ofwel hun hoge gain. Deze zijn weergegeven in tabel 2.2 .

Figuur 2.9: Een ruwe onbewerkte foto van de galaxie IC3328, genomen met het FORS instrument.

20


Tabel 2.2: De lage en hoge gain voor de 4 verschillende circuits uitgedrukt in ADU per elektron. Circuit(kwadrant) 1 2 3 4

2.5.1

Lage gain 0.31 0.32 0.34 0.29

Hoge gain 0.61 0.64 0.68 0.56

Reductie van de foto’s

Bias Eerst splitsen we al onze bias frames op in diegene die gebruik maken van een lage gain en diegene die gebruik maken van een hoge gain om dan voor beide een masterbias te maken. Dit gebeurt door de mediaan te nemen in een bepaalde pixel over alle frames. In MIDAS kan dit gedaan worden door een cataloog aan te maken die een reeks biasframes bevat. Hiervan kan dan in iedere pixel de mediaan berekend worden en deze waarden worden dan samengevoegd in een masterbias. Deze frame zal van alle verder gebruikte foto’s moeten afgetrokken worden. In figuur 2.10 is een bias frame weergegeven. We zien hier ook duidelijk een onderscheid tussen de 4 verschillende uitleescircuits van de camera. In de onderste figuur is er een plot gemaakt van een rij in het derde en vierde kwadrant, rij 500 van de frame en deze is weergegeven door de zwarte lijn. De rode lijn geeft een plot weer van rij 1500 van de frame die het tweede en eerste kwadrant snijdt. Hier is duidelijk te zien dat de bias van de 4 kwadranten lichtjes verschilt. Flatfield Zoals eerder vermeld, wordt het flatfield gebruikt om een correctie door te voeren voor het feit dat niet alle pixels even gevoelig zijn. We beschikken over een aantal flatfields die genomen zijn gedurende de 5 nachten van de observaties en deze kunnen we wederom combineren tot een master flatfield, maar voordat we de mediaan nemen over al deze flatfields gaan we ze eerst wat moeten bewerken. We beginnen met van iedere flatfield de masterbias af te trekken. Omdat de randen van de frame niet altijd even goed zijn, worden deze er afgeknipt (nadat de frames gecorrigeerd zijn met de bias). We zullen dus ook voor onze foto’s van de standaardsterren en van onze galaxie er een zelfde rand van af moeten knippen voordat we de flatfield correctie doorvoeren. Later zullen we onze foto’s van de waargenomen standaardsterren en galaxie delen door het master flatfield. Als we er voor zorgen dat de pixelwaarden van onze flatfield in ADU/elektronen staan, dan zullen onze finale foto’s direct in elektronen/s staan, wat we uiteindelijk ook willen bereiken. Om dit te doen wordt normaal de mediaan F van alle pixels berekend, waarbij F een flux voorstelt die uitgedrukt is in ADU. Deze wordt dan gedeeld door de Gain G van de telescoop en nadat de frame dan gedeeld wordt door f = F/G staat deze in ADU/elektronen. Voor onze foto’s is het iets moeilijker omdat we 4 uitleescircuits hebben. Wat we dus doen in de plaats is in ieder kwadrant de mediaan van de flux berekenen. We noemen deze F 1, F 2, F 3 en F 4. Vervolgens worden ze alle 4 gedeeld door de gain in elk kwadrant, die we kunnen voorstellen door G1, G2, G3 en G4. Dan hebben we voor alle 4 de kwadranten een flux (f 1, f 2, f 3 en f 4)


21

(a)

(b)

Figuur 2.10: (a) Voorbeeld van een bias frame en (b) de pixelwaarden uitgelezen op rij 500 (zwart) en op rij 1500 (rood). De effecten in de bias van de 4 verschillende uitleescircuits zijn op beide figuren duidelijk te zien.

die uitgedrukt is in elektronen. Ten slotte kunnen we de frame delen door het gemiddelde f = (f 1+f 2+f 3+f 4)/4 van deze 4 fluxen en zo hebben we een flatfield die ook uitgedrukt wordt in ADU/elektronen. Bovenstaande procedure wordt toegepast op ieder flatfield en daarna kunnen we, net zoals we bij de biassen gedaan hebben, een cataloog aanmaken van deze flatfields en een master flatfield samenstellen door de mediaan te nemen over alle frames in de cataloog. In figuur 2.11 wordt het masterflat field weergegeven voor foto’s met de lage gain. De donkere vlekken wijzen duidelijk op gebieden met lagere efficiëntie.

22


Figuur 2.11: Foto van het uiteindelijk bekomen masterflat field voor de lage gain. De donkere gebieden zijn gebieden waar de efficiëntie van de CCD beduidend lager ligt dan gemiddeld.

Kalibratie standaardsterren en galaxie Nu we onze masterbias en master flatfield bepaald hebben, kunnen we eens kijken naar de foto’s van de standaardsterren en de galaxie zelf. Net zoals voor het flatfield starten we voor deze foto’s ook met het aftrekken van de bias van de frames. Hierna moet gedeeld worden door het flatfield, maar omdat we het flatfield bijgeknipt hadden, moeten we dit ook doen voor deze foto’s, dus eerst de randen eraf knippen en dan delen door het flatfield. Nu staan de pixelwaarden van deze foto’s in elektronen. Ten slotte delen we de foto’s nog door hun belichtingstijd zodat de pixelwaarden nu uitgedrukt worden in elektronen per seconde. Nu hebben we onze gereduceerde foto’s van de standaardsterren en de galaxie. Maar als we magnituden berekenen, willen we er zeker van zijn dat we enkel het licht van de sterren/galaxie zien en dat de achtergrond uit onze foto’s weg is. Dit kunnen we doen door een polynoom te fitten aan de achtergrond, wat eenvoudig kan in MIDAS door het selecteren van ”lege”gebieden (waar geen sterren of galaxieën te zien zijn) op je foto. Het programma fit dan een eerste of tweede graadspolynoom (zoals je zelf kiest) aan de achtergrond en trekt deze van je foto af. In figuur 2.12 is dezelfde foto afgebeeld als in figuur 2.9, maar waarop alle bovenstaande stappen van de datareductie toegepast zijn.


23

Figuur 2.12: Een volledig gereduceerde foto van IC3328. Dit is dezelfde foto van in figuur 2.9 waarop alle processen van de datareductie uitgevoerd zijn.

2.5.2

Bepaling van het fotometrisch nulpunt

Nu we de gereduceerde afbeeldingen hebben van onze galaxie, kunnen we bijna beginnen met de analyse van de afbeeldingen. Enkel het fotometrisch nulpunt moet nog bepaald worden. Wat we willen weten is hoeveel licht de galaxie uitzendt. Om dit te bepalen zullen we gebruik maken van de standaardsterren die zich in de standaardvelden bevinden en waarvan de magnituden in verschillende kleurenbanden gekend zijn. De magnitude van zo’n ster wordt bepaald door de vergelijking, X X mX = −2.5 log f + mX 0 + kc (V − R) − ke z − AX ,

(2.1)

waarin mX de magnitude van de ster is in de X-band, f het aantal elektronen/sec(/pixel), X mX efficiënt is (Hier gegeven voor 0 het fotometrisch nulpunt in de X-band, kc de kleurco¨ X het (V − R) kleur), ke de extinctiecoëfficiënt, z de luchtmassa en AX de interstellaire absorptie. Eénmaal het fotometrisch nulpunt berekend is, kunnen we dan ook magnituden berekenen van het geobserveerde object. De standaardvelden zijn zodanig gekozen dat de interstellaire absorptie verwaarloosd kan worden. De vergelijking kan dus geschreven worden als, X X mX = −2.5 log f + mX (2.2) 0 + kc (V − R) − ke z. Normaalgezien worden kcX en keX bepaald samen met de waarnemingen, maar in de eerste meetperiode van FORS1 was dit niet het geval. We zullen deze 2 parameters dus samen

24


met het fotometrisch nulpunt fitten aan de waarnemingen van de standaardsterren. Dit is ook de reden waarom er data van 5 verschillende nachten gebruikt werd. Om genoeg data te hebben om deze parameters te fitten hebben we ook standaardvelden nodig, waargenomen bij luchtmassa’s die genoeg van elkaar verschillen (anders zou z als een constante kunnen beschouwd worden tijdens het fitten). Na het bepalen van deze 3 parameters zullen we dan voor een pixel in de foto van onze galaxie de magnitude(/pixel) kunnen uitrekenen en zo een oppervlaktehelderheidsprofiel kunnen maken. Zoals eerder vermeld, gebeurt bovenstaande kalibratie met behulp van standaardvelden. Een voorbeeld van zo een gebruikt veld is weergegeven in figuur 2.13. De linkse figuur is een afbeelding die je op de ESO website kan terugvinden. Op deze foto staan een aantal sterren aangeduid waarvan de magnituden in verschillende banden gekend zijn. De rechtse figuur is de waargenomen foto van hetzelfde veld sterren, waarin zonder veel moeite de aangeduide sterren (A,1,2,3) terug gevonden kunnen worden. In totaal werden er zo’n zestal standaardvelden gebruikt waarin telkens 3 tot 7 sterren gebruikt kunnen worden voor de bepaling van het fotometrisch nulpunt. Er werden 13 foto’s van deze standaardvelden gebruikt, wat goed was voor 53 sterren waaraan we onze 3 parameters uit vergelijking 2.2 kunnen gaan fitten.

(a)

(b)

Figuur 2.13: (a) Een foto van de ESO website waarop de standaardsterren afgebeeld c ESO) en (b) een gereduceerde foto van hetzelfde veld genomen met het FORS staan ( instrument. Nu moeten we eerst voor iedere standaardster de flux bepalen, die uitgedrukt staat in elektronen/seconde. Een probleem is hierbij dat je bij het uitlezen van de CCD geen hogere ADU dan 65536 kan krijgen, wat met behulp van de gain zal zorgen voor een maximaal aantal elektronen per seconde die je kan tellen. Dit zal er voor zorgen dat de helderste sterren gesatureerd zijn. In figuur 2.14 wordt het verschil getoond tussen een gesatureerde ster en een niet-gesatureerde ster. Enkel de sterren die niet gesatureerd zijn kunnen gebruikt worden voor de kalibratie. Als we deze sterren gevonden hebben, dan kunnen we deze uitknippen uit de foto van de standaardvelden en eventueel correcties gaan doorvoeren voor andere objecten die in de buurt van de ster liggen. In figuur 2.15


25

is een foto van een standaardster afgebeeld waarbij een ander object gemaskt werd, zodat enkel het licht van de ster gebruikt wordt bij de bepaling van het fotometrisch nulpunt. Nu we voor iedere ster de flux in elektronen per seconde bepaald hebben, kunnen we een eenvoudig programma schrijven in python dat de parameters kcX , keX en mX 0 uit vergelijking 2.2 bepaalt door gebruik te maken van de optimize.leastsq methode van scipy. De bekomen waarden voor de parameters zijn, kcX = 0.05226, keX = 0.08518 en mX 0 = 27.66096.

(a) gesatureerde ster

(b) niet-geastureerde ster

Figuur 2.14: Plot van de pixelwaarden doorheen een gesatureerde ster (a) en een niet gesatureerde ster (b).

(a) standaardster

(b) ster met gemodifieerde oppervlakte

Figuur 2.15: Foto van een (a) standaardster en (b) een standaardster waar een nabijgelegen object (aan de hemel) gemaskt werd.

26


Tabel 2.3: De flux van de standaardsterren, de V-R magnitude, de V-band magnitude en de luchtmassa z van de foto waarop de ster staat. Al deze sterren kunnen gebruikt worden om een fit te maken aan (2.2). f(e/s) SA 110-1 2.85227e6 141415 99510.4 1.65438e6 PG1633-2 170585 135004 1.19285e6 973684 485351 PG1323-2 420973 693660 392109 SA110-3 2.80313e6 1.63639e6 138052 97722.5 740481 1.60739e6 PG1323-1 414578 692493 387848 PG1633-1 170726 136837 1.21204e6 978186 485504

2.6

V-R

V

z

0.538 0.918 1.066 1.360

12.018 15.693 16.252 13.470

1.1455 1.1455 1.1455 1.1455

-0.093 0.505 0.590 0.618 0.324

14.397 15.256 12.969 13.229 13.691

1.247 1.247 1.247 1.247 1.247

-0.048 0.426 0.395

13.481 13.406 14.003

1.045 1.045 1.045

0.538 0.361 0.918 1.066 0.697 1.360

12.018 12.425 15.693 16.252 13.615 13.470

1.124 1.124 1.124 1.124 1.124 1.124

-0.048 0.426 0.395

13.481 13.406 14.003

1.1555 1.1555 1.1555

-0.093 0.505 0.590 0.618 0.324

14.397 15.256 12.969 13.229 13.691

1.213 1.213 1.213 1.213 1.213

f(e/s) SA110-2 138370 98040.2 738426 1.52093e6 MarkA-1 445571 62777.9 214188 205177 PG2213-1 210422 325700 1.27583e6 138389 PG1633-3 170344 136296 1.19387e6 972590 484189 PG2331-1 86855.1 937161 210739 PG1323-3 419793 687231 399924 PG1323-4 417620 693814 391372

V-R

V

z

0.918 1.066 0.697 1.360

15.693 16.252 13.615 13.470

1.111 1.111 1.111 1.111

-0.115 0.367 0.379 0.587

13.258 15.911 14.540 14.818

1.982 1.982 1.982 1.982

-0.092 0.406 0.427 0.426

14.124 14.178 12.706 15.109

1.48 1.48 1.48 1.48

-0.093 0.505 0.590 0.618 0.324

14.397 15.256 12.969 13.229 13.691

1.231 1.231 1.231 1.231 1.231

-0.012 0.419 0.481

15.182 13.051 14.744

1.297 1.297 1.297

-0.048 0.426 0.395

13.481 13.406 14.003

1.052 1.052 1.052

-0.048 0.426 0.395

13.481 13.406 14.003

1.05 1.05 1.05

Datareductie WFC

Foto’s van galaxieën (of andere objecten) die met een camera van de Hubble ruimtetelescoop getrokken werden zijn terug te vinden in het HST archief (zoals eerder vermeld). Deze foto’s hebben het voordeel dat de datareductie al uitgevoerd is, dus in het archief kan je de gereduceerde foto’s terugvinden, alsook de andere files die gebruikt werden voor deze kalibratie. De basisprincipes voor de kalibratie zijn dezelfde als hierboven beschreven, maar een stuk ingewikkelder om uit te voeren en ook andere niet-triviale correcties voor de detectoren in de HST moeten in rekening gebracht worden. Daarom worden de ruwe foto’s

2.7. GEREDUCEERDE DATA EN OPMERKINGEN

27

al door een pipeline gestuurd die de kalibratie voor ons doet zodat er kant-en-klare foto’s terug te vinden zijn in het archief. Op figuur 2.16 is zo’n foto van IC3328 weergegeven. De volledige foto is iets groter, maar op de figuur wordt enkel het relevante stuk met de door ons onderzochte galaxie getoond. De zwarte band te zien op de figuur is te wijten aan het feit dat de WFC van de Hubble telescoop uit 2 CCD’s bestaat die samen 1 foto vormen.

Figuur 2.16: Foto van de galaxie IC3328, die te downloaden is in he HST archief en die gebruikt kan worden voor de data analyse.

2.7

Gereduceerde data en opmerkingen

Standaardsterren controleren Het is belangrijk om alles wat je gebruikt te controleren, zelfs de data/info die je op de ESO website vindt. Omdat we data gebruiken van een paar nachten hebben we 4 foto’s getrokken in verschillende nachten van het standaardveld PG1323-086. Hiervoor werd een grote variabiliteit gevonden voor ster A in dit veld gedurende deze nachten. Na wat opzoekingswerk bleek dat deze ster eigenlijk een eclipserende dubbelster is met een korte periode zoals in 2009 beschreven werd door Clem en Landolt. 3 jaar later staat deze echter nog steeds op de standaardvelden, voorhanden door ESO, aangeduid. Het is bijgevolg niet mogelijk om deze ster te gebruiken als standaardster voor kalibratie. Verder is de aanduiding op de template van hetzelfde standaardveld de aanduiding van PG 1323086 onduidelijk. De naam van iedere ster staat er rechts van aangeduid op de kaarten, deze staat echter links zodat als er een hoop sterren te zien zijn in het standaardveld, wat ook enige verwarring teweeg kan brengen. Een goede controle van de standaardsterren is dus nodig om een goede kalibratie te hebben!

28


Astero¨ıde gespot Als je een stukje van de hemel fotografeert, dan is er natuurlijk veel meer te zien op je foto dan enkel het object dat je wilt waarnemen. Zo is in figuur 2.17 te zien dat er zich een object beweegt linksonder de galaxie. Na wat opzoekingswerk bleek dit een astero¨ıde te zijn uit de Planeto¨ıdengordel tussen Mars en Jupiter. De baan van de astero¨ıde in ons zonnestelsel wordt weergegeven in figuur 2.18.

(a)

(b)

(c)

(d)

Figuur 2.17: De astero¨ıde beweegt van noordoost naar zuidwest op de foto’s. Foto (a) is eerst getrokken, dan (b) en dan (c). Op (d) is een combinatie van de 3 foto’s te zien en de astero¨ıde doet zich daar voor als een lijn, aangeduid door de rode pijl, en valt makkelijk op tijdens het bekijken van de gereduceerde data. De kleuren op deze foto’s zijn verdeeld op een logaritmische schaal.

c NASA) Figuur 2.18: De baan van de asteroide 2000 RS 36 in ons zonnestelsel. (

Hoofdstuk 3

Analyse van de data 3.1 3.1.1

Model fitten aan galaxie Schijnbare vorm van de galaxie

De schijnbare vorm van de galaxie is de vorm die we waarnemen aan de hemel, dus dit is eigenlijk een 2 dimensionale projectie van de galaxie. Voor elliptische galaxieën kan een schijnbare ellipticiteit gedefinieerd worden als = 1 − b/a, met a en b respectievelijk de schijnbare halve grote en halve kleine as van een isofoot (lijn met constante oppervlaktehelderheid). Een projectie van zo’n elliptische galaxie kan nooit ronder zijn dan de oorspronkelijke galaxie, maar wel platter. Nu is het de bedoeling dat we een model maken van de projectie van IC3328, door de isofoten te fitten, om later hieruit het oppervlaktehelderheidsprofiel te halen. In eerste instantie kunnen ellipsen gefit worden aan deze isofoten, maar er kunnen zich afwijkingen voordoen van deze ellipsen zodat de isofoten meer schijfvormig of doosvormig zijn. Een meer nauwkeurige beschrijving kan dan gegeven worden door de volgende vergelijking, I(θ) = I0 × 1 + C3 cos(3θ) + C4 cos(4θ) + S3 sin(3θ) + S4 sin(4θ) , (3.1) waarbij I(θ) de intensiteit volgt langs een ellips met een bepaalde assenverhouding b/a, I0 de gemiddelde intensiteit is van de isofoot en θ de hoek is, gemeten vanaf de grote as a. De coëfficiënten Ci en Si met i = 1, 2 mogen achterwege gelaten worden omdat deze het centrum van de ellips bepalen en aangeven hoe elliptisch het profiel is en aangezien er al elliptische isofoten gefit worden, zijn deze coëfficiënten nutteloos. De coëfficiënten met i = 3, 4 geven een maat voor de afwijkingen tegenover de perfecte ellips, deze worden dus wel behouden. Zo kunnen de isofoten er doosvormig uitzien of eerder schijfvormig. Verder

Figuur 3.1: De positiehoek gedefinieerd van noord (N) naar oost (E). 29

30

HOOFDSTUK 3. ANALYSE VAN DE DATA

is het ook zo dat de isofoten niet perfect concentrische ellipsen zijn. Dit kan aangeduid worden met de positiehoek, die gedefinieerd is van noord naar oost zoals op figuur 3.1 te zien is en de hoek aanduidt tussen de grote as a en een referentierichting.

3.1.2

Code

Nu kan met behulp van een fortran code gebaseerd op bovenstaande principes een model voor de isofoten gefit worden aan de galaxie. Deze code is dezelfde als gebruikt door De Rijcke et al. (2003) en werd ook door hen geschreven. Dit programma bestaat uit 5 afzonderlijke bestandjes die dienen gerund te worden. Het eerste maakt niet echt deel uit van het programma zelf, maar zorgt ervoor dat de juiste keywoorden in de header van de gebruikte foto’s staan, die als fits-bestand dienen meegegeven te worden aan het programma. De waarden die naar de header van de fits files geschreven worden zijn die van de extinctiecoëfficiënt keX , de belichtingstijd en de luchtmassa z. Vervolgens wordt het programma fit dat de isofoten gaat fitten aan de foto van de galaxie gebruikt. Het centrum van de galaxie wordt meegegeven en het programma start met het fitten van een isofoot voor een bepaalde oppervlaktehelderheid en gaat dan in stapjes ∆x die meegegeven worden in de input file steeds isofoten fitten met een oppervlaktehelderheid die ∆x lager is dan de vorige isofoot. Dit blijft duren tot een isofoot gefit wordt voor de ondergrens in magnitude/arcsec2 die ook meegegeven wordt in de input file. Dit programma schrijft al zijn gegevens in een bestand. Hierin staan onder andere de coördinaten van het centrum van de ellips, de halve grote as a en halve kleine as b, de positiehoek en de Fouriercoëfficiënten die de afwijking van de isofoten op een pure ellips karakteriseren. Het volgende programma dat gebruikt wordt is clean. Dit filtert alle isofoten uit de output file van het programma fit waarvoor de a of b kleiner zijn dan de a en b van de vorige isofoot. Dit is nodig om het 3de programma recon dat alle fotometrische grootheden zoals oppervlaktehelderheid, halve kleine as, positiehoek, ellipticiteit en Fouriercoëfficiënten in functie zet van de halve grote as a met behulp van kubieke splines. Dit gebeurt zodat vrij gladde waarden bekomen worden voor deze grootheden. Al deze grootheden worden dan in aparte data files opgeslagen. Deze data files worden dan gebruikt door het programma dpmd om een gladde versie van de oorspronkelijke foto te maken. Input file Bovenstaande programma’s maken gebruik van een input file. Hieronder is de input file voor de VLT data gegeven, de nummers in de eerste kolom dienen enkel om te weten op welke regel wat moet ingegeven worden. 1 2 3 4 5 6 7 8 9

vcc0856 R.fits vcc0856 R filt.fits mask.fits 3 400 41.7 689 1293 0.2

3.1. MODEL FITTEN AAN GALAXIE 10 11 12 13 14 15 16 17

31

0.05 out.fits 1 0.0852 5 1.4 24.166 27.5

Op de eerste 3 regels staan 3 fits files die moeten meegegeven worden aan de programma’s. De eerste file is een volledig gereduceerde foto van de waargenomen galaxie en de 2de file is een gefilterde foto. Deze wordt gemaakt door in een gebied van 10 pixels rond iedere pixel de mediaan te nemen van alle pixels in dit gebied en deze waarde toe te kennen aan de pixel zelf. Deze file wordt gebruikt omdat, indien je bij de eerste foto naar lagere magnitudes zou gaan om isofoten te fitten aan de oppervlaktehelderheid, er te veel ruis zit op de foto terwijl de ruis onderdrukt wordt in de gefilterde foto zodat er voor de iets lagere magnitudes toch nog isofoten kunnen gefit worden. Deze 2 foto’s van zowel de VLT als HST data worden getoond in figuur 3.2. De 3de file is een mask file. Deze foto’s worden voor beide telescopen weergegeven in figuur 3.3(c) en 3.3(d) en is op een bepaalde plaats ofwel 0, ofwel 1. Het is de bedoeling van deze file om voorgrondsterren en andere galaxieën in de buurt van het geobserveerde object te ’masken’. Als het programma de isofoten fit, zal het geen rekening houden met de gemaskte gebieden. Dus effecten van andere heldere objecten zullen opgegeven worden bij het reconstrueren van de isofoten. In figuur 3.3 worden 3 files weergegeven. Op de bovenste rij staan de oorspronkelijk gereduceerde foto’s, op de tweede rij zijn zoals eerder vermeld de mask files die 0 zijn in de donkere (zwarte) gebieden en 1 zijn in de heldere (witte) gebieden. De derde rij is het product van de eerste 2 foto’s. Bij een vergelijking met de eerste foto is duidelijk dat de helderste objecten gemaskt zijn. Op de 4de regel van de input file staan het aantal foto’s die gebruikt werden om de uiteindelijke foto te maken, wat in ons geval 3 is. De 5de en 6de regel geven respectievelijk de belichtingstijd weer van 1 foto en de achtergrond van 1 foto. Dit is de achtergrond voordat deze gefit en afgetrokken werd op het einde van de datareductie. Deze achtergrond is dus diegene die je nog op je foto hebt nadat je de correctie van de bias en flatfields doorgevoerd hebt en voordat je er een eerste- of tweedegraadspolynoom aan gaat fitten. De 7de en 8ste regel zijn de co¨ ordinaten van het midden van de galaxie, in pixelcoördinaten. De eerste is de x-co¨ ordinaat en de tweede de y-coördinaat. Op de 9de regel staat de pixelschaal, of met andere woorden de grootte van een pixel, wat voor de VLT 0.2 boogseconden is. Op regel 10 staat de stapgrootte dat het programma moet nemen om zijn isofoten te fitten. Dus tussen 2 gefitte isofoten zal een verschil in magnitude van 0.05 zitten. Op de 11de regel staat de naam van de output file die het model van de galaxie bevat na het fitten van de isofoten. Op de 12de regel staat de gain van de telescoop. In ons geval is deze al gecorrigeerd voor de gain en staan de input files al in elektronen per seconde. Deze waarde zal dus verschillen van 1 als de gebruikte foto’s (die op de eerste 3 regels staan in de input file) nog in ADU per pixel staan. De 13de, 14de en 15de regel bevatten respectievelijk de extinctiecoëfficiënt, de read out noise en de luchtmassa van de foto’s, waarin

32


(a) VLT foto

(b) HST foto

(c) gefilterde VLT foto

(d) gefilterde HST foto

Figuur 3.2: De originele gereduceerde foto van de VLT (a) en de HST (b) en de gefilterde foto’s van beide (c) en (d). de extinctiecoëfficiënt deze is bepaald door de eerder vermelde fit aan vergelijking 2.2. Op de 2 laatste regels staan de nulpuntsmagnitude en de magnitude waarbij het programma stopt met het fitten van isofoten, dus de laagste magnitude waarvoor het programma een isofoot zal fitten. De nulpuntsmagnitude wijkt af van mX 0 = 27.661 gegeven door de fit aan vergelijking 2.2 . Dit komt omdat de nulpuntsmagnitude op regel 16 van de input file uitgedrukt is in magnitude/arcsec2 en deze zal dus gelijk zijn aan: m0 = mX 0 + 5 log(pixelgrootte),

(3.2)

waarin de pixelgrootte de grootte van een pixel in boogseconden is, wat voor de FORS detector van de VLT 0.2 boogseconden is. Dit invullen levert dat m0 = 24.166, wat teruggevonden wordt in de input file. Het gebruik van dit nulpunt zorgt ervoor dat de berekende magnitudes ook direct in magnitudes per vierkante boogseconde zullen staan. De ondergrens die weergegeven is op regel 17 staat dus logischerwijs ook uitgedrukt in magnitude/arcsec2 .

3.1. MODEL FITTEN AAN GALAXIE

33

(a) VLT foto

(b) HST foto

(c) VLT mask file

(d) HST mask file

(e) VLT masked

(f) HST maked

Figuur 3.3: De originele, gereduceerde foto’s (a), (b) en de mask files die gebruikt worden om de heldere objecten te verduisteren (c), (d). In (e) en (f) is het product van respectievelijk (a)×(c) en (b)×(d) weergegeven. De kleurverdeling is in deze foto’s op een logaritmische schaal verdeeld.

34

3.2


Resultaten van het model

Na het fitten van het model voor de galaxie, kunnen we de waarden plotten voor alle parameters die het programma berekend heeft. Eerst kunnen we een blik werpen op het gefitte model van de galaxie. De modellen voor beide data (VLT en HST) zijn weergegeven in figuur 3.4 op de middelste rij. De eerste rij geeft de oorspronkelijk gereduceerde foto’s weer en in de laatste rij zijn de foto’s te zien waarbij het model van de oorspronkelijke foto’s afgetrokken is. Zoals op het eerste zicht aan deze foto’s te zien is, lijkt het model geslaagd. Nu kunnen we de andere parameters van de isofoten gaan bekijken. In figuur 3.5 en 3.6 zijn respectievelijk de parameters van de VLT data en HST data geplot in functie van de straal r, die uitgedrukt staat in boogseconden. Zoals in de tweede en derde figuur (van boven naar onder) in de linkerkolom van beide figuren te zien is, is er op het eerste zicht niks abnormaals aan de data. Deze duiden op mooie concentrische ellipsen voor de eerste 30 à 40 boogseconden, maar de plot van de ellipticiteit en de positiehoek vertoont abnormaal gedrag. De coherente wiebelingen duiden op een onderliggende structuur. Die structuur blijkt een spiraalstructuur te zijn zoals aangetoond werd door Jerjen, Kalnajs en Binggeli (2000). Uit de plot kunnen we zien dat het spiraalpatroon zich uitstrekt tot ongeveer 30 à 40 boogseconden. De aandachtige lezer ziet dat het oppervlaktehelderheidsprofiel (bovenste plot in figuur 3.5 in de linkerkolom) een knik vertoont rond deze straal. Dit zal zijn omdat het oppervlaktehelderheidsprofiel uit een component bestaat die afkomstig is van de schijf, die ook de spiraalstructuur bevat, en een component die afkomstig is van de halo van de galaxie. Zoals verder beschreven wordt (sectie 3.3), zullen we een Sérsicprofiel fitten aan dit helderheidsprofiel. Het totale profiel zal dus een som zijn van de bijdragen van de verschillende componenten. Naast de bijdrage van de schijf en de halo zal er ook nog een bijdrage van de nucleus zijn. Dit is te zien aan de piek die voorkomt bij de kleinste stralen (kleiner dan 200 ). Voor de iets grotere stralen (200 tot ongeveer 4000 ) zal de schijf de belangrijkste bijdrage leveren aan het helderheidsprofiel en voor nog grotere stralen zal de halo de overhand nemen. Om de structuur en het profiel van de nucleus aan het licht te brengen zullen we gebruik moeten maken van data met een veel hogere resolutie dan diegene die we kunnen bereiken met de VLT telescoop. De beperkende factor voor deze telescoop is namelijk de seeing van de atmosfeer, ook al was deze uitstekend tijdens de observaties (0.600 ). Dit is niet genoeg om de galaxie op te lossen voor de kleinste stralen. Hier komt de data van de HST van pas. Zoals te zien in figuur 3.6 is er data voor veel kleinere stralen beschikbaar. De pixelgrootte van de HST is 0.05 boogseconden en aangezien de telescoop zich in een baan rond de aarde bevindt buiten de atmosfeer, heeft deze geen last van de seeing van de atmosfeer. Zo zie je dat bij het oppervlaktehelderheidsprofiel dat volgt uit de HST data, de helderheid veel hoger piekt voor kleine stralen. Als de helderheidsprofielen van beide datasets bekeken worden ziet de lezer dat de magnituden van de HST data iets hoger liggen dan de magnituden van de VLT data. Dit kan veroorzaakt zijn door het gebruik van een andere magnitudeschaal en door het feit dat de waarnemingen niet met dezelfde filter gebeurd zijn. De VLT waarnemeningen zijn bij ongeveer 600 nm gedaan terwijl die bij de HST op ongeveer 900 nm gedaan zijn. Nu is het verschil tussen beide profielen een zo goed als constante waarde. Hierdoor kan de data van de HST verschoven worden zodat deze op die van de VLT ligt. Het op deze manier combineren van beide datasets zorgt voor een helderheidsprofiel met hoge resolutie

3.2. RESULTATEN VAN HET MODEL

35

(a) Oorspronkelijke foto VLT

(b) Oorspronkelijke foto HST

(c) Model VLT data

(d) Model HST data

(e) Oorspronkelijke foto met het model (f) Oorspronkelijke foto met het model ervan afgetrokken, VLT ervan afgetrokken, HST

Figuur 3.4: Op de bovenste rij zijn de foto’s van de 2 telescopen weergegeven. Op de tweede rij staan de modellen die uit het programma beschreven in sectie 3.1 volgen. Op de onderste rij staan de oorspronkelijke foto’s met de modellen ervan afgetrokken.

36


(data voor kleine stralen) in de R-band van de galaxie. De fouten, die ook weergegeven zijn op de figuren, zijn bepaald met een bootstrapmethode, die willekeurige datapunten trekt met terugplaatsing, uit de dataset en een fit maakt. Door dit herhaaldelijk te doen, wordt een fout op de data bepaald. Bepaling van de schijnbare magnitude uit de modellen voor de galaxie Uit de modellen die weergegeven zijn in figuur 3.4(c) en 3.4(d) kunnen de schijnbare magnitude in de R-band en die in de F 850LP -band (filter gebruikt op HST) bepaald worden aan de hand van vergelijking 2.1 . Het product van de kleurcoëfficiënt kcX en de kleur en de interstellaire absorptie kunnen verwaarloosd worden. Zowel de kcX als V − R zijn klein en hun bijdrage tegenover de andere termen zal verwaarloosbaar zijn en verder is er in de gezichtslijn met de galaxie zo goed als geen interstellaire absorptie, dus deze term mag ook achterwege gelaten worden. We krijgen dus als vergelijking voor de schijnbare magnitude, X mX = −2.5 log(f ) + mX (3.3) 0 + ke z. Met behulp van het ESO-MIDAS programma kunnen we de flux f , uitgedrukt in elektronen per sec, van de galaxie berekenen uit de modellen. De waarden voor alle parameters zijn weergegeven in tabel 3.1. De nulpuntsmagnitudes zijn de nulpunten op de Vega magnitudeschaal, dus de schijnbare magnitudes zijn ook weergegeven in de Vega magnitudeschaal. De bekomen waarden voor de magnitudes in de R-band en de F 850LP -band zijn respectievelijk 13.15 en 12.47. Tabel 3.1: De waarden van de parameters uit vergelijking 3.3, weergegeven voor de VLT en de HST. parameter f mX 0 keX z mX

VLT 568475 27.66096 0.085 1.4 13.15

HST 56577.2 24.34711 — — 12.47

3.2. RESULTATEN VAN HET MODEL

37

Figuur 3.5: Linkerkolom(van boven naar onder): het oppervlaktehelderheidsprofiel, afwijking van het middelpunt van de isofoot tegenover het middelpunt van de galaxy langs de kleine as, afwijking van het middelpunt van de isofoot tegenover het middelpunt van de galaxy langs de grote as, de positiehoek en het elliptisch type van de isofoten, wat gelijk is aan 10 keer de ellipticiteit. Rechterkolom: De Fouriercoëfficiënten die de afwijkingen tegenover een echte ellips beschrijven. VLT data.

38


Figuur 3.6: Linkerkolom(van boven naar onder): het oppervlaktehelderheidsprofiel, afwijking van het middelpunt van de isofoot tegenover het middelpunt van de galaxy langs de kleine as, afwijking van het middelpunt van de isofoot tegenover het middelpunt van de galaxy langs de grote as, de positiehoek en het elliptisch type van de isofoten, wat gelijk is aan 10 keer de ellipticiteit. Rechterkolom: De Fouriercoëfficiënten die de afwijkingen tegenover een echte ellips beschrijven. HST data.

3.3. FITTEN VAN EEN PROFIEL AAN DE OPPERVLAKTEHELDERHEID

3.3

39

Fitten van een profiel aan de oppervlaktehelderheid

Galaxieën zijn drie-dimensionale objecten, maar het enige wat we kunnen waarnemen is een twee-dimensionale projectie van de galaxie op de hemelsfeer. De intrinsieke grootheden van een galaxie beschrijven de eigenschappen ervan in de drie-dimensionale ruimte. Het oppervlaktehelderheidsprofiel dat we waargenomen hebben is een geprojecteerde grootheid, wat een waarneembare grootheid is aan het tweedimensionale hemeloppervlak. Het is dus de bedoeling om het helderheidsprofiel zoals bepaald in voorgaande sectie te deprojecteren zodat we de stellaire emissiviteit bekomen, welke uitdrukt hoeveel radiatief vermogen uitgestraald wordt door de sterren per kubieke meter per golflengte per ruimtehoek. Deze grootheid zal via de massalichtkrachtverhouding gerelateerd zijn aan de dichtheidsverdeling van de sterren. In deze sectie beginnen we met een profiel te beschrijven dat de waargenomen intensiteit (of oppervlaktehelderheid), wat dus de projectie is van de emissiviteit, moet beschrijven. Hierna zullen we het effect van de point-spread functie bekijken op de waargenomen data, om ten slotte een profiel te fitten aan de waargenomen intensiteit.

3.3.1

Profielen voor de intensiteit

Eén van de populairste modellen om oppervlaktehelderheden van elliptische galaxieën te fitten is de wet van Vaucouleur. Deze wet gaat uit van het feit dat de oppervlaktehelderheid een functie is van r1/4 en kan dus geschreven worden in de vorm, µ(r) = A + Br1/4 ,

(3.4)

waarin µ de oppervlaktehelderheid is uitgedrukt in magnitude, A en B te fitten constanten zijn en r de straal tot het middelpunt van de galaxie voorstelt. Een equivalente voorstelling hiervan is om het profiel uit te drukken als de intensiteit in plaats van de oppervlaktehelderheid. Deze voorstelling kan geschreven worden als volgt, r I(r) = I0 exp − ( )1/4 , (3.5) r0 waarin I(r) de intensiteit op straal r voorstelt en I0 e−1 de intensiteit voorstelt op de schaalengte r0 . Beiden zijn gerelateerd op de volgende manier, µ(r) = m0 − 2.5 log(I(r)).

(3.6)

Vergelijking 3.5 hierin invullen levert vergelijking 3.4 met A = m0 − 2.5log(I0 ) en B = −1/4 −2.5 log(e)r0 op. Voor sommige galaxieën is bovenstaande beschrijving heel goed, maar voor elliptische dwerggalaxieën blijkt een exponentiële beschrijving voor de intensiteit een beter resultaat op te leveren, r I(r) = I0 exp(− ). (3.7) r0 Bij een exponentieel model voor de intensiteit zal de oppervlaktehelderheid lineair afhangen van de straal r, wat ook te zien is in de plot van de oppervlaktehelderheid in figuur 3.5. De schijf weegt door van ongeveer 2 tot 35 boogseconden. Bij deze stralen zien we een lineair verloop van het oppervlaktehelderheidsprofiel en om het intensiteitsprofiel te

40


fitten kan dus een exponentiële gebruikt worden. Hetzelfde geldt voor de bijdrage van de halo, die doorweegt vanaf ongeveer 40 boogseconden. Deze 2 componenten kunnen we dus fitten door gebruik te maken van de som van 2 exponentiëlen. Nu kunnen we bovenstaande profielen nog uitbreiden tot een model met een algemene macht 1/n, r (3.8) I(r) = I0 exp − ( )1/n , r0 wat het Sérsicmodel is. Dit model heeft 3 parameters om te fitten aan het profiel. Dit profiel zal gebruikt worden om het oppervlaktehelderheidsprofiel van de nucleus te fitten, die duidelijk de overhand neemt voor de kleinste stralen (kleiner dan 2 boogseconden). In figuur 3.7 zijn 4 Sérsicprofielen weergeven met verschillende waarden voor de exponent 1/n. Hoe groter de exponent, hoe vlugger het profiel naar beneden zal vallen voor grotere stralen. Voor de parameters m0 , I0 en r0 werden respectievelijk de waarden 18, 20 en 2 gebruikt. Deze figuur werd gemaakt met behulp van python.

Figuur 3.7: Plot van de oppervlaktehelderheid van Sérsicprofielen waarvoor in vergelijking 3.7 verschillende waarden voor n gebruikt werden. n = 1 levert een exponentieel profiel op, terwijl n = 4 een Vaucouleursprofiel oplevert.

3.3.2

De point-spread functie

Zelfs in een optisch perfect instrument zal het licht van een puntbron niet gefocust kunnen worden in een perfect punt. De lichtgolven convergeren naar de focus van het optisch systeem en komen in bijna perfecte fase aan in de omliggende regio van de geometrische focus. Ze komen niet toe in een perfect punt, maar zullen een smalle schijf vormen, die de diffractieschijf genoemd wordt. Deze wordt weergegeven in figuur 3.8. Door het golfpatroon van licht zal een deel van het licht buiten de schijf terecht komen en een diffractiepatroon vormen, wat ook te zien is op de figuur.


41

Figuur 3.8: De voorstelling van de diffractieschijf van een optisch systeem. De ringen naast c wikipedia ) de middelste schijf worden veroorzaakt door het golfkarakter van licht. (

De distributie van het licht van een puntbron wordt dus weergegeven zoals op de figuur en zal beschreven worden door de point-spread function (PSF). De PSF kan een perfect diffractiepatroon zijn of kan onscherp zijn door optische aberraties of de seeing van de atmosfeer. De PSF beschrijft dus hoe de telescoop het beeld van een perfecte puntbron weergeeft. De hoekdiameter van de diffractieschijf geeft dus weer hoe groot de puntbron zal weergegeven worden en wordt gegeven door,

θ(radialen) = 2.44

λ , D

(3.9)

waarin λ de golflengte van het licht is en D de diameter van de telescoop. Voor de VLT (8.2m) op 660 nm en de HST (2.4m) op 900.8 nm wordt θ respectievelijk gegeven door 0.0400 en 0.18900 . Vergelijken we dit met de pixelgrootte van beide telescopen, 0.200 voor de VLT en 0.0500 voor de HST, dan zien we dat de diffractieschijf volledig binnen de pixelgrootte van de telescoop valt voor de VLT. Hierdoor zou je verwachten dat je een puntbron ook als een puntbron zal waarnemen. Dit is echter niet het geval omdat we rekening moeten houden met de seeing van de atmosfeer, welke bij de observaties ongeveer 0.600 bedroeg. Het probleem van de atmosfeer valt natuurlijk weg voor de HST. De PSF zal voor de HST dan ook door het optisch systeem bepaald worden en de binnenste piek van het diffractiepatroon zal een diameter hebben bepaald door formule 3.9 . In figuur 3.9 is te zien dat de diffractieschijf goed benaderd kan worden door een Gaussische functie. Deze Gaussische functie zal gekarakteriseerd worden door een standaardafwijking σ, die een maat zal zijn voor hoe uitgesmeerd de puntbron er zal uitzien nadat het licht ervan door het optisch systeem gepasseerd is. De bepaling van σ kan gedaan worden door een Gaussisch profiel te fitten aan de sterren die zichtbaar zijn op de foto die van IC3328 getrokken is. Dit kan eenvoudig gedaan worden met MIDAS dat de FWHM (Full width half maximum) van de gefitte Gaussische functie weergeeft.

42


Figuur 3.9: vaste lijn: Een radiale doorsnede van de diffractieschijf. stippellijn: Een c wikipedia ) Gaussische benadering voor de radiale doorsnede. (

3.3.3

Bepaling standaardafwijking van de PSF

De standaardafwijking van de PSF kan bepaald worden door een Gaussisch profiel te fitten aan een puntbron waargenomen met de telescoop. Gelukkig zijn er een heleboel puntbronnen voor handen om de PSF van een telescoop te bepalen. Alle sterren die naast de waargenomen galaxie ook in het gezichtsveld van de telescoop liggen, kunnen als puntbron fungeren. Dit is niet zo evident want puntbronnen op laboratoriumschaal zijn lang niet zo eenvoudig om te maken, een ledlampje van 1mm groot zou op een afstand van meer dan 5 km moeten staan om door een 8.2m telescoop als een puntbron waargenomen te worden, maar de gigantische afstanden tot de sterren in het gezichtveld zorgen er voor dat deze daar ideaal voor zijn.

Figuur 3.10: De pixelwaarden langsheen een doorsnede van de ster. De rode curve is een Gaussische functie gefit aan het profiel van de ster.


43

Tabel 3.2: De waarden voor de FWHM van de gefitte Gaussische verdelingen uitgedrukt in pixels. ster 1 2 3 4 5 6 7 8 gemi σgemi gem σgem

FWHM x-as FWHM y-as 3.5306 3.4678 3.4307 3.3852 2.6473 2.4622 3.5683 3.4637 2.5656 2.5759 3.0803 3.4052 2.5996 2.5372 3.0794 3.2254 3.0627 3.0653 0.1778 0.2066 3.0640 0.1794

De sterren zullen door de PSF uitgesmeerd worden op de detector van de camera zoals beschreven door de PSF. Uit figuur 3.9 bleek al dat een Gaussische functie hiervoor een goede benadering levert. Op figuur 3.10 zijn de pixelwaarden van een doorsnede van een ster te zien. De rode lijn is een Gaussische fit aan het profiel van de ster. Dit gebeurt met het ESO-MIDAS programma dat de FWHM (Full width at half maximum) weergeeft. Dit is de waarde waarop de verdeling de helft van zijn maximale waarde bereikt. We kunnen dus als profiel voor de PSF een 2 dimensionale genormeerde Gaussische functie voorstellen van de vorm, 1 y 2 1 x2 √ P SF (r) = √ , (3.10) + exp − 2σx2 2σy2 2πσx 2πσy p waarin r = x2 + y 2 de afstand is tot het centrum van het profiel en σx en σy de standaardafwijkingen van het Gaussisch profiel zijn. In tabel 3.2 zijn de FWHM, uitgedrukt in pixels, van 8 sterren weergegeven in de x- en de y-richting. Het gemiddelde van beide richtingen is weergegeven en de fout erop. Deze gemiddelde waarden liggen heel goed bij elkaar. Het feit dat deze waarden voor de FWHM in beide richtingen gelijk zijn, zorgt ervoor dat de vorm van de PSF iets makkelijker wordt. De PSF zal enkel afhangen van de straal r en niet van de richtingshoek θ, waarin r en θ poolcoördinaten zijn. De PSF kan nog geschreven worden als, P SF (r) =

1 r 2 exp − , 2πσ 2 2σ 2

(3.11)

waarin σ = σx = σy is. De gemiddelde FWHM van beide richtingen is 3.06. Nu moeten we dit nog gaan omrekenen naar σ. De FWFM van de Gaussische verdeling wordt gegeven door, (F W HM/2)2 1 02 1 exp − = exp − , (3.12) 2 × 2πσ 2 2σ 2 2πσ 2 2σ 2 waaruit volgt, p F W HM = 2 2 ln(2)σ,

(3.13)

44


of nog, F W HM σ= p . 2 2 ln(2)

(3.14)

Met de gemiddelde waarde uit tabel 3.2 kan nu de standaardafwijking van het Gaussisch profiel voor de PSF berekend worden. Hiervoor wordt een waarde σ = 1.30 bekomen. De straal in vergelijking 3.11 wordt uitgedrukt in boogseconden, dus moet σ ook in boogseconden staan. Rekening houdend met een pixelgrootte van 0.0500 vinden we σ = 0.065 boogseconden. Nu de standaardafwijking en daarmee de PSF gekend is, kunnen we het effect van de PSF op de intensiteit bekijken.

3.3.4

Convolutie van de PSF en de intensiteit

De PSF zoals hierboven beschreven, zal als resultaat hebben dat het ware intensiteitsprofiel van de galaxie zal afwijken van hetgeen waargenomen wordt door de telescoop. Beide profielen zullen door volgende convolutie aan elkaar gerelateerd zijn, Z f (r) = P SF (r0 )I(r + r0 )dr0 , (3.15) waarin f het waargenomen profiel is en I het echte profiel is. Zowel de PSF als de helderheidsprofielen zijn onafhankelijk van de hoek θ en dus enkel functie van de afstand r, waarbij θ en r poolco¨ ordinaten zijn. De vectoren r en r0 zijn weergegeven in figuur 3.11.

Figuur 3.11: Schematische voorstelling van de vectoren r, r0 en de hoek θ tussen de 2 vectoren. Het centrum van de galaxie valt samen met de oorsprong van het assenstelsel en is daarvoor logischerwijs r = 0. r is de vector die de plaats tegenover het middelpunt van de galaxie aanduidt, terwijl de vector r0 de plaats van een punt aanduidt en in bovenstaande integraal de PSF in dat punt zal weergeven tegenover r. De hoek tussen beide vectoren is θ. We kunnen voorgaande formule dus herschrijven als, Z ∞ Z 2π 0 0 f (r) = r dr dθP SF (|r0 |)I(|r + r0 |), (3.16) 0

0


45

waarin dr0 vervangen wordt door r0 dr0 dθ. Het uitwerken van de variabele van de functie I geeft als resultaat, ∞

Z

0

r dr

f (r) =

0

Z

2π

dθP SF (r0 )I

p r02 + r2 + 2rr0 cos(θ) .

(3.17)

0

0

Als we nu de uitdrukkingen voor de P SF (3.11) en I (3.8) invullen krijgen we, I0 f (r) = 2πσ 2

Z

∞

0

r dr 0

0

Z 0

2π

p 02 r + r2 + 2rr0 cos(θ) 1/n r02 dθ exp − 2 exp − . (3.18) 2σ r0

In figuur 3.12 wordt het effect van de convolutie tussen het werkelijke profiel van de galaxie en het waargenomen profiel getoond. De rode curve in de figuur is het werkelijke profiel I terwijl de groene curve het waargenomen profiel f is. Het effect is enkel beduidend op kleine stralen van het centrum van de galaxie. Het zal ervoor zorgen dat de centrale piek van het intensiteitsprofiel afgevlakt wordt. Deze wordt platter (eerst gaat de groene curve onder de rode) en breder (daarna gaat de groene curve boven de rode) en op grotere afstanden is er geen verschil meer (beide curves naderen naar elkaar). Het feit dat de piek niet alleen platter maar ook breder wordt, komt doordat de totale lichtkracht behouden blijft bij de convolutie, of met andere woorden, de oppervlakte onder beide curves blijft dezelfde. De waarden gebruikt om de plot te maken in de figuur zijn r0 = 0.0004, I0 = 141, n = 5 en m0 = 20. De figuur werd gemaakt door gebruik te maken van een pythonscript. Om de integraal uit te rekenen werd gebruik gemaakt van de scipy methode integrate.quad.

Figuur 3.12: Het effect van een convolutie op een oppervlaktehelderheidsprofiel. De piek wordt smaller en breder.

46

3.3.5


Het fitten van de data

Het oppervlaktehelderheidsprofiel van beide datasets bekomen zoals beschreven in sectie 3.2, is weergegeven in figuur 3.13. De rode datapunten zijn de datapunten van de HST terwijl de groene die van de VLT voorstellen. Op het eerste zicht is te zien dat de resolutie van de HST veel groter is dan die van de VLT. De lagere resolutie van de VLT leidt net als de PSF tot een platter en breder profiel voor de piek. Op figuur 3.13(a) is te zien dat beide datasets door ongeveer een constante van elkaar verschillen. Zoals eerder vermeld, komt dit vooral door het gebruik van verschillende filters in een ander golflengtegebied. We willen de intensiteit in de R-band bepalen en dus kunnen we het profiel van de HST verschuiven, wat weergegeven is in figuur 3.13(b). Hierop is ook te zien dat beide profielen enkel voor kleine stralen significant verschillen van elkaar, wat logisch is door de hogere resolutie van de HST. Om uit deze 2 datasets 1 set punten te halen die een zo goed mogelijke weergave van het profiel weergeven, selecteren we voor een straal kleiner dan 2 boogseconden de data van de HST, die de centrale piek van het profiel veel beter beschrijft en voor grotere stralen de data van de VLT, waarvoor we data hebben tot op ongeveer 80 boogseconden van het centrum. Deze geselecteerde data is weergegeven in figuur 3.13(c). Het is nu de bedoeling om aan deze data een som van profielen te fitten. Voor het profiel van de nucleus van de galaxie stellen we om te beginnen een Sérsicprofiel voorop en voor zowel de schijf als de halo een exponentieel profiel, r 1/nN,0 IN (r) = IN,0 exp − ( , ) rN,0 IS (r) = IS,0 exp(−

IH (r) = IH,0 exp(−

r rS,0

(3.19)

),

(3.20)

),

(3.21)

r rH,0

waarin N, S en H respectievelijk staan voor nucleus, schijf en halo. In totaal moeten dus 7 parameters gefit worden, de centrale intensiteiten van de 3 profielen, IN,0 , IS,0 , IH,0 , de 3 schaallengtes van de profielen, rN,0 , rS,0 , rH,0 en nN,0 die de exponent is in het Sérsicprofiel van de nucleus. Zoals in figuur 3.12 al te zien was, zal de convolutie tussen het werkelijke profiel en het waargenomen profiel enkel een effect hebben op de kleinere stralen. We hoeven dus enkel de convolutie van de PSF met de nucleus te beschouwen. Je zou een convolutie van allemaal kunnen nemen ook, maar het computationeel werk voor de fit procedure mag niet onderschat worden en kan een hoop tijd in beslag nemen. Daarom proberen we de te fitten functie zo eenvoudig mogelijk te houden (we fitten al 7 parameters en de convolutie zorgt ervoor dat voor ieder datapunt een dubbelintegraal uigerekend moet worden!). De convolutie van de PSF met de nucleus ziet er dus uit zoals vergelijking 3.17, ∞

p r02 + r2 + 2rr0 cos(θ) 1/nN,0 r dr − rN,0 0 0 (3.22) Het is niet nodig om van 0 tot +∞ te integreren, maar een aantal a keer sigma volstaat, IN,0 fN (r) = 2πσ 2

Z

0

0

Z

2π

r02 dθ exp − 2 exp 2σ


(a)

47

(b)

(c)

Figuur 3.13: (a) Het oppervlaktehelderheidsprofiel bekomen uit de data van beide telescopen en (b) dezelfde data waarbij de datapunten van de HST verschoven werden naar de datapunten van de R-band van de VLT. (c) De datapunten van beide datasets gecombineerd tot 1 profiel dat zowel een goede beschrijving geeft voor kleine stralen als voor grote stralen.

omdat de PSF naar nul nadert voor grote waarden van r0 .

p r02 + r2 + 2rr0 cos(θ) 1/nN,0 r dr − rN,0 0 0 (3.23) Omdat de PSF een Gaussische verdeling is zou je op het eerste zicht verwachten dat 3 à 4 keer sigma genoeg is. Dit is echter niet zo omdat het Sérsicprofiel een exponentiële afhankelijkheid kent en voor stralen r − aσr0 nog zwaar zal doorwegen. Hiermee moeten we dus rekening houden en a moet van de grootte 15-25 zijn. Het waargenomen profiel IN,0 fN (r) = 2πσ 2

Z

+aσ

0

0

Z

2π

r02 dθ exp − 2 exp 2σ

48


ziet er dan als volgt uit, ftot (r) = fN (r) + IS (r) + IH (r) p Z 2π Z +aσ r02 + r2 + 2rr0 cos(θ) 1/nN,0 IN,0 r02 0 0 dθ exp − 2 exp − r dr = 2πσ 2 0 2σ rN,0 0 r ) + IS,0 exp(− rS,0 r + IH,0 exp(− ). (3.24) rH,0 We kunnen dit nu nog schrijven als, ftot (r) = IN,0 × r02 dθ exp − 2 exp r dr 2σ 0 0 ! r r + kS,0 exp(− ) + kH,0 exp(− ) rS,0 rH,0 1 2πσ 2

Z

+aσ

0

0

Z

2π

p r02 + r2 + 2rr0 cos(θ) 1/nN,0 − rN,0

= IN,0 × Ftot (r),

(3.25)

waarin, IS,0 = kS,0 × IN,0 ,

(3.26)

IH,0 = kH,0 × IN,0 .

(3.27)

De datapunten uit figuur 3.13(c) staan echter in oppervlaktehelderheden. We hebben dus, µ(r) = m0 − 2.5 log(ftot (r)) = m0 − 2.5 log(IN,0 ) − 2.5 log(Ftot (r)) = µ0 − 2.5 log(Ftot (r)),

(3.28)

waarbij we in de laatste stap m0 en IN,0 tot 1 parameter µ0 gecombineerd hebben. Later zullen we IN,0 uit µ0 bepalen. De te fitten parameters zijn dus, µ0 , kS,0 , kH,0 , rN,0 , rS,0 , rH,0 en nN,0 , wat nog steeds 7 parameters zijn. Dit is de vergelijking die we aan de datapunten in figuur 3.13(c) moeten fitten. Om het computationeel werk nog iets draaglijker te maken voor de gebruikte PC, fitten we eerst het profiel IS (r) + IH (r) van de schijf en de halo aan de datapunten waarvan r groter is dan 2 boogseconden en fN (r) aan de datapunten met r kleiner dan 2 boogseconden. Hierna doen we de fitting procedure opnieuw met de 3 componenten samen en maken we gebruik van alle datapunten, maar de schaallengtes van de schijf en de halo worden vast gehouden, zodat er nog slechts 5 te fitten parameters overblijven. De reden dat de schaallengte van de nucleus niet vast gehouden wordt is omdat de datapunten voor de kleinere stralen moeilijker aan het profiel te fitten zijn, dus laten we daar iets meer vrijheid om toch een goede fit te bekomen. De datapunten op de langere afstanden zijn al zeer goed gefit door het profiel van de schijf en de halo en hun schaallengtes kunnen dus vast gehouden worden. Ook zal het profiel van de nucleus helemaal geen effect hebben op stralen groter dan 5 boogseconden, terwijl vooral het profiel van de schijf wel nog een bijdrage zal leveren rond 1 à 2 boogseconden,


49

Figuur 3.14: De data wordt voorgesteld door de rode punten, de groene volle lijn stelt de gefitte convolutie voor van de Sérsicprofielen en de groene stippellijn stelt het werkelijke oppervlaktehelderheidsprofiel voor. Het gefitte model is het profiel waarbij de nucleus door 1 Sérsicprofiel beschreven wordt. dus de schaallengte van de nucleus zal hierdoor aangepast worden. Op figuur 3.14 is het resultaat van deze fit te zien. Deze fit werd gemaakt met de methodes optimize.leastsq en integrate.quad van het scipy pakket van python. De volle lijn stelt de fit van ftot (r), dat dus de convolutie tussen de PSF en het werkelijke profiel is, aan de datapunten voor en de stippelijn stelt het werkelijke profiel Itot (r) = IN (r) + IS (r) + IH (r) voor. We zien dat de binnenkant van de nucleus wel goed gifit wordt door het profiel, maar dat de buitenkant (van de nucleus) het laat afweten. Een oplossing voor dit profiel is om het intensiteitsprofiel van de nucleus anders (en met meer parameters) te gaan voorstellen. We kunnen bijvoorbeeld de som van 2 Sérsicprofielen nemen, r 1/nN,0 r 1/nN,1 IN (r) = IN,0 exp − ( + IN,1 exp − ( , (3.29) ) ) rN,0 rN,1 terwijl we de modellen voor de profielen van de schijf en de halo kunnen behouden, r IS (r) = IS,0 exp(− ), (3.30) rS,0 IH (r) = IH,0 exp(−

r rH,0

).

(3.31)

We kunnen nu het waargenomen profiel schrijven als, ftot (r) = fN (r) + IS (r) + IH (r) p Z +aσ Z 2π 02 r02 + r2 + 2rr0 cos(θ) 1/nN,0 1 r 0 0 = r dr dθ exp − I exp − N,0 2πσ 2 0 2σ 2 rN,0 0 p 02 2 0 r + r + 2rr cos(θ) 1/nN,1 + IN,1 exp − rN,1 r r + I0,S exp(− ) + I0,H exp(− ), (3.32) r0,S r0,H

50


aangepast wat kan geschreven worden als, ftot (r) = IN,1 × p Z 2π Z +aσ r02 + r2 + 2rr0 cos(θ) 1/nN,0 r02 1 0 0 dθ exp − k exp − r dr N,0 2πσ 2 0 2σ 2 rN,0 0 p r02 + r2 + 2rr0 cos(θ) 1/nN,1 + exp − rN,1 ! r r + k0,S exp(− ) + k0,H exp(− ) r0,S r0,H = IN,1 × Ftot (r).

(3.33)

Opnieuw stellen we, IN,0 = kN,0 × IN,1 ,

(3.34)

IS,0 = kS,0 × IN,1 ,

(3.35)

IH,0 = kH,0 × IN,1 ,

(3.36)

µ(r) = µ0 − 2.5 log(Ftot (r)),

(3.37)

en

waarin µ0 = m0 − 2.5 log(IN,1 ). De te fitten parameters zijn nu, µ0 , kN,0 , kS,0 , kH,0 , rN,0 , rN,1 , rS,0 , rH,0 , nN,0 en nN,1 . Net zoals bij de fit waar het profiel van de nucleus uit slechts één Sérsicprofiel bestond, gaan we nu eerst de schijf en de halo aan de data met stralen groter dan 2 boogseconden fitten en de 2 Sérsicprofielen aan de data met stralen kleiner dan 2 boogseconden om dan de fit opnieuw te doen voor alle data en met de schaallengtes van de schijf en de halo die vast gehouden worden. In figuur 3.15 is het resultaat te zien van deze fit. Ook deze fit werd gemaakt met de methodes optimize.leastsq en integrate.quad van het scipy pakket van python. In tegenstelling met het eerste profiel dat gebruikt werd, worden nu wel alle datapunten op een goede manier gefit aan het profiel. De volle lijn is wederom het waargenomen profiel, terwijl de stippellijn het werkelijke profiel is. De gefitte parameters worden weergegeven in tabel 3.3. De parameter kN,1 is geen gefitte parameter maar is hier toch vermeld omdat de onderlinge verhouding van de k’s weergegeven in de tabel gelijk zijn aan die van de centrale intensiteiten van de 4 gebruikte profielen, maar om deze centrale intensiteiten te weten te komen moet gebruik gemaakt worden van de verhoudingen onderling en de waarde µ0 Bepaling van de centrale intensiteiten en de schaallengtes De waarden voor de schaallengtes in tabel 3.3 staan nog uitgedrukt in boogseconden. Verder zijn enkel de verhoudingen van de centrale intensiteiten weergegeven. Met behulp van µ0 kan IN,1 berekend worden en zo ook de 3 andere centrale intensiteiten. De omrekening van de schaallengtes kan makkelijk gedaan worden als de afstand tot de galaxy gekend is. Deze is terug te vinden in de Nasa/Ipac Extragalactic Database (NED) waar een waarde


51

Figuur 3.15: De data wordt voorgesteld aan de rode punten, de groene volle lijn stelt de gefitte convolutie voor van de Sérsicprofielen en de groene stippellijn stelt het werkelijke oppervlaktehelderheidsprofiel voor. Het gefitte model is het model waarbij de nucleus uit de som van 2 Sérsicprofielen bestaat. Tabel 3.3: De gefitte parameters voor het profiel voorgesteld door vergelijkingen 3.32 en 3.36. parameter µ0 rN,0 rN,1 rS,0 rH,0 kN,0 kN,1 kS,0 kH,0 nN,0 nN,1

gefitte waarde 17.71 0.00538500 0.242100 6.97800 19.3700 141.4 1 0.1144 0.01125 1.964 1.126

van D = 16.375 ± 2.157 Mpc gegeven wordt. Nu kunnen we de schaallengtes omzetten van boogseconden naar kpc met behulp van, R = D tan(θ),

(3.38)

met R de grootte van de schaallengte uitgedrukt in pc en θ de grootte van de schaallengte uitgedrukt in boogseconden. Om de centrale intensiteiten te vinden, moeten we eerst de oppervlaktehelderheid µR uitrekenen die overeenkomt met 1 L ,R pc−2 , wetende dat de absolute R-band magnitude van de zon M ,R = 4.42 is. We starten van een object dat als magnitude m heeft en een ruimtehoek van p2 arcsec2 omspant, lichtkracht LX heeft en

52


een afstand dx van ons verwijderd is. De uitdrukking voor de magnitude per arcsec2 µX is nu, LX ! µX = M ,R − 2.5 log

4πp2 d2x L ,R 4πd20

.

(3.39)

We kunnen nu LX = D2 IX stellen met D de diameter van het object en IX uitgedrukt in L ,R /pc2 . Gebruik makend van het feit dat d0 per definitie gelijk is aan 10 pc, kan bovenstaande vergelijking herschreven worden als, µX = M ,R − 2.5 log(

IX 10D ) − 5 log( ). 2 (L ,R /pc ) p dX

(3.40)

Nu kunnen we α = D/dX stellen waarbij α de grootte is van het object uitgedrukt in radialen. Verder willen we de oppervlaktehelderheid weten die overeenkomt met IX = 1L ,R /pc2 . We krijgen dan, α µR = M ,R − 5 − 5 log( ) p 360 ) π × 60 × 60 = 25.992 mag/arcsec2 .

= M ,R − 5 − 5 log(

(3.41)

Nu hebben we voor onze galaxie een centrale oppervlaktehelderheid van, µ(0) = µ0 − 2.5 log(kN,0 + kN,1 + kS,0 + kH,0 ) = 12.32 mag/arcsec2

(3.42)

De centrale oppervlaktehelderheid kan nu ook uitgedrukt worden in functie van de centrale intensiteit, IN,1 (kN,0 + kN,1 + kS,0 + kH,0 ) µ(0) − µR = −2.5 log( ), (3.43) L ,R /pc2 waaruit IN,1 kan berekend worden. Hieruit volgt een waarde IN,1 = 2.0545e9 L ,R /kpc2 . De andere waarden zijn weergegeven in tabel 3.4 en in figuur 3.16 zijn de intensiteitsprofielen van de nucleus, de schijf en de halo geplot. Deze figuur werd gemaakt met behulp van python.

3.3.6

De totale intensiteit en de schijnbare magnitude

De totale lichtkracht in de R-band kunnen we berekenen door de Sérsicprofielen van de vorm 3.8 te integreren, Z 2π Z ∞ LR = dφ rdrI(r) 0 0 Z ∞ r 1/n = 2πI0 exp − ( ) rdr. (3.44) r0 0 We kunnen nu de substitutie u = (r/r0 )1/n doorvoeren en we krijgen dan, Z ∞ 2 LR = 2πI0 nr0 exp(−u)u2n−1 du. 0

(3.45)


53

Tabel 3.4: De gefitte parameters voor het profiel voorgesteld door vergelijkingen 3.32 en 3.36. parameter rN,0 rN,1 rS,0 rH,0 IN,0 IN,1 IS,0 IH,0 nN,0 nN,1

gefitte waarde 4.275 × 10−4 kpc 0.01922 kpc 0.5540 kpc 1.538 kpc 2.905 × 1011 L ,R /kpc2 2.055 × 109 L ,R /kpc2 2.35 × 108 L ,R /kpc2 2.311 × 107 0.01125 L ,R /kpc2 1.964 1.126

Figuur 3.16: De 3 bijdragen aan de totale intensiteit. De rode curve stelt de bijdrage van de nucleus voor, de groene curve de bijdrage van de schijf en de blauwe curve de bijdrage van de halo. Voor de nucleus is een som van de 2 Sérsicprofielen gebruikt en voor zowel de schijf als de halo een exponentieel model.

De integraal kan geschreven worden als een gammafunctie Γ(2n) en dit levert een resultaat voor de totale lichtkracht op van,

LR = 2πI0 nr02 Γ(2n).

(3.46)

54


De totale lichtkracht van de galaxie kan nu uitgedrukt worden in termen van de gefitte parameters van de verschillende componenten, wat als resultaat, 2 2 LR,galaxie = 2π IN,0 nN,0 rN,0 Γ(2nN,0 ) + IN,1 nN,1 rN,1 Γ(2nN,1 ) 2 2 + IS,0 rS,0 Γ(2) + IH,0 rH,0 Γ(2) = 8.14 × 108 L

(3.47)

oplevert en waarbij rekening gehouden werd met het feit dat voor de schijf en de halo een exponentieel profiel gebruikt werd, waarvoor logischerwijs nS,0 = nH,0 = 1 geldt. Uit de lichtkracht kan nu de schijnbare magnitude van de galaxie bepaald worden met behulp van, ! mR,galaxie − MR, = −2.5 log

LR,galaxie D2 LR, d20

,

(3.48)

waarin D de afstand tot de galaxie is,d0 = 10pc en MR, de absolute magnitude van de zon in de R-band is. Zoals eerder vermeld wordt MR, = 4.42 en D = 16.375M pc gebruikt. Met deze waarden vinden we een schijnbare magnitude mR,galaxie = 13.21 ± 0.66 terug in de R-band. De fout hierop werd bepaald door rekening te houden met de fout op de afstand tot de galaxie, die de fout op de schijnbare magnitude zal domineren (tegenover de fouten op de andere gebruikte waarden in de formule). Deze waarde is consistent met de gevonden waarde in sectie 3.2 waar mR,galaxie = 13.15 gevonden werd en met de waarde 13.17 die door Jerjen et al (2000) gegeven wordt.

3.4 3.4.1

Het dichtheidsprofiel De massa-lichtkrachtverhouding

Om uit de gedeprojecteerde intensiteit de dichtheid te bepalen moeten we gebruik maken van de massa-lichtkrachtverhouding. Deze werd bepaald door dr. M. Koleva door een fit te maken aan het spectrum van IC3328 (bekomen uit de SDSS database) uitgaande van modellen voor de sterpopulatie. Deze fit werd gemaakt met ULySS, een software pakket geschreven door Koleva et al (2009) en is terug te vinden op http://ulyss.univ-lyon1.fr/ . Dit pakket fit het volledig spectrum en minimaliseert het verschil tussen het geobserveerde spectrum en een grid van modellen. In ons geval werden SSP’s (Single stellar populations) gebruikt. De fit van het spectrum is weergegeven in figuur 3.17. De beste fit is getoond in blauw en het spectrum is weergegeven in het zwart. De gemaskte pixels zijn geplot in het rood en de lijn in cyaan is een polynoom gebruikt in de fit om te corrigeren voor errors in de fluxkalibratie of onzekerheden in de extinctie. Omdat we niet zeker weten hoe de stellaire evolutie er uit ziet in de galaxie werden er 2 modellen gebruikt. Pegase-HR (Le Borgne et al., 2004) en Vazdekis (Vazdekis et al., 2010) werden gebruikt als modellen. hiermee werden de waarden voor de leeftijd en de metalliciteit bepaald zoals weergeven in tabel 3.5. Uit de best gefitte SPP parameters kunnen we dan de massa-lichtkracht verhouding bepalen. De stellaire populatie geeft ons de massa van de sterren weer voor een gegeven leeftijd en metalliciteit, en als we de lichtkracht weten, kunnen we de massalichtkrachtverhoudingen bepalen, welke ook weergegeven is in tabel 3.5. Het verschil in

3.4. HET DICHTHEIDSPROFIEL

55

Figuur 3.17: De fit aan het spectrum van de galaxie. De blauwe lijn stelt de fit voor en de zwarte lijn stelt het waargenomen spectrum voor.

Tabel 3.5: De leeftijd, metalliciteit en de massa-lichtkrachtverhouding in de R-band voor beide modelen. leeftijd [F e/H] YR

Pegase-HR 4 Gyr −0.46 2.4

Vazdekis 3 Gyr −0.41 1.2

deze waarden kan verklaard worden door het gebruik van hoge- en lage-massasterren in de stellaire populatie. De werkelijke waarde voor de massa-lichtkrachtverhouding zal dus waarschijnlijk ergens tussen deze 2 uitersten liggen. Om in hoofdstuk 4 een distributiefunctie te fitten zullen we de waarde tussen 1.2 en 2.4 laten variëren. Een beste fit voor de distributiefunctie werd daar bekomen met een waarde YR = 1.6. Voor de figuren die in de komende secties gemaakt worden zullen we dus ook altijd deze waarde gebruiken voor de massa-lichtkrachtverhouding.

3.4.2

Dichtheidsprofiel van de sferische componenten

Deprojectie van het intensiteitsprofiel Als we een galaxie waarnemen is het logisch dat het vermogen van het totaal aantal fotonen dat vanuit een bepaalde richting op de detector valt de som is van het vermogen dat de sterren uitzenden langs de gezichtslijn. Het vermogen dat door de sterren uitgestraald wordt per kubieke meter is de emissiviteit, wat een intrinsieke grootheid is. Wat we zien is het vermogen per vierkante meter, genaamd intensiteit, wat de projectie is van de emissiviteit langs de gezichtslijn. Als we een coördinatenstelsel definiëren, de z-as langs de gezichtslijn leggen en de x-as en y-as in het vlak van de hemel en de oorsprong van het assenstelsel laten samenvallen met de oorsprong van de galaxie, kunnen we de projectie

56


als volgt schrijven, Z

∞

IR (x, y) =

jR (x, y, z)dz.

(3.49)

−∞

De kromming aan de hemel kunnen we verwaarlozen omdat de afstand tot de galaxie veel groter is dan de diameter ervan. Het subscript R duidt er op dat we in de R-band werken en dus de intensiteit en emissiviteit ook in die band aan elkaar relateren. Hetzelfde kan gedaan worden voor andere banden of voor specifieke golflengten. Voor de nucleus en de halo van de galaxie kunnen we een sferische symmetrie veronderstellen, wat ervoor zal zorgen dat de intensiteit en de emissiviteit enkel afhangt van de straal. Vanuit deze veronderstelling kunnen we bovenstaande formule herschrijven als, Z ∞ IR (r) = jR (r)dz, (3.50) −∞

met, r=

p p R 2 + z 2 = x2 + y 2 + z 2 .

(3.51)

waarin R de geprojecteerde straal voorstelt. Aangezien de emissiviteit enkel van de straal afhangt, kunnen we de integraal herschrijven als, Z ∞ IR (r) = 2 jR (r)dz. (3.52) 0

We kunnen nu de integratieveranderlijke veranderen in r, en met behulp van, r=

p p rdr R2 + z 2 ⇒ z = r2 − R2 ⇒ dz = √ , 2 r − R2

krijgen we Z

(3.53)

∞

j (r)rdr √R . (3.54) r 2 − R2 R Deze vergelijking maakt het mogelijk om de intensiteit te berekenen uit de emissiviteit. Praktisch gezien is het altijd de intensiteit van een galaxie die gemeten wordt en hebben we dus nood aan het inverse proces, een deprojectie van de intensiteit om zo het uitgestraalde vermogen per kubieke meter te bepalen. In het geval van een sferisch symmetrische emissiviteit is het mogelijk om de deprojectie exact uit te voeren. Vergelijking 3.54 is een Abelse integraalvergelijking, waarvoor de oplossing gegeven wordt door, Z 1 ∞ dIR dR jR (r) = − . (3.55) (R) √ π r dR R2 − r 2 IR (r) = 2

Uitdrukking voor de dichtheid We kunnen de emissiviteit die uitgedrukt is in vermogen per kubieke meter linken aan de massadichtheid, uitgedrukt in kilogram per kubieke meter met behulp van de massalichtkrachtverhouding YR = M/LR in de R-band, die gelijk gesteld werd aan 1.4 M /L ,R op de volgende manier voorstellen: jR (r) =

1 LR ρ(r) = ρ(r). YR M

(3.56)

3.4. HET DICHTHEIDSPROFIEL

57

Nu hebben we voor de sferische componenten een intensiteitsprofiel dat verloopt volgens, I(R) = I0 exp(−(R/r0 )n ).

(3.57)

Hiervan hebben we als afgeleide, dI(R) Rn−1 exp(−(R/r0 )n ). = I0 n − n dR r0

(3.58)

Dit invullen in vergelijking 3.55, wat we op zijn beurt in vergelijking 3.56 invullen, levert de volgende uitdrukking op voor de dichtheid, ρ(r) =

I0 YR n π

Z r

∞

Rn−1 dR exp(−(R/r0 )n ) √ . r0n R2 − r 2

(3.59)

Tenslotte kunnen we nog overgaan op een dimensieloze variabele, wat het computationeel werk een stuk eenvoudiger en overzichtelijker zal maken. Hiervoor voeren we x = r/r0 in en noteren we X = R/r0 , wat nu de integratieveranderlijke is in onderstaande vergelijking, ρ(x) = =

Z I0 YR n ∞ n−1 dX X exp(−X n ) √ πr0 x X 2 − x2 I0 YR n J(x). πr0

(3.60)

In figuur 3.18 zijn de dichtheidsprofielen van de nucleus en de halo te zien, die uitgerekend werden met behulp van het python pakket scipy en de methode integrate.quad. Om de figuur te bekomen werd vergelijking 3.60 ge¨ıntegreerd.

Figuur 3.18: Het dichtheidsprofiel van de nucleus (rood) en de halo (blauw).

58

3.4.3


Oppervlaktedichtheidsprofiel van de schijf

De oppervlaktedichtheid van de schijf kan zonder veel problemen rechtstreeks uit de intensiteit afgeleid worden. Door deze met de massa-lichtkrachtverhouding te vermenigvuldigen krijgen we direct de oppervlaktedichtheid, Σ(r) = YR IR (r) = YR I0 exp(−

r ). r0

(3.61)

Ook hier kunnen we op de dimensieloze variabele x = r0 overgaan, Σ(x) = YR I0 exp(−x)

(3.62)

In figuyr 3.19 is een plot van de oppervlaktedichtheid weergegeven. Deze figuur werd gemaakt met behulp van python.

Figuur 3.19: De oppervlaktedichtheid van de schijf.

3.4.4

De totale massa van de galaxie

De bijdrage tot de massa voor de sferische componenten kan berekend worden door, Z M

2π

=

Z

π

dφ 0

Z sin(θ)dθ

0

+∞

ρ(r)r2 dr

0

+∞

Z

ρ(r)r2 dr

= 4π 0

=

4πr03

Z 0

+∞

ρ(x)x2 dx

(3.63)

3.5. DE POTENTIAAL

59

en voor de component van de schijf is deze gegeven door, Z +∞ Z 2π Σ(r)rdr dφ M = 0 0 Z +∞ Σ(r)rdr = 2π 0 Z +∞ 2 = 2πr0 Σ(x)xdx

(3.64)

0

We kunnen nu numeriek deze integralen oplossen met voor ρ(x) en Σ(x) de uitdrukkingen gegeven door 3.60 en 3.62. De waarden voor de massa’s van de verschillende componenten en die van de totale galaxie zijn weergegeven in tabel 3.6. We kunnen met behulp van de massa-lichtkrachtverhouding ook de massa berekenen uit de lichtkracht berekend uit de intensiteit, die gegeven wordt door vergelijking 3.47. Met deze totale lichtkracht bekomen we een massa 1.30 × 109 M , wat hetzelfde is als wat bekomen werd uit de integraal over de dichtheid en die weergegeven is in de tabel. Hieruit kunnen we dus besluiten dat er in de numerieke uitrekening van de dichtheid geen computationele fouten geslopen zijn. Tabel 3.6: De massa van de verschillende componenten en de totale massa van de galaxie. Component Nucleus Schijf Halo Totaal

3.5 3.5.1

Massa (M ) 1.56 × 107 7.32 × 108 5.55 × 108 1.30 × 109

massafractie 0.012 0.562 0.426 1

De potentiaal De bijdrage van de materie

Om de bijdragen van de materie aan de potentiaal te berekenen, vertrekken we van de Poisson vergelijking, gegeven door, ∆Φ(r) = 4πGρ(r),

(3.65)

waarin φ(r) de potentiaal is, die wegens de sferische massaverdeling enkel van de straal r zal afhangen. We zullen nu bindingspotentiaal Ψ(r) definieren als volgt, Ψ(r) = −Φ(r),

(3.66)

zodat Ψ(r) overal positief zal zijn. We kunnen de Poisson vergelijking nu herschrijven als, ∆Ψ(r) = −4πGρ(r).

(3.67)

De potentiaal van de sferische componenten We kunnen vergelijking 3.67 nu schrijven met de Laplaciaan uitgedrukt in bolcoördinaten. Dit leidt voor de sferische potentiaal tot, 1 d 2 dΨ r = −4πGρ(r). (3.68) r2 dr dr

60


Deze vergelijking een eerste maal integreren levert, Z r 2 dΨ(r) 2 dΨ(r) ρ(r0 )r02 dr0 = −4πG r − lim r r→0 dr dr 0

(3.69)

Verder is het zeer aannemelijk om een kracht te veronderstellen die verdwijnt in het centrum, zodat de limiet van de tweede term gelijk is aan 0. We kunnen deze vergelijking dan verder integreren en bekomen de volgende uitdrukking voor de potentiaal, Z +∞ Z t dt Ψ(∞) − Ψ(r) = Ψ(r) = 4πG ρ(r0 )r02 dr0 , (3.70) 2 t r 0 waarin gebruik gemaakt werd van het feit dat de potentiaal verdwijnt voor grote stralen. Deze dubbelintegraal kan herleid worden tot een enkelvoudige integraal door een partiële integratie uit te voeren, +∞ Z Z +∞ 4πG u −1 0 02 0 Ψ(r) = − ρ(r )r dr − 4πG ρ(u)u2 du u 0 u r r Z Z +∞ 4πG r 0 02 0 = ρ(r )r dr + 4πG ρ(r0 )r0 dr0 . (3.71) r 0 r Ook hier kunnen we terug overgaan op de dimensieloze variabele x = r/r0. We krijgen dan de volgende uitdrukking voor de potentiaal, Z Z +∞ 4πGr02 x 0 02 0 2 ρ(x )x dx + 4πGr0 ρ(x0 )x0 dx0 . (3.72) Ψ(r) = x 0 x Hierin kunnen we dan de uitdrukking voor ρ(r) (3.60) substitueren en dan krijgen we, Ψ(r) = 4Gr0 I0 YR n

1 Z

waarin,

x

r

0

02

0

J(x )x dx +

0

∞

Z

X n−1 exp(−X n ) √

J(x) =

Z

x

+∞

J(x0 )x0 dx0 ,

(3.73)

x

dX . X 2 − x2

(3.74)

De potentiaal van de schijf In tegenstelling tot de potentiaal van de sferische componenten, heeft de potentiaal van een dunne schijf, met als oppervlaktedichtheid een profiel van de vorm, Σ = Σ0 exp(−(r/r0 )),

(3.75)

wel een analytische oplossing, zoals gegeven wordt door Binney en Tremaine (1988), Ψ(r) = −πGΣ0 r I0 (y)K1 (y) − I1 (y)K0 (y) , (3.76) met y = r/(2r0 ) en In en Kn gemodifieerde Besselfuncties van respectievelijk de eerste en de tweede soort.

3.5. DE POTENTIAAL

61

Figuur 3.20: De bijdrage tot de potentiaal van de 3 componenten van de galaxie, de nucleus, de schijf en de halo. Nu we voor alle componenten een uitdrukking hebben voor de potentiaal, kunnen we deze ook berekenen en deze is weergegeven in figuur 3.20. Deze figuur werd bekomen door in python vergelijking 3.73 te integreren voor de sferische componenten met behulp van de methode integrate.quad van het scipy pakket. Voor de potentiaal van de schijf werd vergelijking 3.76 gebruikt. Het testen van de potentiaal Omdat we de potentiaal volledig numeriek hebben uitgerekend en het numeriek rekenwerk van een computer ook niet feilloos is, kunnen we nu de dichtheid bepalen uit de potentiaal met behulp van de Poisson vergelijking (3.67). Uit vergelijking 3.68 volgt dan, 2 dΨ d2 Ψ + 2 = −4πGρ(r), r dr dr

(3.77)

waaruit we ρ(r) kunnen halen, ρ(r) =

2 dΨ d2 Ψ + 2 /(−4πG) r dr dr

(3.78)

De afgeleiden van de potentiaal kunnen berekend worden door in python een spline te trekken door de potentiaal met behulp van de methode interpolate.UnivariateSpline van het scipy pakket. In figuur 3.21 is de dichtheid die berekend werd op 2 verschillende manieren weergegeven. De volle lijn toont de dichtheid zoals die berekend werd uit de potentiaal met behulp van de vergelijking van Poisson en de stippellijn heeft de potentiaal voor zoals berekend uit vergelijking 3.60.

62


Figuur 3.21: De dichtheid berekend uit de potentiaal via vergelijking 3.78 (volle lijnen) en de dichtheid berekend uit de gedeprojecteerde intensiteiten via vergelijking 3.60 (x punten). De zwarte lijn en punten stelt de totale dichtheid van de nucleus en de halo voor. We zien dat beide leden van de vergelijking hiervoor overeen komen. De afwijkingen voor de nucleus zijn te wijten aan computationele fouten, maar zijn niet van belang omdat de bijdrage van de nucleus al 6 grootte ordes kleiner is dan de totale potentiaal op de plaats waar de fouten optreden.

3.5.2

De bijdrage van de donkere materie

Om de donkere materie halo te beschrijven waarin de dwerggalaxie ingebed is, nemen we een profiel zoals voorgesteld door Burkert (1995), ρDM (r) =

ρ0 r03 . (r + r0 )(r2 + r02 )

(3.79)

Hierin zijn r0 en ρ0 met elkaar verbonden door de relatie, ρ0 = 4.5 × 10−2 (r0 kpc−1 )−2/3 M pc−3 .

(3.80)

Omdat ρ0 geen rechtstreeks waarneembare grootheid is, werd door Burket de relatie tussen v0 en r0 opgesteld, met v0 de bijdrage tot de cirkelsnelheid door ρ0 op r0 , v0 = 17.7(r0 kpc−1 )2/3 km/s.

(3.81)

Dus als we een waarde voor v0 meegeven, kan de vergelijking voor de dichtheid opgelost worden. We kunnen nu de uitrdukking van de donkere materiedichtheid in vergelijking

3.5. DE POTENTIAAL

63

3.72 stoppen, wat de volgende uitdrukking oplevert, Ψ(r) = 4πGr02 ρ0

1 x

Z 0

x

x02 dx0 + (x0 + 1)(x02 + 1)

Z x

+∞

! x0 dx0 , (x0 + 1)(x02 + 1)

die analytisch op te lossen is tot, r0 r0 2 Ψ(r) = πGr0 ρ0 log(r02 + r2 )( − 1) + 2( + 1) log(r0 + r) r r r0 r0 − 2( + 1) arctan(r/r0 ) − 4 log(r0 ) + π . r r

(3.82)

(3.83)

De plot voor de dichtheid en de potentiaal van de donkere materie is weergegeven in figuur 3.22. Deze figuur werd opnieuw gemaat met behulp van python. De parameters voor de donkere materiedichtheid werden gekozen zodanig dat de massa aan donkere materie ongeveer 90% van de totale massa uitmaakt (dus 10% van de totale massa is afkomstg van de sterren in de galaxie).

64


Figuur 3.22: Plot voor de dichtheid (boven) en de potentiaal (onder) van de donkere materie.

3.5. DE POTENTIAAL

65

De cirkelsnelheid De cirkelsnelheid van een ster is gerelateerd aan de middelpuntzoekende kracht die op de ster inwerkt, die het gevolg is van de gravitationele potentiaal, ∂2r ∂Φ ∂Ψ mvc2 =m 2 =m = −m . r ∂t ∂r ∂r

(3.84)

Hierin werd opnieuw gebruik gemaakt van het feit dat Φ = −Ψ de gravitatiepoteniaal is. We kunnen de cirkelsnelheid nu uitdrukken als, r ∂Ψ vc = −r . (3.85) ∂r In figuur3.23 is de cirkelsnelheid in de potentiaal weergegeven, de rode en de zwarte curve stellen respectievelijk de bijdrage voor van de materie en de donkere materie. De groene curve stelt de cirkelsnelheid voor in de totale potentiaal. Deze figuur werd ook opnieuw bekomen door gebruik te maken van python en de module interpolate.UnivariateSpline.

Figuur 3.23: De cirkelsnelheid in de totale potentiaal (groene lijn) en de bijdragen van de materie (rode lijn) en de donkere materie (zwarte lijn).

66


Hoofdstuk 4

De distributiefunctie Vooraleer we met de stabiliteitsanalyse van de galaxie kunnen beginnen moeten we eerst een model hebben voor de schijf van de galaxie. Dit model zal gefit worden aan de eerder beschreven potentiaal en dichtheid van de galaxie en de schijf. Verder zal ook data over de snelheid en de snelheidsdispersies gebruikt worden om een zo goed mogelijk model voor de schijf te reproduceren. Hierna kan er in hoofdstuk 5 overgegaan worden tot de stabiliteitsanalyse, die gebruik zal maken van dit model voor de schijf. In dit hoofdstuk zullen we eerst kijken naar de beschikbare kinematische data die voorhanden is en de gemeten waarden voor de snelheid en dispersies langs de gezichtslijn bekijken. Hierna kunnen we overgaan naar een beschrijving van de distributiefunctie en haar momenten en bekijken hoe die in het progamma, dat een model voor de schijf fit, ge¨ımplementeerd is.

4.1

Kinematische data

Naast de potentiaal en de dichtheid hebben we ook nog informatie nodig over de snelheid en de snelheidsdispersie van de sterren. De gebruikte data voor deze snelheden en dispersies is de data verzameld door Prugniel en Simien (2002) en is weergegeven in tabel 4.1 en figuur 4.1.

Figuur 4.1: Links: De snelheid langs de gezichtslijn. Rechts: De dispersie van de snelheid langs de gezichtslijn. 67

68

HOOFDSTUK 4. DE DISTRIBUTIEFUNCTIE

Tabel 4.1: De data voor de snelheid en de snelheidsdispersie, met de fout erop. straal(kpc) -1.882 -1.318 -0.897 -0.611 -0.4287 -0.2858 -0.1905 -0.0953 -0.0079 0.0873 0.1826 0.27786 0.42076 0.6034 0.88915 1.3099 1.8656

vl.o.s. -27 -25 -9 -16 -6 -8 -1 -7 1 -1 2 0 8 8 8 7 4

σvl.o.s. 16 6 4 3 3 3 3 2 2 2 2 2 3 3 4 6 9

σl.o.s. 24 33 31 39 34 34 40 35 39 31 34 25 38 33 40 38 26

σσl.o.s. 25 8 5 4 4 4 4 3 3 3 3 3 4 4 5 8 15

We kunnen enkel de snelheid en de snelheidsdispersie langs de gezichtslijn (l.o.s. = line of sight) waarnemen. We kunnen de waargenomen snelheid schrijven als, 2 = vφ2 sin2 (i) + vz2 cos2 (i). vl.o.s.

(4.1)

Uit de analyse gedaan in hoofdstuk 3 bleek dat de schijf domineert tot ongeveer 30 à 40 boogseconden van het centrum van de galaxie, wat overeenkomt met ongeveer 2.3 à 3 kpc. Dit wil zeggen dat alle kinematische data die we hebben afkomstig zal zijn van de schijf. Hierdoor kunnen we veronderstellen dat de z-component van de snelheid geen bijdrage zal leveren aan vl.o.s. en hebben we dus een relatie tussen de snelheid langs de gezichtslijn en de gemiddelde rotatiesnelheid vφ , vl.o.s. = vφ sin(i),

(4.2)

waarin i de inclinatie van de galaxie is en i = 0◦ een face-on galaxie en i = 90◦ een edge-on galaxie voorstelt. De inclinatie werd door Lisker en Fuchs (2009) bepaald uit de spiraalstructuur. Zij vonden een inclinatie van i = 30◦ . Analoog als voor de snelheid kunnen we de waargenomen dispersie schrijven als 2 σl.o.s. = σφ2 sin2 (i) + σz2 cos2 (i).

(4.3)

Er bestaat een relatie tussen de dispersie in de z- en de r-richting van de vorm σz = cσr , met c een constante factor. c = 0.8 wordt voorgesteld door Lisker en Fuchs (2009), maar ook voor de waarden c = 0 en c = 0.5 werd een model gefit, maar met 0.8 werd de beste fit voor de distributiefunctie (sectie 4.4.3) bekomen. We kunnen nu bovenstaande vergelijking verder uitwerken tot, 2 = σφ2 sin2 (i) + c2 σr2 cos2 (i). (4.4) σl.o.s.

4.1. KINEMATISCHE DATA

69

Nu kunnen we uit de epicykelbenadering (Bijlage A) de relatie tussen σφ en σr bepalen, σφ2 σr2

=

κ2 , 4Ω2

(4.5)

met, dΩ2 + 4Ω2 . dr Bovenstaande vergelijkingen combineren levert op, κ2 = r

σφ2 σr2

(4.6)

r dΩ2 4Ω2 dr r3 d vc2 = 1+ 2 4vc dr r2 r dvc 1 − = 1+ 2v dr 2 c 1 r dvc = 1+ , 2 vc dr = 1+

(4.7)

waarin vc de cirkelsnelheid in de potentiaal voorstelt. Vergelijking 4.7 invullen in 4.4 levert dan op, 2 = σφ2 sin2 (i) + σl.o.s.

= σφ2

2 c2 σphi

1 r dvc 1+ 2 vc dr

cos2 (i)

! 2 c cos2 (i) . sin2 (i) + 1 r dvc 1+ 2 vc dr

(4.8)

Nu hebben we een relatie gevonden die σφ rechtstreeks aan de waargenomen σl.o.s. linkt, gebruik makende van de inclinatie en de potentiaal, waaruit de cirkelsnelheid en zijn afgeleiden berekend kunnen worden door middel van een spline te trekken door de potentiaal in python. Met behulp van vergelijkingen 4.2, 4.7 en 4.8 kunnen we nu respectievelijk vφ , σr en σφ bepalen uit de gemeten waarden langs de gezichtslijn. Deze zijn allemaal weergegeven in figuur 4.2. In deze figuren staat de data van de snelheden en de dispersies uitgedrukt voor positieve stralen. Hiervoor werden de datapunten met negatieve stralen en snelheden omgezet naar datapunten met positieve stralen en snelheden. Daarna werd er uitgemiddeld over de 2 snelheden die bij eenzelfde straal behoren. Het resultaat is te zien op de figuur.

70


Figuur 4.2: Boven: De gemiddelde rotatiesnelheid van de sterren. Midden: De snelheidsdispersie in de φ-richting. Onder: De snelheidsdispersie in de r-richting.

4.2. DE DISTRIBUTIEFUNCTIE EN HAAR MOMENTEN

4.2 4.2.1

71

De distributiefunctie en haar momenten De distributiefunctie

Als we de dynamische evolutie van een galaxie willen bestuderen, kunnen we dit doen door voor alle sterren in de galaxie de posities ri en de snelheden vi te kennen en door met behulp van de onderlinge interacties tussen de sterren de evolutie van de galaxie te beschrijven. Door het enorme aantal sterren die galaxieën bevatten is dit N-lichamen vraagstuk niet makkelijk op te lossen. N-lichamen simulaties gebruiken dan ook veel minder sterren dan er effectief aanwezig zijn in een galaxie. De exacte posities en snelheden van alle sterren afzonderlijk kennen is echter niet steeds nodig. We kunnen deze ook op een statistische manier beschrijven. Hiertoe definiëren we de distributiefunctie F (r, v, t) in de 6 dimensionale faseruimte zodat, dM = F (r, v, t)drdv,

(4.9)

de massa aan sterren weergeeft op een tijdstip t tussen r en r + dr en een snelheid heeft tussen v en v + dv. De totale massa van de galaxie kan dan weergegeven worden door de integraal te nemen over de faseruimte, Z F (x, t)dx. (4.10) M= Ω(t)

Hierin is x een zes-dimensionale vector die de plaats- en snelheidscoördinaten bevat en Ω(t) de faseruimte waarover ge¨ıntegreerd wordt. Verder kunnen we een infinitisimaal volume element tussen r en r + dr beschouwen waarin deeltjes zitten met een snelheid tussen v en v + dv. Hiermee beschrijven we een zes-dimensionaal volume element in de faseruimte. We zullen er verder vanuit gaan dat de deeltjes onderling niet botsen, wat een plausibele veronderstellling is door de immense afstanden tussen de sterren in een galaxie. Als resultaat van het krachtveld F waarin de sterren zich bevinden zal de snelheid van elke ster in het zojuist beschreven volume element veranderen van v naar v + Fdt en de positie van r naar r + vdt. De massa aan sterren F (r, v, t) is dus gelijk aan de massa in sterren F (r + vdt, v + Fdt, t + dt), of met andere woorden, F (r + vdt, v + Fdt, t + dt) − F (r, v, t) = 0.

(4.11)

Het linkerlid stelt de verandering in F voor gedurende een tijdsinterval dt. We kunnen dit dus nog schrijven als, dF ∂F ∂F ∂F (r, v, t) = +v − ∇Φ = 0, dt ∂t ∂r ∂v

(4.12)

∂r ∂v waarin logischerwijs v = en −∇Φ = gesteld werd, met Φ de potentiaal. Vergelijking ∂t ∂t 4.12 wordt de vergelijking van Boltzmann genoemd.

4.2.2

De momenten

De door ons gebruikte distributiefunctie wordt gebruikt om een schijf te modelleren. Hiervoor kan de schijf voorgesteld worden in twee dimensies. De distributiefunctie zal dan

72


slechts in een vier-dimensionale faseruimte leven, met 2 plaatscoördinaten r, φ en 2 snelheidscoördinaten vr , vφ . De bovenstaande benadering van de distributie blijft echter volledig geldig. De momenten van de distributiefucntie kunnen als volgt worden gedefinieerd, Z µ2i,j (r) =

F (r, v)vr2i vφj dv,

(4.13)

waarin r = (r, φ) en v = (vr , vφ ) en r en φ poolcoördinaten voorstellen. Het is eenvoudig in te zien dat de dichtheid gegeven wordt door i = j = 0 zodat, Z ρ(r) = µ0,0 (r) =

F (r, v)dv.

(4.14)

We weten uit voorgaande hoofdstukken dat we met een sferische dichtheidsverdeling te maken hebben, deze zal dus enkel van r afhangen en niet van φ. We kunnen dus schrijven, Z ρ(r) = µ0,0 (r) =

F (r, v)dv,

(4.15)

waarbij de distributiefunctie nu ook niet meer van de φ-coördinaat afhangt. Dit houdt in dat we de andere momenten enkel moeten kennen met een afhankelijkheid van de rcoördinaat, zoals de data die in figuur 4.2 weergegeven is. Buiten de dichtheid zal het model voor de distributiefunctie nog aan 3 andere momenten gefit worden ook, µ1,0 (r) = ρvr2 (r) = ρσr2 (r) =

Z

F (r, v)vr2 dv

(4.16)

Z µ0,1 (r) = ρvφ (r) = F (r, v)vφ dv Z µ0,2 (r) = ρvφ2 (r) = ρ σφ2 (r) + vφ 2 (r) = F (r, v)vφ2 dv.

(4.17) (4.18)

Hierin is σi2 = vi2 − vi 2 en is logischerwijs vr = 0.

4.3

Integraal van de beweging en theorema van Jeans

Hieronder worden de constanten en integralen van beweging samen met het theorema van Jeans uitgelegd, zoals grotendeels teruggevonden kan worden in de cursus dynamica van galaxieën van prof. dr. H. Dejonghe.

4.3.1

Constante van beweging

Een constante van beweging is een functie van de plaatscoördinaten r, de snelheidscoördinaten v en de tijd t, die constant is langs een oplossing van de bewegingsvergelijkingen. Dus elke functie C(r, v, t) die constant is voor de baan is een constante van beweging. We kunnen ons afvragen hoeveel constanten van beweging er kunnen zijn. Het antwoord is maximaal 5 en dit kan aangetoond worden door onderstaande redenering.

4.4. MODEL VOOR DE SCHIJF

73

Stel dat je 6 onafhankelijke constanten van beweging hebt (onafhankelijk van de tijd). De constanten van beweging zijn dus enkel functie van de plaats- en snelheidscoördinaten. Uit deze 6 vergelijkingen zouden dan de 3 plaatscoördinaten x,y,z en de 3 snelheidscoördinaten vx ,vy ,vz bepaald kunnen worden. Dit houdt in dat er geen beweging meer mogelijk zou zijn. Een baan in de zes-dimensionale faseruimte wordt dus bepaald door 5 onafhankelijke voorwaarden. 6 onafhankelijke voorwaarden bepalen in de faseruimte 1 punt of meerdere ge¨ısoleerde punten en als deze posities tijdsonafhankelijk zijn is dit geen baan.

4.3.2

Integraal van de beweging

Een integraal van beweging is een constante van beweging die de tijd niet expliciet bevat. We kunnen dit nog uitdrukken als, dI (r, v, t) = 0. (4.19) dt Uit voorgaande redenering weten we dat er dus maximaal 5 onafhankelijke integralen van beweging zijn. Verder is een willekeurige functie van integralen van de beweging ook een integraal van de beweging. Als we de baan van een deeltje niet expliciet kennen en we weten niets over eventuele integralen van de beweging, dan kan het deeltje overal in de ruimte aangetroffen worden. Kennen we bijvoorbeeld wel de bindingsenergie van het deeltje, dan ligt er een beperking in het volume waarin het deeltje kan aangetroffen worden. De kennis van de bindingsenergie beperkt dus het toegelaten volume in de ruimte. De bindingsenergie wordt dan een isolerende integraal van beweging genoemd. Deze beperkt het volume waarin een deeltje kan teruggevonden worden. Let wel dat niet alle integralen van de beweging isolerende integralen zijn. Bijgevolg leert zo’n integraal ons niets bij over het deeltje.

4.3.3

Theorema van Jeans

Het theorema van Jeans luidt: • 1) Als de distributiefunctie F louter een functie is van de constanten van de beweging, dan is F oplossing van de vergelijking van Boltzmann (4.12). • 2) Als F voldoet aan de vergelijking van Boltzmann, dan is F een constante van de beweging. • 3) Men kan de oplossingen van de stationaire vergelijking van Boltzmann beperken tot functies van de isolerende integralen.

4.4

Model voor de schijf

Om een model van de schijf te maken wordt het python programma pyQP gebruikt dat ontwikkeld werd door mijn promotor prof. dr. S. De Rijcke naar een idee van Dejonghe (1989). Dit programma gebruikt kwadratische programmering om een zo goed mogelijke fit te maken van de distributiefunctie aan de beschikbare data. Hieronder wordt eerst de distributiefunctie uitgedrukt in termen van de integralen van beweging en daarna wordt gekeken hoe het programma een fit maakt aan de distributiefunctie aan de hand van de beschikbare data.

74

4.4.1


De distributiefunctie in functie van de integralen van beweging

De stelling van Jeans zorgt ervoor dat we de distributiefunctie in functie van de energie E en het draaimoment langs de z-as Jz , wat allebei integralen van de beweging zijn , kunnen schrijven. We hebben dus, F (r, v, t)drdv = F (E, Jz )dV,

(4.20)

waarin V het volume van de vier-dimensionale faseruimte is. Verder zijn de integralen van beweging gegeven door, 1 E = Ψ(r) − (vr2 + vφ2 ), 2 Jz = rvφ .

(4.21) (4.22)

Beiden voldoen aan vergelijking 4.19 en zijn dus integralen van de beweging.

4.4.2

Vorm van de distributiefunctie

We moeten proberen om aan de hand van de beschikbare data, die bestaat uit de dichtheid, de potentiaal en de kinematische data van de snelheden en de dispersies, een distributiefunctie te fitten. Dit is geen gemakkelijke zaak, omdat we a priori geen idee hebben hoe deze distributiefunctie er zal uitzien. Zoals eerder vermeld, kan de distributiefunctie geschreven worden in functie van de gekende integralen van beweging. Verder zullen we de werkelijke distibutiefunctie proberen te benaderen door deze te schrijven als een som van verschillende componenten, X F (E, J) = ci Fi (E, Jz ). (4.23) i

De componenten worden gegeven door, C(E − E0,i )α (Jz − J0,i (E))β als E ≥ E0,i en Jz ≥ J0,i (E) Fi (E, J) = 0 als E < E0,i of Jz < J0,i (E) met, C=

Ψ(0)α , Jmax (rcut )β

(4.24)

(4.25)

wat een normalisatiefactor voor de distributiefunctie is, en J0,i = f × Jmax (E).

(4.26)

Jmax (E) is het maximale draaimoment bij een energie E en f is een factor tussen 0 en 1. Als f = 1 zullen enkel cirkelbanen toegestaan zijn, terwijl voor een kleinere f steeds meer banen met grotere ellipticiteiten zullen voorkomen. Verder is E0 de minimale energie, of met andere woorden de energie van een cirkelbaan, bij een bepaald angulair momentum. E0 is dan de energie waarop 1 f Jmax (E) 2 f (E) = Ψ(rcut ) − − E = 0. (4.27) 2 rcut


75

Met deze distributiefuncties kunnen dan de momenten, Z µ2i,j (r) = F (E, Jz )vr2i vφj dv

(4.28)

berekend worden, welke later gefit kunnen worden aan de datapunten. Bovenstaande vorm van de distributiefunctie is gerelateerd aan 4 parameters α, β, f en rcut . rcut is de cut-off bij een bepaalde straal van een component van de distributiefunctie, of met andere woorden, de maximale straal die bereikt kan worden voor sterren in deze component van de distributiefunctie. f zal beschrijven hoe radiaal, cirkulair de banen van de sterren zijn. Voor grote f (dicht bij 1) zal J0,i dicht bij Jmax (E) liggen en zullen bijgevolg enkel quasi cirkelvormige banen mogelijk zijn, terwijl voor kleine f (dicht bij 0), J0,i klein zal zijn en dus steeds meer radiale banen mogelijk worden. Deze vorm van de distributiefunctie was de initiële implementatie in het programma pyQP, maar houdt echter geen rekening met tegenroterende sterren. We kunnen de implementatie veranderen, zodat deze hier wel rekening mee houdt. Dit kan gebeuren door nog een extra parameter h in te voeren die −1 of 1 kan zijn en die aanduidt of de sterren op een tegenroterende baan (h = −1) of een normaal roterende (h = 1) baan bewegen. We kunnen dit als volgt in de distributiefunctie implementeren, C(E − E0,i )α (|Jz | − J0,i (E))β als E ≥ E0,i , |Jz | ≥ J0,i (E) en hJz ≥ 0 Fi (E, J) = 0 als E < E0,i of |Jz | < J0,i (E) of hJz < 0 (4.29) en in de uitdrukking voor de momenten, j µ2i,j α,β,rcut ,f,h = (h) µ2i,j .

(4.30)

Verder blijven vergelijkingen 4.23, 4.25, 4.26 en 4.27 geldig. We hebben nu enkele parameters die we moeten combineren om de ide component van de distributiefunctie te maken. Deze parameters zijn α, β, f en rcut en h. De gebruikte waarden voor deze parameters zijn gegeven in tabel 4.2. Verder wordt er bij een waarde β = 0 enkel f = 0 genomen om een component van de distributiefunctie te maken. De parameters die weergegeven zijn in de tabel leveren dan 468 combinaties op, wat logischerwijs leidt tot 468 componenten Fi (dus 0 ≤ i ≤ 467). Tabel 4.2: De parameters en hun waarden gebruikt om de componenten van de distributiefunctie te maken. parameter rcut α β f h

4.4.3

gebruikte waarden 1 , 2, 3, 4, 5, 6, 7, 8, 10 2, 5 0, 2, 5 0, 0.05, 0.2, 0.5, 0.8, 0.98 -1, 1

De fit aan de data

Nu we een model hebben voor de distributiefunctie is het de bedoeling dat we dit fitten aan de beschikbare data. Dit gebeurt door de momenten voor de 468 verschillende componenten uit te rekenen voor de stralen r, op welke we over data van de dichtheid ρ, de

76


snelheid vφ , en de snelheidsdispersies σr en σφ beschikken. De momenten die gefit worden zijn, Z Z (4.31) µ0,0 = ρ = dvr dvφ F (E, Jz ), Z Z (4.32) µ1,0 = ρσr2 = dvr dvφ F (E, Jz )vr2 , Z Z µ0,2 = ρ(vφ2 + σφ2 ) = dvr dvφ F (E, Jz )vφ2 , (4.33) Z Z µ0,1 = ρvφ = 2π dvr dvφ F (E, Jz )vφ . (4.34) Verder zal ook de potentiaal van de galaxie en de dichtheid van de schijf gebruikt worden. De potentiaal wordt meegegeven als een reeks datapunten, die kunnen berekend worden zoals in hoofdstuk 3 beschreven. Hierna zal het programma pyQP er een spline door trekken. Verder geven we aan het programma een maximale straal rmax mee waarop de potentiaal en de dichtheid op 0 zal vallen. Dit gebeurt niet op een discrete manier, maar beide functies worden getapet naar 0. Dit wil zeggen dat we de oorspronkelijke functie/spline behouden tot op 95% van rmax en hierna de functie op een gladde manier naar 0 laten naderen op rmax . Noteren we nu de momenten in 1 lijst en noemen we die µj met 1 ≥ j ≥ N en N het aantal datapunten en waarop j duidt op de straal rj . Verder wordt in de lijst ook opgeslagen om welk moment het gaat, (0, 0), (1, 0), (0, 2) en (0, 1), zodat we slechts 1 lijst van datapunten meer hebben, maar deze hoeven dus niet noodzakelijk hetzelfde moment te beschrijven. We kunnen nu de variabele 2 X X 2 χ = wj µ j − ci µi,j (4.35) j

i

definiëren. Hierin zijn µj de datapunten (Het zijn in principe de momenten van de datapunten, maar deze kunnen direct uit de datapunten geconstrueerd worden, want de dichtheid is ook gekend in de datapunten) en µi,j de momenten van de ide component van de distributiefunctie op een straal rj . De gewichtsfactoren worden in het programma gegeven door, 100 wj = 2 . (4.36) rj Verder moeten we ook nog rekening houden met het feit dat de distributiefunctie niet negatief kan worden. Dit gebeurt door de voorwaarde, X F = ci Fi (E, Jz ) ≥ 0 (4.37) i

op te leggen. We beschikken over 468 verschillende componenten om een fit te maken aan de distributiefunctie. We zullen deze echter niet allemaal nodig hebben. Wat gedaan wordt is eerst voor alle componenten een fit maken en de beste hiervan overhouden. Hierna worden alle overige componenten met de eerste component gecombineerd en diegene die samen met


77

de eerste dan de beste fit heeft, wordt ook behouden. Hierna worden deze 2 componenten gecombineerd met de 466 overgebleven componenten om zo opnieuw de beste fit met 3 componenten te vinden. Zo gaat dit verder tot χ2 niet significant meer kleiner wordt. In figuren 4.3(a) en 4.4(a) zijn respectievelijk de fit van de distributiefunctie aan de datapunten en de distributiefunctie zelf weergegeven waarbij geen gebruik gemaakt werd van tegenroterende sterren. De bekomen χ2 uit vergelijking 4.35 bedraagt hier 597, terwijl deze bij de fit voor de distributiefunctie waar wel rekening gehouden werd met tegenroterende sterren 206 bedraagt. In figuur 4.3(b) en 4.4(b) zijn respectievelijk de fit van de distributiefunctie aan de datapunten en de distributiefunctie zelf weergegeven voor het model dat gebruik maakt van tegenroterende sterren. Het is duidelijk te zien dat de tegenroterende sterren nodig zijn om de snelheden in de φ-richting te verlagen en de dispersies σφ te verhogen.

78


(a) Voor een model zonder tegenroterende sterren

(b) Voor een model met tegenroterende sterren

Figuur 4.3: Per figuur hebben we: Linksboven: De datapunten van de snelheidsdispersie in de radiale richting met hun fout en de fit eraan (volle lijn). Rechtsboven: De datapunten van de snelheidsdispersie in de φ-richting met hun fout erop en de fit eraan (volle lijn). Linksonder: de datapunten van de rotatiesnelheden van de sterren en hun fout en de fit eraan (volle lijn). Rechtsonder: de datapunten van de dichtheid en de fit eraan (volle lijn).


79

(a) Voor een model zonder tegenroterende sterren

(b) Voor een model met tegenroterende sterren

Figuur 4.4: De gefitte distributiefunctie aan de datapunten uit figuur 4.3(a) (a) en 4.3(b) (b). Op de x-as staat het pericentrum r− en op de y-as staat het apocentrum r+ , wat respectievelijk de stralen zijn waarop de baan het dichtst en het verst verwijderd is van het centrum van de galaxie. Banen met een negatieve waarde voor het pericentrum duiden op banen met tegenroterende sterren. De kleurschaal geeft de waarde voor de distributiefunctie weer in de eenheden M (kpc)−2 (km/s)−2 .

80


Hoofdstuk 5

Stabiliteitsanalyse In de volgende behandeling van de schijf van de galaxie zullen we veronderstellen dat de spiraalarmen dichtheidsgolven zijn. Dit idee werd voor de eerste keer voorop gesteld door Lin en Shu (1964). Hun zoektocht ging initieel uit naar niet groeiende stationaire patronen. Hieronder volgt een beschrijving over het vinden van de eigenmodes van de sterschijf. Dit zijn patronen die zichzelf versterken, en dus exponentieel zullen groeien. De groei van deze patronen zullen echter terug afgeremd worden, omdat de dynamische component te heet zal worden, waardoor de respons op de geperturbeerde potentiaal terug zal afnemen. Om deze benadering te doen wordt uitgegaan van een galaxie die bestaat uit een donkere materiehalo, een bult en een axiaalsymmetrische schijf. Hier zal echter enkel de theorie beschreven worden van hoe op zoek gegaan wordt naar de eigenmodes van de schijf. De code pyStab (Dury et al., 2008) houdt geen rekening met tegenroterende sterren, dus kan op dit moment nog niet toegepast worden op de door ons verzamelde gegevens. Deze code kan eventueel op een later tijdstip nog uitgebreid worden en dan kan de stabiliteitsanalyse alsnog uitgevoerd worden.

5.1

De eerste-orde Boltzmannvergelijking

We zullen eerst het effect bekijken van een lichte perturbatie op de bindingspotentiaal Ψ(r), die we definieerden als Ψ(r) = −Φ(r) met Φ(r) de gravitatiepotentiaal, op de distributiefunctie. We kunnen de pertubatie op de bindingspotentiaal beschrijven als, Ψ(r, t) = Ψ0 (r) + Ψp (r, t),

(5.1)

waarin Ψ0 (r) de ongestoorde potentiaal is en Ψp (r, t) een lichte perturbatie op de potentiaal voorstelt. De perturbatie van de potentiaal zal een effect hebben op de sterren in de schijf. De distributie zal dus met andere woorden een respons vertonen op de geperturbeerde potentiaal. We kunnen deze dan schrijven als, F (r, v, t) = F0 (E, Jz ) + Fr (r, v, t),

(5.2)

waarin F0 (E, Jz ) de ongestoorde distributiefunctie voorstelt en Fr (r, v, t) de respons van de distributie op de geperturbeerde potentiaal voorstelt. We kunnen nu de Boltzmannver81

82

HOOFDSTUK 5. STABILITEITSANALYSE

gelijking 4.12 schrijven als, ∂ ∂ F0 (E, Jz ) + Fr (r, v, t) + v F0 (E, Jz ) + Fr (r, v, t) ∂t ∂r ∂ + ∇ Ψ0 (r) + Ψp (r, t) F0 (E, Jz ) + Fr (r, v, t) ∂v = 0. (5.3)

dF (r, v, t) = dt

De nulde-orde Boltmannvergelijking is vanzelfsprekend voldaan en om de vergelijkingen verder niet nodig ingewikkeld te maken zullen we de afhankelijkheid van de distributiefunctie en de potentiaal van r, v en t niet meer expliciet schrijven. De eerste-orde Boltzmannvergelijking wordt dan, ∂Ψp ∂F0 ∂Ψ ∂Fr ∂Fr ∂Fr +v + + = 0, ∂t ∂r ∂r ∂v ∂r ∂v

(5.4)

wat nog kan geschreven worden als, ∂Ψp ∂F0 dFr =− . dt ∂r ∂v

5.2

(5.5)

De respons van de distributiefunctie

Nu we de eerste-orde Boltzmannvergelijking bepaald hebben, kunnen we hieruit een uitdrukking voor de respons van de distributiefunctie op de geperturbeerde potentiaal bepalen. Verder kunnen we schrijven, ∂F0 ∂E ∂F0 ∂Jz ∂F0 = + . ∂v ∂v ∂E ∂v ∂Jz

(5.6)

Met behulp van vergelijkingen 4.21 en 4.22 kunnen we schrijven, ∂E ∂v ∂Jz ∂v

= −v

(5.7)

= reφ .

(5.8)

We kunnen dan schrijven, −

∂Ψp ∂F0 ∂r ∂v

= =

∂Ψp ∂F0 ∂Ψp ∂F0 v − ∂r ∂E ∂φ ∂Jz ∂F0 dΨp ∂Ψp ∂F0 ∂Ψp − − . ∂E dt ∂t ∂Jz ∂φ

(5.9)

Nu kan een willekeurige storing van de potentiaal steeds opgebouwd worden als iets van de vorm, Ψp (r, φ, t) = Ψp (r)ei(mφ−ωt) , (5.10)

5.2. DE RESPONS VAN DE DISTRIBUTIEFUNCTIE

83

Wat consistent is met de Poissonvergelijking in cilindercoördinaten, ! ∂ 1 ∂2 ∂2 1 ∂ r + 2 2 + 2 Ψp (r)ei(mφ−ωt) ∆Ψp (r, φ, t) = r ∂r ∂r r ∂φ ∂z ! ∂Ψp (r) 1 ∂ m2 = r − 2 Ψp (r) ei(mφ−ωt) r ∂r ∂r r = −4πGρr (r)ei(mφ−ωt) .

(5.11)

Hierin stelt m het symmetriegetal van de perturbatie voor en komt overeen met het aantal spiraalarmen. Ω = <{ω/m} stelt de rotatiesnelheid voor en ={ω} stelt de groeisnelheid van het patroon voor. Als ω niet-reëel is, zal dit zorgen voor een patroon dat exponentieel kan gaan groeien. Verder stelt ρr (r) de respons van de dichtheid op de geperturbeerde potentiaal voor, die ook uit de respons van de distributiefunctie kan berekend worden, zoals gedaan zal worden in sectie 5.3. Voor een potentiaal van de vorm 5.10 hebben we, ∂Ψp ∂t ∂Ψp ∂φ

= −iωΨp ,

(5.12)

= imΨp .

(5.13)

Nu kunnen we de eerste-orde Boltzmannvergelijking schrijven als, dFr ∂F0 dΨp ∂F0 ∂F0 = +i ω −m Ψp . dt ∂E dt ∂E ∂Jz

(5.14)

We kunnen bovenstaande vergelijking nu integreren naar de tijd (in het linkerlid staat een volledige tijdsafgeleide) langs een ongestoorde baan. We krijgen dan als oplossing de volgende integraal, Fr (r, φ, vr , vφ , t) =

∂F0 Ψp (r)ei(mφ−ωt) ∂E Z t ∂F0 ∂F0 0 0 +i ω −m Ψp (r(t0 ))ei(mφ(t )−ωt ) dt0 , (5.15) ∂E ∂Jz −∞

welke berekend wordt langs een ongestoorde baan die op tijdstip t door het punt (r, φ, vr , vφ ) in de faseruimte gaat. Deze uitdrukking gaat uit van een storing die 0 was op t = −∞ en groeide in de loop van de tijd. Analoog zou ook een integraal tussen t en t = +∞ opgesteld kunnen worden om een storing te beschrijven die afneemt in de tijd. Voeren we nu volgende substitutie door, t0 = t + t00 0

(5.16) 00

φ(t ) = φ(t) + φ(t ),

(5.17)

dan kunnen we de uitdrukking voor de geperturbeerde distributiefunctie nog schrijven als, met r00 = r(t00 ) en φ00 = φ(t00 ), Fr (r, φ, vr , vφ , t) =

∂F0 Ψp (r)ei(mφ−ωt) ∂E Z ∂F0 ∂F0 i(mφ−ωt) 0 00 00 +i ω −m e Ψp (r00 )ei(mφ −ωt ) dt00(5.18) ∂E ∂Jz −∞

84


waaruit volgt dat we de distributiefunctie kunnen schrijven als, Fr (r, φ, vr , vφ , t) = Fr (r, vr , vφ )ei(mφ−ωt) ,

(5.19)

met logischerwijs, Fr (r, vr , vφ ) =

∂F0 Ψp (r) ∂E Z 0 ∂F0 ∂F0 00 00 Ψp (r00 )ei(mφ −ωt ) dt00 . +i ω −m ∂E ∂Jz −∞

(5.20)

Verder volgt uit de epicykelbenadering (A.24) dat we de hoek-coördinaat kunnen schrijven als, φ(t00 ) = Ωt00 + φp (t00 ), (5.21) waarin het lineair deel de rotatie van de ster voorstelt rond het galaxiecentrum met hoeksnelheid Ω en waarbij φp (t) een periodiek deel is van de straal met periode κ, wat het vooren achterlopen op de eenparige rotatie door de radiale beweging gecombineerd met het behoud van draaimoment weergeeft. We kunnen nu eerst de integraal in 5.20 herschrijven als, Z 0 Z 0 00 00 00 i(mφ00 −ωt00 ) 00 Ψp (r )e dt = Φ1 (r)eimφp (t ) ei(mΩ−ω)t dt00 −∞

−∞ 0

Z

00

I(t00 )ei(mΩ−ω)t dt00 ,

=

(5.22)

−∞

met logischerwijs,

00

I(t00 ) = Ψp (r(t00 ))eimφp (t ) .

(5.23)

We kunnen dit nu als een Fourier-reeks ontwikkelen, X 00 I(t00 ) = Il eilκt ,

(5.24)

l

met, Il = =

Z 1 T I(t)e−ilκt dt, T 0 Z 1 T Ψp (r(t))ei(mφp (t)−lκt) dt T 0

Dit invullen in 5.22 levert op, Z 0 Z Ψp (r)ei(mφ−ωt) dt = −∞

0

X

−∞

=

X

00

00

Il eilκt ei(mΩ−ω)t dt00

l

Z

0

Il

00

ei(mΩ−ω+lκ)t dt00

−∞

l

= −i

(5.25)

X l

Il . mΩ + lκ − ω

(5.26)

5.3. DE RESPONS VAN DE DICHTHEID

85

Verder volgt uit 5.23 en 5.24 dat voor t00 = 0 of dus voor de eerste term in 5.20, X Ψp (r(t)) = Ψp (r(t00 = 0)) = I(t00 = 0) = Il .

(5.27)

l

We kunnen de geperturbeerde distributiefunctie dus herschrijven als, +∞ X

F1 (r, vr , vφ ) =

Il

∂F0 0 (lκ + mΩ) ∂F ∂E − m ∂Jz

lκ + mΩ − ω

l=−∞

5.3

.

(5.28)

De respons van de dichtheid

Uit de Poissonvergelijking 5.11 bleek dat de dichtheid kon geschreven worden als, ρr (r, φ, t) = ρr (r)ei(mφ−ωt) .

(5.29)

Deze kan ook nog geschreven worden aan de hand van de respons van de distributiefunctie op de gestoorde potentiaal, Z Z ρr (r) = dvr dvφ Fr (r, vr , vφ ). (5.30) Vergelijking 5.28 hierin invullen en de sommatie laten lopen over een eindig aantal waarden voor l levert op, ρr (r) =

lX max

Z Z Il

∂F0 0 (lκ + mΩ) ∂F ∂E − m ∂Jz

l=−lmax

lκ + mΩ − ω

dvr dvφ .

(5.31)

We kunnen nu overgaan op een ander coördinatensysteem (v, α) dat gerelateerd is aan (vr , vφ ) door, vr = vcos(α),

(5.32)

vφ = vsin(α),

(5.33)

zodat we de respons van de dichtheid kunnen herschrijven als, ρr (r) =

lX max l=−lmax

Z

vesc (r)

Z vdv

0

2π

Il 0

Al (r, v, α) dα. pl (r, v, α) − ω

(5.34)

Hierin stellen Al en pl functies voor van de overgebleven coördinaten (r, v, α). Als we enkel positieve draaimomenten zouden hebben, dus geen tegenroterende sterren, kunnen we α laten lopen van 0 tot π. Voor onze distributiefunctie is echter gebruik gemaakt van tegenroterende sterren, dus zal α van 0 tot 2π moeten lopen. We hebben echter nog steeds met een twee-dimensionale integraal te maken, wat computationeel niet zo eenvoudig te berekenen is. Hieronder zullen we bespreken hoe de integraal als een één-dimensionale integraal van de vorm, Z pmax W (r, p) dp, (5.35) ρr (r) = p−ω pmin

86


kan uitgerekend worden. Dit kan vlugger uitgerekend worden dan de twee-dimensionale integraal. Dit zal handiger zijn als voor verschillende ω’s de geperturbeerde dichtheid moet uitgerekend worden. Het is echter niet zo simpel om ω in 1 van de integralen analytisch naar buiten te brengen. Dit probleem valt echter wel numeriek op te lossen. De uitdrukking in 5.34 heeft een rechthoekig domein [0, vesc ] × [0, 2π]. We kunnen dit nu verdelen in kleinere gebieden met oppervlakte ∆v × ∆α en we verondertsellen dat het integrandum in al deze gebieden constant is. We kunnen dit integrandum verder schrijven als, Al (r, v, α)v∆v∆α . (5.36) l (r, v, α) − ω De twee-dimenionale integraal moet uitgerekend worden door de bijdrage van al deze gebieden, voor elke l uit de sommatie, te sommeren. Om naar de één-dimensionale integraal over te gaan kunnen we bovenstaande uitdrukking gaan schrijven als, ∆p ∆v∆α Al (r, v, α)v , pl (r, v, α) − ω ∆p

(5.37)

die een integraal over de posities van de polen p is. We houden nu bij een bepaalde straal r, W (p) bij als een rij van fucntiewaarden voor een eindig aantal waarden p tussen pmin en pmax in stappen van ∆p. Initieel worden alle waarden W (p) op 0 gezet om daarna voor elk gebied ∆v × ∆α en elke l de p-waarde te zoeken uit de rij van W (p)-waarden die het dichts bij positie P ligt. We tellen bij deze waarde W (p) dan Al (r, v, α)v∆v∆α op. Als we zo alle gebieden afgelopen hebben moet er enkel nog door ∆p gedeeld worden. Op deze manier is de twee-dimensionale integraal 5.34 omgevormd tot een één-dimensionale integraal 5.35.

5.4 5.4.1

Zelf-consistente storingen Eigenmodes van de schijf

Als een storing in een galaxieschijf zelf-consistent is en er geen uitwendige invloeden zijn die op de galaxie inwerken, zal de respons van de massadichtheid aanleiding geven tot een respons van de potentiaal. We hebben eerder de respons van de dichtheid berekend in voorgaande sectie. We kunnen nu de respons van de potentiaal berekenen aan de hand van de Poissonvergelijking, ∇Ψr (r, φ, t) = −4πGρr (r, φ, t),

(5.38)

die de respons zal zijn op de dichtheidsrespons. Deze potentiaal kan ook geschreven worden als, Ψr (r, φ, t) = Ψr (r)ei(mφ−ωt) . (5.39) Als de respons van de potentiaal gelijk is aan de oorspronkelijke perturbatie op de potentiaal, Ψr (r) = Ψp (r), (5.40) voor een bepaalde ω, die de eigentrillingen of eigenmodes van de schijf genoemd worden, dan zal dit aanleiding geven tot het versterken van de originele perturbatie. Hierdoor zal de amplitude van de perturbatie exponentieel gaan groeien. We zullen hieronder een methode beschrijven om deze eigenmodes van de schijf te vinden.

5.4. ZELF-CONSISTENTE STORINGEN

5.4.2

87

Familie van potentialen en dichtheden

De respons van de massadichtheid en de potentiaal zijn verbonden door middel van de Poissonvergelijking, " # ∂Ψr (r) m2 1 1 ∂ r − 2 Ψr (r) . ρr (r) = − 4πG r ∂r ∂r r

(5.41)

We kunnen nu een familie van potentiaal-dichtheidsparen vooropstellen die aan deze vergelijking voldoen. Deze familie zal dan dienen als basis om een willekeurige respons in te ontwikkelen. Beide hebben dezelfde φ en t afhankelijkheid en zijn te schrijven als, Ψm,j (r, φ, t) = Ψm,j (r)ei(mφ−ωt) ,

(5.42)

i(mφ−ωt)

(5.43)

ρm,j (r, φ, t) = ρm,j (r)e

.

De geperturbeerde potentiaal kan nu ook in deze basis geschreven worden als, Ψp (r, φ, t) =

X

aj Ψm,j (r) ei(mφ−ωt) .

(5.44)

j

Als we de potentiaal ontwikkelen in deze basisfuncties, is het van belang dat we met een eindig aantal basisfuncties de potentiaal zo goed mogelijk kunnen beschrijven. Hiervoor is het best om een familie van basisfuncties te hebben die goed lijken op de trillingsmodes van het systeem. Stel nu dat je een infinitesimale verplaatsing hebt in een ongestoorde schijf in een richting in het vlak van de schijf. In dit voorbeeld is de verplaatsing infinitesimaal klein en kunnen we dus lineairiseren. De respons van de dichtheid kan dan gevonden worden door de ongestoorde massadichtheid af te leiden naar bijvoorbeeld de x-coördinaat. We weten dat x = <(reiφ ) en we krijgen dan, ρr (r, φ) =

∂ρ0 dρ0 (r) = (r)eiφ , ∂x dr

(5.45)

wat overeenkomt met een m = 1 storing met ω = 0. Analoog aan dit eenvoudige voorbeeld en rekening houdend met het feit dat de eerste m − 1 afgeleiden naar de straal 0 moeten zijn in het centrum, wat door de Poissonvergelijking opgelegd wordt, kunnen we volgende keuzes voor de massadichtheden voorop stellen, ( ρm,j (r) =

r ρ0 (r) cos(πj rmax ) als m = 0 en k ≥ 0 dρ r m−1 0 −r dr (r) cos(πj rmax ) als m ≥ 0 .

(5.46)

De potentiaal horende bij deze dichtheden kan dan gevonden worden door, Z Ψm,j (r) = G

ρm,j (r) 0 dr . |r − r0 |

(5.47)

88

5.4.3


Zoeken naar eigenmodes van de schijf

Wanneer we een familie van potentialen en massadichtheden hebben, die een basis vormen om de respons van de potentiaal en massadichtheid in te ontwikkelen, kan de respons van de massadichtheid geschreven worden als, ρr,m,l (r, φ, t) =

jX max

cj,l (ω)ρm,j (r)ei(mφ−ωt) .

(5.48)

j=0

Het subscript l wordt gebruikt om de geperturbeerde potentiaal aan te duiden. In bovenstaande vergelijking stelt ρr,l (r, φ) de respons voor op de lde geperturbeerde potentiaal uit de basis. De coëfficienten cj,l zullen namelijk verschillen voor een andere storende potentiaal. Verder is cj,l ook afhankelijk van de gebruikte ω. We kunnen dit ook nog schrijven als, jX max ρr,m,l (r) = cj,l (ω)ρm,j (r), (5.49) j=0

waar de φ en t afhankelijkheden geëlimineerd zijn. Anderzijds kan ρr,m,l,i = ρr,l (ri ) ook berekend worden uit vergelijking 5.35 en met behulp van de kleinste kwadratenmethode die volgende uitdrukking minimaliseert, !2 jX n max X χ2 = ρr,m,l,i − cj,l (ω)ρm,j (ri ) (5.50) i=0

j=0

kunnen de coëfficiënten ck,j (ω) bepaald worden. Nu kunnen we de respons van de dichtheid op de geperturbeerde potentiaal met behulp van 5.44 nog schrijven als, ρr (r) =

lX max max jX

cj,l (ω)aj ρm,j (r).

(5.51)

j=0

l=0

De bijhorende respons van de potentiaal kan dan geschreven worden als, Ψr (r) =

lX max max jX

cj,l (ω)aj Ψm,j (r).

(5.52)

j=0

l=0

We weten nu dat voor de eigenmodes van de schijf de voorwaarde 5.40 voldaan is. De vergelijkingen 5.44 en 5.52 in deze voorwaarde invullen levert op, jX max

cj,l (ω)aj = al ,

(5.53)

j=0

of in matrixvorm, C(ω)A = A.

(5.54)

De eigenmodes ω van de schijf zijn dus diegene waarvoor aan vergelijking 5.54 voldaan is. Dit komt neer op de eigenvectoren zoeken van de matrix C(ω) horend bij de eigenwaarde λ = 1. Deze matrix zal echter niet steeds λ = 1 als eigenwaarde hebben en er moet dus gezocht worden naar de ω waarvoor vergelijking 5.54 geldig is.

Hoofdstuk 6

Besluit We zullen deze masterproef afsluiten met de belangrijkste resultaten hier op een rijtje te zetten. We begonnen in hoofdstuk 2 met de datareductie te doen. Hierna kon in hoofdstuk 3 een model gefit worden aan de galaxie, waaruit het oppervlaktehelderheidsprofiel bepaald kon worden. Hieraan kon dan een profiel gefit worden, waarbij rekening werd gehouden met de bijdragen van de nucleus, de schijf en de halo. De bijdragen van deze 3 componenten tot het oppervlaktehelderheidsprofiel van de galaxie zijn dus gekend. De totale lichtkracht van de galaxie kan hieruit berekend worden en een waarde van 8.14 × 108 L werd hiervoor gevonden. Ook de massa-lichkrachtverhouding in de R-band werd bepaald door een fit te maken aan het spectrum van IC3328. Hiervoor werd een waarde tussen 1.2 en 2.4 gevonden. Uit deze fit volgen ook de leeftijd van de sterren die ongeveer 3 tot 4 Gyr is en de metalliciteit die tussen -0.46 en -0.41 ligt. De fit aan de distributiefunctie in hoofdstuk 4 leverde het beste resultaat met een massa-lichtkrachtverhouding van 1.6 in de R-band. Met deze waarde kan de massa van de dwerggalaxie bepaald worden. Hiervoor werd een waarde 1.30 × 109 M gevonden. Omdat we bij de fit aan het oppervlaktehelderheidsprofiel de 3 componenten afzonderlijk beschouwd hebben, konden we de bijdrage van elk van deze componenten aan de totale massa bepalen. De schijf bleek een bijdrage van 56.2% te leveren terwijl de nucleus en de halo verantwoordelijk zijn voor een bijdrage van respectievelijk 1.2% en 42.6%. Uit de bepaalde dichtheden kon verder ook de potentiaal berekend worden door gebruik te maken van de vergelijking van Poisson. In hoofdstuk 4 kon dan met behulp van kinematische data en de informatie over de dichtheid en de potentiaal, die eerder verzameld werd, een fit gemaakt worden aan de distributiefunctie van de galaxie. Hieruit beek dat enkel een goed resultaat bekomen kon worden als er ook tegenroterende sterren aanwezig zijn in de galaxie. Verder werd in hoofdstuk 5 de theorie besproken over hoe op zoek gegaan kan worden naar eigenmodes van de schijf. De stabiliteitscode voorhanden houdt echter geen rekening met een distributiefunctie die ook tegenroterende sterren beschrijft. Om de stabiliteitsanalyse op een later tijdstip toch nog uit te voeren, zou de code aangepast moeten worden.

89

90

HOOFDSTUK 6. BESLUIT

Bijlage A

De epicykelbenadering Hieronder zullen we de epicykelbenadering beschrijven. De afleiding is grotendeels dezelfde als deze die terug te vinden is in Binney & Tremaine (1989). In schijfvormige galaxieën bevolken de meeste sterren quasi-cirkulaire banen. Voor de schijf van onze galaxie zal dit niet anders zijn. De algemen bewegingsvergelijking voor een ster heeft de vorm, d2 r = −∇Φ(R, z), dt2

(A.1)

waarin Φ(R, z) de potentiaal is die symmetrisch is rond het z = 0 vlak en R, φ en z cilinderco¨ ordinaten voorstellen. We hebben, r = ReR + zez , ∂Φ ∂Φ ∇Φ = eR + ez . ∂R ∂z

(A.2) (A.3)

In cilinderco¨ ordinaten kunnen we schrijven, d2 r ¨ φ + zëz . ¨ − Rφ˙ 2 )eR + (2R˙ φ˙ + Rφ)e = (R dt2

(A.4)

Bovenstaande vergelijkingen combineren levert op, ¨ − Rφ˙ 2 = − ∂Φ , R ∂R d 2 ˙ ¨ ˙ 2R˙ φ + Rφ = 0 ⇒ R φ = 0, dt ∂Φ . z¨ = − ∂z

(A.5) (A.6) (A.7)

Vergelijking A.6 drukt het behoud van draaimoment van de z-component uit. We kunnen dus schrijven, ˙ Lz = R2 φ. (A.8) Hiermee kunnen we vergelijkingen A.5 en A.7 herschrijven als, ¨ = − ∂Φef f , R ∂R ∂Φef f z¨ = − , ∂z 91

(A.9) (A.10)

92

BIJLAGE A. DE EPICYKELBENADERING

waarin, L2z . 2R2 bereikt een minimum als, Φef f ≡ Φ(R, z) +

De effectieve potentiaal Φef f

∂Φef f ∂R ∂Φef f ∂z

= 0=

∂Φ L2 − z3 ∂R R

= 0

(A.11)

(A.12) (A.13)

De tweede vergelijking is overal in het equatoriaalvlak z = 0 voldaan en de eerste vergelijking heeft voor elke Lz een straal Rg als oplossing. Op deze straal krijgen we, Rg φ˙ 2 =

∂Φ (Rg , 0), ∂R

(A.14)

die een cirkelbaan definieert met cirkelsnelheid φ˙ op straal Rg . De effectieve potentiaal kan nu als een Taylorreeks ontwikkeld worden rond de cirkelbaan, ∂Φef f ∂Φef f (Rg , 0) + (z − 0) (Rg , 0) ∂R ∂z (R − Rg )2 ∂ 2 Φef f (z − 0)2 ∂ 2 Φef f + (R , 0) + (Rg , 0) g 2 ∂R2 2 ∂z 2 + O{3}. (A.15)

Φef f (R, z) ≈ Φef f (Rg , 0) + (R − Rg )

Verder kunnen we respectievelijk de epicykelfrequentie en de vertikale frequentie definiëren, κ2 (Rg ) = ν 2 (Rg ) =

∂ 2 Φef f (Rg , 0) ∂R2 ∂ 2 Φef f (Rg , 0) ∂z 2

(A.16) (A.17)

en kunnen we uitdrukking A.15 vereenvoudigen door bovenstaande uitdrukkingen in te vullen en rekening te houden met het feit dat A.12 en A.13 voldaan zijn voor R = Rg . Stellen we verder x = R − Rg , dan krijgen we, Φef f (R, z) ≈ Φef f (Rg , 0) +

x2 2 (z − 0)2 2 κ + ν 2 2

(A.18)

Deze vergelijking invullen in A.9 en A.10 levert, x ¨ = −κ2 x ⇒ R = Rg + AR cos(κt + αR ) 2

z¨ = −ν z ⇒ z = Az cos(νt + αz )

(A.19) (A.20)

op, waarbij AR en Az de radiale en vertikale amplitude en αR en αz willekeurige fases zijn. De baan van de ster kan dus gezien worden als een quasi cirkelvormige baan die bestaat uit een cirkelvormige baan met de superpositie in vertikale en radiale richtingen van de 2 oscillaties gegeven door vergelijkingen A.19 en A.20 . De cirkelbaan wordt beschreven door een gidscentrum met straal Rg met een constante hoeksnelheid Ω2 (Rg ) =

1 ∂Φ L2 (Rg , 0) = z4 . Rg ∂R Rg

(A.21)

93 De epicykelfrequentie wordt dan gegeven door, ∂Φef f ∂ ∂Φ L2z 2 κ (Rg ) = (Rg , 0) = − ∂R2 ∂R ∂R R3 (Rg ,0) L2z ∂ dΩ2 2 RΩ (R) − 3 = = Rg (Rg ) + 4Ω2 (Rg ). ∂R R (Rg ,0) dR

(A.22)

De hoeksnelheid kan geschreven worden als, Lz Lz = 2 2 R Rg (1 + x/Rg )2 Lz 2x ≈ 1− R2 Rg AR cos(κt + αR ), = Ωg − 2Ωg Rg

φ˙ =

(A.23)

waaruit volgt, φ = φ0 + Ωg t − 2Ωg

AR sin(κt + αR ). κRg

(A.24)

We hebben eerder al de x- en z-co¨ ordinaat bepaald, x = AR cos(κt + αR )

(A.25)

z = Az cos(νt + αz )

(A.26)

en nu kunnen we uit A.24 de y-coördinaat bepalen die loodrecht staat op de x- en zcoördinaat, 2Ωg AR sin(κt + αR ) κ = −Ay sin(κt + αR ).

y = −

(A.27)

De beweging in de z-richting is onafhankelijk van de beweging in de x- en y-richting. De ster zal op een ellips, genaamd de epicykelbeweging, bewegen in het xy-vlak rond het gidscentrum. De lengten van de assen van de ellips hebben de verhouding, AR κ = . Ay 2Ωg

(A.28)

We beschouwen nu de beweging van een ster die beweegt op een epicykelbaan gezien door een waarnemer die meebeweegt met het gidscentrum van de ster. Op verschillende tijdstippen tijdens de baan meet de waarnemer de maximumsnelheden x˙ en y˙ in de xen y-richting die respectievelijk κAR en κAy zijn. Deze waarnemingen zijn gerelateerd aan de galactische potentiaal. Dit is echter vrij omslachtig omdat de levensduur van de waarnemer een pak kleiner is dan de epicykelperiode. We kunnen dus een resultaat proberen te bekomen door uit te middelen over verschillende sterren die enkel in hun epicykelfase αR van elkaar verschillen. We krijgen dan, y˙ 2 x˙ 2

=

(vφ − vc (Rg ))2 2 vR

=

4Ω2g . κ2

(A.29)

94

BIJLAGE A. DE EPICYKELBENADERING

Dit is echter geen praktische procedure aangezien we de straal Rg van de sterren niet kennen. We kunnen enkel de snelheid op een straal R meten, dus het enige wat we kunnen meten is, vR (R) en vφ (R) − vc (R) op een straal R tot het middelpunt van de galaxie. We kunnen nu voor beide metingen een uitdrukking berekenen, vφ (R) − vc (R) = R(φ˙ − Ω) = R(φ˙ − Ωg + Ωg − Ω) dΩ ˙ ' R (φ − Ωg ) − x dR R 2Ωg dΩ − , = −Rx Rg dR R

(A.30)

waarin we in de laatste stap gebruik gemaakt hebben van vergelijking A.23 . Voor kleine x kunnen we Ωg /Rg vervangen door Ω/R. Bovenstaande vergelijking kan dan herschreven worden als, dΩ vφ (R) − vc (R) = −x 2Ω + R dR R dΩ x 2 4Ω + 2RΩ = − 2Ω dR 2 R x dΩ = − 4Ω2 + R 2Ω dR R = −

xκ2 (R) , 2Ω

(A.31)

waarin we in de laatste stap vergelijking A.22 gebruikt hebben. Tenslotte hebben we dus, (vφ (R) − vc (R))2 =

κ4 (R) 2 x . 4Ω2

(A.32)

We hebben in de radiale richting, 2 =x vR ˙ 2 = κ(R)2 x2 = κ2 (R)x2 .

(A.33)

We krijgen dan voor de verhoudingen van de dispersies, σφ2 2 σR

=

(vφ (R) − vc (R))2 2 vR

=

κ2 . 4Ω2

(A.34)

Bibliografie [1] Baes M., Galactische Sterrenkunde, Cursus, Gent, UGent, 2009-2010. [2] Baes M., Astrophysical Simulations, Cursus, Gent, UGent, 2010-2011. [3] Binney J. & Tremaine S., Galactic Dynamics, Princeton, NJ: Princeton University Press, 1987. [4] Burkert A., The structure of dark matter halos in dwarf galaxies, The Astrophysical Journal, 447:L25-28, 1995. [5] Clem J., Landolt A., The eclipsing binary PG1323-086A, Publications of the Astronomical Society of the Pacific, Vol. 122, No. 887, 27-28, 2010. [6] Dejonghe H., A quadratic programming technique for modeling gravitating systems, The Astrophysical Journal, 343: 113-124, 1989. [7] Dejonghe H., Dynamica van galaxieën, Cursus, Gent, UGent, 2010-2011. [8] De Rijcke S., Dejonghe H., Zeilinger W., Hau G., Embeddd disks in Fornax dwarf elliptical galaxies, A&A 400, 119-125, 2003. [9] De Rijcke S., Fysica van galaxieën, Cursus, Gent, UGent, 2011-2012. [10] Dury V., Stabiliteit van schijfvormige sterrenstelsels, Thesis, Gent, UGent, 20042005. [11] Dury V., De Rijcke S., Debattista V., Dejonghe H., Global m=1 instabilities and lopsidedness in disc galaxies, Mon. Not. R. Astron. Soc. 387, 2-12, 2008. [12] Jerjen H., Kalnajs A., Binggeli B., IC3328: A dwarf elliptical galaxy with spiral structure, A&A 358, 845-849, 2000. [13] Koleva M., Prugniel Ph., Bouchard A., Wu Y., ULySS: a full spectrum fitting package, A&A 501, 1269-1279, 2009. [14] Le Borgne D. et al., Evolutionary synthesis of galaxies at high spectral resolution with the code PEGASE-HR, A&A 425, 881-897, 2004. [15] Lin C., Shu F., On the spiral structure of disk galaxies, Astrophysical Journal, vol. 140, 646, 1964. 95

96

BIBLIOGRAFIE

[16] Lisker T., Fuchs B., On the nature of IC 3328, an early-type dwarf galaxy with week spiral structure, A&A 501, 429-435, 2009. [17] Moore, B., Lake G., Katz N., Morphological transformation from galaxy harassment, Astophysical Journal, 495, 139-151, 1998. [18] Nakos. T. & Baes M., Observational techniques in astronomy, Cursus, Gent, UGent, 2009-2010. [19] Salucci P., Burkert A., Dark matter scaling relations, The Astrophysical Journal, 537:L9-L12, 2000. [20] Simien F., Prugniel Ph., Kinematical data on early-type galaxies, A&A 384, 371382, 2002. [21] Vazdekis A. et al., Ecolutionary stellar population synthesis with MILES - I. The base models and a new line index system, MNRAS 404, 1639-1671, 2010.

Datareductie en stabiliteitsanalyse van IC3328

Recommend Documents