Thursday, December 30, 2010

Brazil to Join ESO

Uz najakou dobu se o tom susklo. A tesne pred koncem roku se to skutecne stalo --- Brazile se stala celkove 15. a zaroven prvni ne-evropskou clenskou zemi Evropske Jizni Observatore (European Southern Observatory -- ESO) -- oficialni verejne prohlaseni ESO.

Order and Progress

Friday, December 3, 2010

XSPEC

... je balík na fitovanie spektier vo vysokých energiách (alternatívou k nemu môže byť napr. spex, ale o ňom možno až dakedy nabudúce) a tiež sme o ňom hovorili na brainstormingu. Tu je nejaký stručný prehľad toho, čo sa hovorilo. Opäť najlepším zdrojom informácií je manuál. XSPEC samotný má v sebe mocný príkaz help.
XSEPC> help plot
XSPEC> help model apec
Prvvý príklad vypíše nápovedu k príkazu plot (duh!), druhý odhalí, čo presne popisuje model apec, aké sú jeho parametre a podobne.

Inštalácia

XSPEC je súčasťou HEASoftu, návod je teda tu.

Monkeying around, data

S XSPECom si, ako s mnohými inými, dopisujeme v príkazovom riadku, čo je cool lebo sa dajú jednoducho vyrábať fitovacie skripty. Balík pozná obrovské množstvo príkazov, výhodou je aj to, že stačí napísať prvých pár znakov, aby sa príkaz dal interpretovať (napr. namiesto setplot energy stačí setp e).

Spektrá do XSPECu zadávame príkazom
XSPEC> data meno_spektra
Často ale používame viac spektrier naraz, v tom prípade ich treba roztriediť do skupín - tie rozlišujeme pomocou čísiel oddelených dvojbodkou, prvé označuje data set (spektrá budú fitnuté rovnakým modelom, druhé plot set. V praxi teda
XSPEC> da 1:1 spektrum1 1:2 spektrum2 2:3 spektrum3
kde prvé dve spektrá popisuje model s rovnakými parapetrami, pre tretie sa parametre menia.
FITS súbor so spektrom by mal v headeri mať nastavené keywords BACKFILE, RESPFILE a ANCRFILE, aby si tieto súbory vedel XSPEC načítať spolu so spektrom. V prípade, že tam vôbec nie sú, odmietne spoluprácu úplne, ak sú nastavené na "none", treba súbory načítať príkazmi backgrnd, response a arf, ktoré sa k danému spektru priradia pomocou rovnakých číselných identifikátorov ako pri načítavaní spektier.

V ďalšom kroku treba vybrať časti spektra, ktoré chceme fitovať (vo veľmi vysokých alebo veľmi nízkych energiách nemusí byť kalibrácia detektoru úplne istá a podobne). Na toto slúžia príkazy ignore a notice.
XSPEC> ig bad
XSPEC> ig 1:12-27
XSPEC> no *:**-0.4
XSPEC> ig 2:7.-**
V prvom príklade ignorujeme všetky časti spektra označené ako bad. V ostatných označuje číslo pred dvojbodkou data sety (hviezdička znamená všetky data sety). Za dvojbodkou je určený interval - celé čísla XSPEC interpretuje ako channely, reálne sú energie v keV. ** znamenajú "kam to až ide", čiže najvyššiu/najnižšiu možnú energiu/channel.

Monkeying around, model

Spektrum modelujeme príkazom model. XSPEC obsahuje množstvo rôznych modelov (vypísať sa dajú spustením model bez parametrov), ktoré sú rozdelené do niekoľkých kategórií - aditívne, multiplikatívne atď. Model sa definuje pohodlným spôsobom analogickým k zápisu rovnice. Príkladom môže byť model známeho objektu MCG-6-30-15 (AGN s relativisticky rozšírenou čiarou železa vo vnútornej časti akrečného disku). Jeho spektrum obsahuje aditívne prvky - blackbody zložku (bbody), powerlaw (powerlaw) a čiaru s Laorovým profilom (laor), to celé, z dôvodu, že žiarenie na ceste k nám prechádza prostredím s absorpciou, prenásobené absorpčnou zložkou (phabs). Zápis modelu potom vyzerá
XSPEC> mo phabs*(bbody+po+laor)
XSPEC sa ďalej spýta na prvotné odhady hodnôt parametrov a to toľkokrát, koľko máme data setov. Parametre môžeme určiť niekoľkými číslami oddelenými čiarkami (odhad hodnoty, veľkosť kroku pri iterácií, maximálna a minimálna hodnota a podobne). V prípade, že chceme ponechať parametru defautlnú hodnotu, píšeme namiesto čísla /, ak chceme defaultné hodnoty pre všetky parametre od aktuálneho ďalej, píšeme /*.

Po zadaní prvotných odhadov XSPEC vypíše tabuľku s modelom a tiež hodnotu štatistiky, ktorá popisuje, ako dobre model popisuje spektrum. Parametre môžu mať rozne vlastnosti. V ďalšom celé čísla pomenovávajú jednotlivé parametre, reálne čísla sú konštanty. Parametre možno fixovať na danej hodnote, ktorá sa pri fitovaní nebede meniť (freeze 3 zafixuje parameter č.3) a opäť uvoľňovať (thaw). Môžeme im priraďovať nové hodoty a taktiež ich medzi sebou zväzovať a rozväzovať.
XSPEC> newpar 4 1.5
XSPEC> newpar 5 = 0.75*4
XSPEC> untie 5
XSPEC> newpar 10 1e-3,-1
V prvom príklade priradíme parametru 4 novú hodotu. Parameter 5 k nemu následne uviažeme s tým, že bude mať hodnotu troch štvrtím parametru 4, čo sa berie do úvahy aj pri fitovaní. V teťom príklade sme si to rozmysleli a parameter 5 opäť osamostatnili. Posledný príklad ukazuje, ako v jednom kroku dať parametru novú hodnotu a zároveň ho zafixovať (v princípe by sme miesto -1 mohli použíť akékoľvek záporné číslo).

S týmto všetkým sa môžeme pustiť do fitovania. Existuje niekoľko typov štatistík, ktoré pri tom môžeme použiť, medzi ktorými prepíname príkazom statistic. Pre zdroje s nízkou povrchovou jasnosťou (plošné zdroje) dáva C-štatistika lepšie výsledky. Ak je v spektre dôležité aj chemické zloženie pozorovaného porstredia, ktoré môžeme prirovnávať k slnečnému, dá sa príkazom abund prepínať medzi jednotlivými modelmi slnečného zloženia. Samotné fitovanie spustíme príkazom
XSPEC> fit
(dajú sa k nemu pridať aj nejaké parametre). Defaultné nastavenie sa po každých 10 iteráciách spýta, či má pokračovať, tohto sa dá zbaviť príkazom query yes, po ktorom si na všetečné otázky XSPEC odpovie sám. Chyby nafitovaných parametrov uvedené v tabuľke modelu sú len orientačné, spočíta ich príkaz error, ktorému je dobré zadať level of confidence (v sigma) a nevyhnutné zadať, ktoré parametre nás zaujímajú. Čiže
XSPEC> err 1. 3,4,6
spočíta 1 sigma chyby pre parametre 3,4 a 6.

Mokneying around, plot

Príkazy na vykresľovanie v XSPECu pochádzajú z balíka PGPLOT. Nastavené parametre obrázkov sa zachovávajú až kým nie sú zmenené, to znamená, že ak graf nakreslíte nejakým šialeným príkazom, stačí pri jeho prekreslení napísať plot a všetky šialené parametre sa zachovajú.

Príklady:
XSPEC> cpd /xw
XSPEC> setp e
XSPEC> pl ld res
XSPEC> pl
Na začiatok treba nastaviť, kam sa má obrázok kresiť. Slúži na to príkaz cpd (change plot device), /xw je pre XWindows, /cps kreslí to postscript súboru (pred lomeno sa dá napísať meno output súboru, cpd mojespektrum.ps/cps). setplot udáva vlastnosti obrázku, e je pre energy, čiže miesto channelov bude xová os v energiách. pl ld res hovorí PGPLOTU, že má kresliť v log-log a v spodnom paneli obrázku vykresliť reziduá (rozdiel dát a modelu). Posledný pl urobí presne to, čo príkaz o riadok vyššie, ale už nemusíme písať parametre.

Monkeying around, zvyšok
Na záver niekoľko technických poznámok. Všetko, čo do sa do XSPECu napíše sa dá uložiť do logu príkazom log meno_logu, logovanie sa vypína log none. Príkazom show vieme zobraziť rôzne informácie - paramtere modelu (sho pa), voľné parametre (sho fre) a podobne. Ukladanie modelu je možné pomocou save model.

Funguje to aj naopak - do XSPEC vieme čítať príkazy, modely atď z textových súborov. XSPEC má rád súbory s koncovkou .xcm, samotné čítanie funguje pomocou @
XSPEC> @subor.xcm
Conclusion

Toľko stučný úvod do XSPECu, opäť sú otázky a komentáre vítané. Príkazy pri ich prvom objavení som sa snažil písať v ich celej forme, no ďalej som používal skrátené verzie - verím, že to nebolo veľmi mätúce :P

Thursday, December 2, 2010

XMM-SAS

(...hlavne pre tých, čo to v pondelok nedali na brainstorming ale aj všetkym ostatným, ktorým to mohlo v pondelok pripadať trochu zmätené...)

XMM-SAS (XMM-Newton Science Analysis System) je balík procedúr a kadečoho iného, čo nám pomáha odhaľovať tajomstvá oblohy tak, ako ich vidí XMM-Newton. Všetko, čo tu píšem, je v tej či onej podobe nájditeľné na XMM-Newton stránkach, v kadejakých manuáloch a tak, ten najužitočnejší je tu. Dôležité je uvedomiť si, že analýza XMM dát, ktorú som robil sa venovala plošným zdrojom (ako sa preloží extended source?) a tá funguje trochu inak ako analýza bodových zdrojov - oveľa väčší dôraz sa kladie na modelovanie pozadia (ktoré plošný zdroj so slabou plošnou jasnosťou oveľa viac ovplyvňuje) a odpadajú problémy ako napríklad pileup a podobne. Preto treba pri vlastnej analýze byť na pozore.

Inštalácia

...vyžaduje predinštalovaný sofrware a umiestnenie (nalinkovanie) perlových bináriek tam, kde ich software očakáva, inak je pomerne jednoduchá a priamočiara (teda len do chvíle, keď sa rozhodnete SAS sami skompilovať, čo nie je až tak triviálne a samotní tvorcovia od toho tak trochu odrádzajú). Existuje niekoľko predkompilovaných verzií, treba si len vybrať svoju najobľúbenejšiu a stiahnuť . Inštalácia samotná sa skladá z dvoch krokov - stiahnutý archív treba rozbaliť,
$ tar xvf "meno_archívu"
a spustiť
$ ./install.sh
Nakoniec treba zadefinovať nejaké premenné, o čo sa stará skript, ktorý treba nasourcovať pred každým spustením (napr. pridaním kódu do .bashrc).
$ . /cesta/xmmsas_20100423_1801/setsas.sh
Detaily inštalácie sú detailne tu.

Kalibračné data, ktoré sú k analýze tiež nevyhnutné, sú stiahnuteľné tu alebo najrýchlejšie pomocou príklazu
$ rsync -v -a --delete --delete-after --force --include='*.CCF' --exclude='*/' xmm.esac.esa.int::XMM_RED_CCF /ccf/
kde posledné /ccf/ treba nahradiť adresárom, kde chceme súbory uložiť. SAS sa o nich dozvie pomocou premennej SAS_CCFPATH, ktorú nastavíme ručne tak, aby naň ukazovala.

XMM-Newton

...má na svojej palube niekoľko inštrumentov, OM (Optical Monitor) ani RGS1 a RGS2 (Reflective Grating Spectrometer) nás v tejto chvíli nebudú zaujímať. Sústrediť sa budeme na kameru s ultimátnym menom EPIC (European Photon Imaging Camera) a jej tri detektory, MOS1, MOS2 a pn.

Monkeying around, initialization

Dáta pre analýzu sú dostupné na XSA stránkach, kde nájdeme niekoľko ciest, ako ich stiahnuť. Vo výsledku sa nám v počítači objaví jeden alebo viac archívov typu xxxxxxxxxx.tar.gz nazvaných podľa identifikátora pozorovania (ObsID). Po jeho kompletnom rozbalení (samotný archív v sebe skrýva ešte jeden .TAR, bacha na to), najlepšie do samostatného adresára, aby sme si nepomiešali surové dáta s výsledkami - napr. adresár s názvom odf - definujeme premennú SAS_ODF, ktorá ukazuje práve na tento adresár
$ export SAS_ODF=/cesta/odf/
Teraz možeme spustiť
$ cifbuild
ktorý vytvorí súbor ccf.cif (Current Calibration File). Najlepšie je spustiť ho priamo v odf/ adresári. Na ccf.cif ukazuje ďalšia premenná
$ export SAS_CCF=/cesta/ccf.cif 
(napr. $ export SAS_CCF=$SAS_ODF/ccf.cif, ak to robíme podľa toho, čo píšem vyššie). Nakoniec spustíme
$ odfingest
ktorý prekope celé pozorovanie a zapíše, čo všetko obsahuje, do .SAS súboru. Pomocou tohto súboru sa ostatné procedúry orientujú v pozorovaní a ak ho nenájdu, budú toho názoru, že pozorovanie neobsahuje žiadne dáta.

Monkeying around, data products

Teraz sme pripravení z pozorovania vytiahnuť najsurovejšie eventfily, teda súbory so všetkým, čo jednotlivé detektory napozorovali. Na tomto mieste je dobré si - opäť z dôvodu zachovania poriadku - vytvoriť podadresáre pre každý detektor.
$ mkdir mos
$ mkdir pn
V adresári mos/ teda spustíme
$ emchain
(bez parametrov). V adresári sme týmto vytvorili veľa obskurne nazvaných súborov, pre nás najdôležitejšie sú *M1*MIEVLI* a *M2*MIEVLI*, ktoré si môžeme poyrieť napr v ds9.Pre pn sa o niečo podobné stará
$ epchain
no treba si dávať pozor - pri rôznych pozorovaniach zaznamenáva pn detektor určitý zlomok OOT (Out Of Time) eventov, ktoré treba od pozorovania odčítať. Tu je treba pozrieť sa do SAS a tiež epchain manuálov pre správne parametre a spôsoby odčítania. Jedným zo spôsobov odčítania je použitie utility farith z heasoftu.

Monkeying around, filtering

Vysoká obežná dráha XMM-Newtona spôsobuje, že detektory sú ovplyvňované slnečným vetrom. Detaily o konkrétnych časticiach sú v manuáloch, pre nás je teraz podstatné, že signál zo slnečného vetra nevieme rozumne oddeliť od pozorovania a musíme sa zbaviť časových intervalov, kedy bolo slnečné pozadie zvýšené. Čiernou skrinkou na celý proces sú procedúry
$ mos-filter
$ pn-filter
no tieto veci sa dajú robiť ručne.

Po prvýkrát sa stretneme s procedúrou evselect, ktorá slúži na prehľadávanie a výber eventfilov. Existuje nekonečne veľa parametrov, ktoré táto procedúra pozná, preto stojí za to zoznámiť sa s manuálom. V tomto prípade chceme urobiť svetelnú krivku vo vysokých energiách s nejakým rozumným zbinovaním času. Pre MOS detektory je to najlepšie nejak takto
$ evselect table=meno_eventfilu withrateset=yes rateset=meno_sv_krivky timecolumn=TIME \
timebinsize=100 maketimecolumn=yes \
expression='#XMMEA_EM && (PI in [10000:12000]) && (PATTERN<=12)'
expression udáva charakteristiky eventov, čo nás zaujímajú. #XMMEA_EM vyberá len good eventy, PI určuje interval energii v eV, PATTERN určuje, aké typy eventov sa majú brať do úvahy (fotóny môžu na detektore vytvoriť rôzne obrazce - nie je to limitované na jeden fotón/jeden pixel - pričom kalibrované sú všetky po 12ty). Pre pn použijeme analogický príkaz, binovanie bude ale nižšie (pokus-omyl), v expression bude XMMEA_EP, energia 10-14keV a pattern len po 4.

Svetelná krivka sa dá zobraziť napr. v fplote.

Filtrovanie samotné prebieha určením GTI (Good Time Intervals), teda časových intervalov, ktoré chceme používať. Procedúra tabgtigen vie vyrobiť GTI súbor, ktorý môžeme použit v evselecte. Príkaz vyzerá
$ tabgtigen table=meno_sv_krivky gtiset=meno_gti_súboru \
expression=COUNTS.gt.0.909262.and.COUNTS.lt.23
pričom spodnú a hornú hranicu v expression môžeme určite viacerými spôsobmi. Jednoduchší ale menej vedecký spôsob je stanoviť hodnotu priamo zo svetelnej krivky tak, aby približne konštantné pozadie zostalo zachované ale flary zmizli. Pomocou procedúry fhisto ale tiež vieme zostrojiť zo svetelnej krivky histogram
$ fhisto meno_sv_krivky meno_histogramu COUNTS 1
Pokiaľ sme použili správnu veľkosť časových binov pri zosrojovaní sv. krivky, mali by sme dostať Poissonovské rozloženie. Toto môžeme s prižmúrením oka preložiť Gaussovou krivkou (v fplote príkaz mo gaus) a brať do úvahy napr. úseky s pozadím plusmínus 3 sigma. Samotné filtrovanie pomocou evselect potom vyzerá
$ evselect table=meno_eventfilu withfilteredset=true destruct=yes keepfilteroutput=true \
expression="GTI(meno_gti,TIME)" filteredset=filtrovany_evetfile
Po tomto všetkom treba rovnako odfiltrovať úseky s vysokým pozadím v nízkych energiách 0.3-1keV.

Monkeying around, spectra

S odfiltrovanými eventfilami môžeme konečne začať robiť nejakú vedu. evselect nám môže vytvoriť obrázky (withimageset=true) a podobne, nás ale budú zaujímať spektrá.

Na začiatok je dobré uvedomiť si, spektrum čoho vlastne chceme robiť. Najčastejšie asi pôjde o nejaký bodový zdroj (AGN, supernova,...), no môže tak isto ísť o horúci plyn ICM v kopách galaxií prípadne niečo úplne iné (spektrá z rôznych častí jednej galaxie a pod.). eveselectom musíme teda správne vybrať eventy z eventfilu. Manuál popisuje pomerne značné množstvo vecí, ktoré sa dajú napísať do expression - môžeme vykrajovať kruh, anulus, obdĺžnik a podobne, môžeme vykrajovať všetko okrem daného kruhu, náš útvar môžeme určovať v rôznych druhoch súradníc. V tomto prípade sa teda všeobecne nedá napísať nič a závisí na konkrétnom prípade. Podstatnou poznámkou ale je, že pomocou ewavelet vieme v pozorovaní vyhľadať bodové zdroje.

Predpokladajme teda, že náš eventfile už obsahuje iba dáta, ktoré nás zaujímajú. Pre samotnú spektrálnu analýzu potrebujeme 4 súbory - samotné spektrum, spektrum pozadia, response maticu (RMF) a anciliary response file (ARF).

Spektrum vyrábame pomocou evselect s parametrom withspectrumset=true, teda napr.
$ evselect table=meno_eventfilu withspectrumset=true withspecranges=true \
energycolumn=PI specchannelmin=0 specchannelmax=11999 spectralbinsize=15 \
updateexposure=yes spectrumset=meno_spektra \
expression="(FLAG==0) && (PATTERN<=12) \ && (PI in [300:12000])"
pre MOS a analogicky pre pn (specchannnelmax=20479, PATTERN<=4, PI in [300:14000]).

Spektrum pozadia je špecifické pre daný problém. Pre bodové zdroje stačí vyrobiť spektrum z inej časti detektora a označiť ho ako pozadie, pre plošné zdroje to už ale také jednoduché nie je - dá sa vyrobiť naškálovaním closed-filter observations (niečo ako dark framy) prípadne z oblastí čipu ležiacixh mimo field of view.

Response matica hovorí, ako sa channely (charakteristika eventov, ktoré vyjdu z detektoru) prekladajú do energií. Vyrába sa procedúrou rmfgen, ktorá potrebuje detektorovú mapu. Je to jednoducho obrázok v detektorových súradniciach, ktorý vyrobíme pomocou evselect s withimageset=true xcolumn=DETX ycolumn=DETY, zvyšné parametre obdobne ako pri normálnej výrobe obrázku. Response maticu potom vytvoríme poocou
$ rmfgen spectrumset=meno_spektra rmfset=meno_rmf detmaptype=dataset \
detmaparray=meno_detektorovej_mapy extendedsource=yes \
badpixlocation=povodny_eventfile
Ak sa nejedná o plošný zdroj, treba upraviť príslušný parameter.

Anciliary response file udáva efektívnu plochu detektoru pre rôzne energie. Vyrába sa pomocou arfgen.
$ arfgen spectrumset=meno_spektra arfset=meno_arf \
withrmfset=yes rmfset=meno_rmf \
extendedsource=yes badpixlocation=meno_eventfile \
detmaptype=dataset detmaparray=meno_detektorovej_mapy
Pre neplošné zdroje treba znovu upraviť príslušný parameter.

Úplne nakoniec je vhodné upraviť header hlavného spektra a do keywords BACKFILE, RESPFILE a ANCRFILE zapísať príslušné súbory. Procedúrou grppha je tiež vhodné zoskupinovať spektrum tak, aby v kažom channeli bol aspoň jeden count - fitovacie programy ako XSPEC to majú hrozne rady :)

Celkom nakoniec je dobré spomenúť ďalšie čierne skrinky, tieto sú ale určené pre plošné zdroje. Procedúry mos-spectra, mos_back, pn-spectra a pn_back, podrobne popísané v cookbooku od Steva Snowdena, vytvoria spektrá, vymodelujú spektrá pozadia, vytvoria RMF a ARF a obrázky sami od seba.

Conclusion

Týchto pár poznámok si nerobí nárok na úplnosť ale niekomu môže pomôcť. V najbližšom čase by som rád urobil niečo podobné pre XSPEC, tým pádom by sme tu mali všetko, čo som hovoril na brainstormingu. Otázky a komentáre sú samozrejme vítané.