Modelling of the "Pi of the Sky" detector
UUniversity of WarsawFaculty of Physics
Modelling of the “Pi of the Sky”detector
Lech Wiktor Piotrowski
PhD thesis written under supervisionof prof. dr hab. Aleksander Filip ŻarneckiWarsaw, 2011 a r X i v : . [ a s t r o - ph . I M ] O c t bstract The ultimate goal of the “Pi of the Sky” apparatus is observation of optical flashes ofastronomical origin and other light sources variable on short timescales, down to tensof seconds. We search mainly for optical emission of Gamma Ray Bursts, but also forvariable stars, novae, blazars, etc. This task requires an accurate measurement of thesource’s brightness (and it’s variability), which is difficult as “Pi of the Sky” single camerahas a large field of view of about ◦ × ◦ . This causes a significant deformation of apoint spread function (PSF), reducing quality of brightness and position measurementwith standard photometric and astrometric algorithms. Improvement requires a carefulstudy and modelling of PSF, which is the main topic of the presented thesis. A dedicatedlaboratory setup has been created for obtaining isolated, high quality profiles, which inturn were used as the input for mathematical models. Two different models are shown:diffractive, simulating light propagation through lenses and effective, modelling the PSFshape in the image plane.The effective model, based on PSF parametrization with selected Zernike polynomialsdescribes the data well and was used in photometry and astrometry analysis of the framesfrom the “Pi of the Sky” prototype working in Chile. No improvement compared tostandard algorithms was observed in brightness measurements, however more than factorof 2 improvement in astrometry accuracy was reached for bright stars. Additionally, themodel was used to recalculate limits on the optical precursor to GRB080319B – a limithigher by . m compared to previous calculations has been obtained.The PSF model was also used to develop a dedicated tool to generate Monte Carlosamples of images corresponding to the “Pi of the Sky” observations. The simulator allowsfor a detailed reproduction of the frame as seen by our cameras, taking into account PSF,electronics and mechanical fluctuations, etc. A comparison of photometry performed onreal and simulated data resulted in very similar results, proving the simulator a worthytool for future “Pi of the Sky” hardware and software development. ontents hapter 1Introduction The sky has fascinated mankind from the very beginning. It’s main inhabitants - Sun,Moon and stars - inspired most ancient cultures to think about their place in the Universeand the passage of time, for behaviour of these celestial bodies was regular and predictable.However, as soon as man started to look into the sky with something more than a pureamazement, he found among the stars some strange ones, not following the rules of themost. These were planets. From time to time the wild blue yonder was raided by someother visitors, such as comets, or much more often – meteors. Those dwellers of heavensshaped the knowledge about the world above our heads for millennia.During passage of centuries, the sky, recognized for a long time as immutable anda rather peaceful place slowly occurred to be a home to many violent events. Amongthese were those more spectacular, such as supernovae or less destructive such as novaeor flare stars. However, quite recently, it were the Gamma Ray Bursts that shocked theastronomers and introduced a set of new ideas to this branch of science.GRBs are perhaps the most energetic phenomena known to man, emitting enormousamounts of energy mainly in gamma-rays, but also in other wavelengths of electromag-netic radiation. More than 4 decades of observation since discovery in 1969 by a pair ofVela satellites showed, that they are also one of the most remote type of objects, appear-ing always at cosmological distances. A few outbursts are present on the sky per week,but their position is totally random. Thus detection can be performed only by satellitegamma-ray detectors covering very large fractions of the sky. To perform simultaneousobservations also in optical band, similar approach has to be introduced to on-groundoptical observatories. Up to now most of the telescopes focused on observing a very smallfraction of the sky with a big magnification and thus were nearly useless for such simul-taneous observations. “Pi of the Sky” introduces a completely new approach – constantmonitoring of a very large part of the sky, aimed at detection of bright transients comingfrom random directions.This work was prepared as a part of the “Pi of the Sky” project. Constant sky monitor-ing requires cameras with a very large field of view (for astronomical standards). Lensesrequired to obtain such parameters introduce large deformations to the image, increas-ing with the distance from the frame centre (optical axis). Study of these deformations,described by the point spread function (PSF) of the detector, is the main subject of thisthesis.The first chapter describes short time-scale variable phenomena on the sky, focusingon Gamma Ray Bursts – the main field of interest of the “Pi of the Sky” experiment.The general idea of the experiment is presented in the second chapter, along with theconstruction of the prototype and the concept of the final system. Analysis software, main5cientific results and uncertainties considered in the data analysis are briefly discussed.The fourth chapter shows results of the laboratory measurements of the point spreadfunction as well as of some other characteristics of the detector. Attempts to reproducethese results with mathematical models are presented in the subsequent two chapters: inchapter 5 via physical modelling of the light passage through lenses and in chapter 6 viaparametrization of the PSF shape on the frame.The final outcome of this work are the applications of the resulting detector responsemodel. The first being the development of a new algorithm for measuring stars (or otherobjects) brightness and position on the frame, as well as search for optical precursor toGRB080319B, presented in chapter 7. The other application, shown in chapter 8, is theprocedure of a realistic simulator of the “Pi of the Sky” frame, featuring deformed starsshapes, signal fluctuations and numerous other features, which can be used for futurehardware and software development .All results presented in chapter 4 to 8 were obtained by the author of this thesis. Simulation with Monte Carlo techniques becomes a common tool for improving our understandingof the complex detector systems. hapter 2Short time-scale phenomena in the sky The first that we know of, to question the immutability of the sky and to write aboutan appearance of a new resident of the firmament was a Greek astronomer Hipparchus(although such events are reported to be observed by Chinese about 1500 BC). In IIcentury BC Hipparchus noticed a “New Star” – gr. Nova (sec. 2.2.1) – that he was surehe had not observed before. That discovery, a shock to the Greek community, whichthough of a sky as of a constant sphere, probably pushed Hipparchus to the creation ofone of the first star catalogues, to allow future generations to discover other new stars.The apparent immutability of the heavens was also questioned by the Chinese. In185 AD they observed a supernova (sec. 2.2.2) SN 185[1], which remained in the skyfor 8 months . Although we hardly know anything about the astronomic observationsin the earlier times, this supernova is believed to be the first recorded by the mankind.Following that other supernovae were rarely discovered, every few hundreds years or so,like SN 1006[3] – the brightest one so far, SN 1054[4] (the one that formed Crab nebula),SN 1572[5] and SN 1604[6], the last one being the latest seen in the Milky Way.The SN 1572 discovered by Tycho Brahe had forced reinvention of the Hipparchusidea, that the night sky beyond planets and Moon is not immutable. This had begun arevolution in European astronomical thinking. Shortly afterwards, in year 1584, an ideathat other stars are other Suns placed far from the Earth was reintroduced by GiordanoBruno, and given physical basis in the next century by Isaac Newton, who explainedabsence of the gravitational pull on the Solar System by equal distribution of stars in thesky.It was not long before the first discovery of a variability of a known star. In year1638 Johannes Holwarda found that the star Omicron Ceti (later named Mira) had a 11months pulsation cycle[7]. Next, variations of the luminosity of star Algol were observedby Geminiano Montanari in 1667. By 1786 ten variable stars were known and the numberstarted to increase rapidly after 1850, induced by an improvement of observational tech-niques. Along with all those discoveries, certain understanding of the nature of celestialbodies was born and is still growing, boosted by invention of quantum physics in the firsthalf of XX century.Other milestones in observation of astronomical objects, including variable stars andother short time-scale phenomena, was introduction of photography in late XIX centuryand CCD imaging in late XX century. However, stepping into XX century, one cannot The first documented Chinese discovery of a star, that could be a supernova comes from 352 BC[2].Some historical sources claim, that suspected supernova was first seen in 2241 BC. Although it was a supernova, at this time it was called “nova”, and ambiguity remained up to year1930.
The main inspiration for the “Pi of the Sky” project, described in more detail in chap.3, are the so-called Gamma Ray Bursts (GRBs) – cosmic explosions radiating mostly ingamma-rays (thus the name), on timescales of milliseconds to tens of seconds[8]. Someof their energy is also radiated in optical band, as well as radio and X-rays – otherwisethey would not be of any interest to the optical experiments. They still remain one ofthe most fascinating and challenging subjects in modern astrophysics, and their mysterybegun with their discovery.
In times when many, if not most of the discoveries in science are born from the jointeffort of experiment and theory, and thus are in some way expected, surprises hold aspecial value and nearly always focus attention of scientific world for at least a short time.This was the case of Gamma Ray Bursts, which were the first extra-terrestrial (excludingSun) sources of gamma rays discovered by man. More so, the discovery was completelyunexpected.In year 1963 United States of America launched the first pair of military satellitesVela equipped with gamma-ray detectors. Their task was to monitor the obedience of thenuclear tests ban treaty by other countries. The treaty banned both the terrestrial andextra-terrestrial tests, therefore Vela satellites were placed on an orbit of about 120000km above the Earth surface – high enough to detect possible explosions on the dark sideof the Moon. Four years passed till the first detection of gamma ray explosion on nd ofJuly, 1967 (fig. 2.1). The outburst contained no known signatures of a nuclear explosion.It’s detection was possible with the 3rd generations of Vela satellites, launched in 1965,which were the first generation to store such events for further analysis.However, th and th generation of Vela satellites, launched in 1969, were requiredto prove that the source of gamma ray bursts was not localized on Sun, Moon or theEarth’s orbit. In the period 1969-1972 satellites detected 16 gamma ray bursts (GRBs)and results were published in Astrophysical Journal Letters in 1973 [9]. Additionally, datafrom Soviet Konus satellites confirming the detection of GRBs were published in 1974.These discoveries begun the era of GRB science, still one of the most intriguing branchesof astrophysics.The main tool for gamma ray bursts research up to 1991 was the InterPlanetaryNetwork (IPN)[10] formed in 1976 . At the beginning it consisted of existing spacecrafts IPN is working up to this day, although spacecrafts forming it have changed. igure 2.1: Historical lightcurve of the first GRB ever detected. The detection was performedby the Vela satellites[9]. equipped with gamma ray detectors, already studying Sun and solar planets. The networkallowed triangulation of GRB coordinates with uncertainty of a few arcminutes. This wasenough to prove that GRBs are not counterparts of any of the known space sources ofradiation, like, for example, X-ray sources. However IPN could not provide any insightsto the distance between Earth and outbursts and thus the origin and nature of GRBsremained a speculation for the next 15 years.Part of the scientific community shared an idea that GRBs origin at cosmologicaldistances. However, very high intensity of the radiation brought most scientists to believethat the phenomenon occurs much closer, at Galactic distances or, at the beginning, evenat the Solar System boundaries. The weight shifted across the years between supportersof both theories, often due to misleading observations.One of the earliest hypothesis supporting Galactic origin was an idea of a neutronstar bombarded by some kind of an object, probably a comet. A lack of any commonpattern in GRBs lightcurves and an estimation of number of neutron stars and comets inthe Milky Way made this a preferred explanation for some period of time. A significantargument for the supporters of cosmological distances came on the th of March 1979 withthe discovery of a burst coincident in position with Supernova N49 in the Large MagellanicCloud[11]. The distance of about ly from the Earth and a very big burst luminositymade the comet-neutron star explanation unacceptable. However, the supporters of theGalactic origin came with two explanations undermining the significance of this detection– the position coincidence between GRB and Supernova could be purely accidental or thedetected outburst was not a GRB, but an other class of phenomena.Following years brought two strong evidences of short distances to GRBs. The firstwas a discovery of − keV absorption lines in several dozen detections performedin years 1979-81. These lines were interpreted as a cyclotron frequencies of electrons inmagnetic field of neutron stars. The calculated flux density of a magnetic field producingrequired effect was about · gauss – typical for Galactic neutron stars and far toosmall for extra-galactic sources. In year 1981 Bradley Schaefer from MIT discovered onan archival photography from 1928 an optical flash in the place of a burst from 1978.In 1984 two more such coincidences were found. That lead to an idea, that GRBs area repetitive phenomenon with an average outburst frequency 1.1 per year, with source9ifetime of more than 50 years. Thus GRBs could not be destructive for their sources,which means their energy is too low for extragalactic origin.In 80ties there were just two arguments supporting extragalactic origin of GRBs, easilyundermined and thus of much less significance. The first was an isotropic distribution ofbursts on the celestial sphere, and thus no correlation with the Milky Way. However, suchcorrelation would also not be observed if bursts were produced in the Galactic halo. Thesecond argument was a shortage of weak GRBs.If one considers the burst intensity L ∝ R (where R is the distance towards the burst)and the number of GRBs N ∝ R , assuming their constant density in the Universe, thenwe expect N ∝ L − , which is a so-called power-law. Balloon experiments performed in1982 showed a shortage of weak bursts and thus a significant deviation to this powerlaw. If the phenomena took place only in much younger Universe, thus on cosmologicaldistances, such a deviation would be expected. Per contra this could also be explainedby a not high enough sensitivity of detectors, which could not detect longer bursts withsmaller amplitude, expected from the larger distances.Nonetheless the idea of Galaxy origin was also shattered. The spectral lines occurred tobe an effect of deconvolution of the signal. The idea of GRBs as a repetitive phenomenonalso lost its impact, for there were still no detections of optical flashes accompanyinggamma ray bursts. Interestingly, GRB from the Large Magellanic Cloud truly turned outto be a different kind of phenomenon – a so called Soft Gamma Repeater (SGR) .This uncertain situation held until 1991, when the satellite Compton Gamma RayObservatory (CGRO) equipped with the BATSE detector[12] was launched. BATSE wasan instrument dedicated to GRB observation and had much better sensitivity and angularresolution then any other experiment up to that day. It detected about 1 burst per day with position uncertainty of ◦ − ◦ . This detection rate together with high enoughsensitivity lead to a definite confirmation of the deviation from the power-law. Thatfact ultimately shifted weight towards the cosmological origin of gamma ray bursts. Alsoconfirmed, with high statistics, was the isotropic distribution of bursts on the sky (see fig.2.2).The definite confirmation of the cosmological distance hypothesis however, was yet tocome. The main role in this event played the Dutch-Italian satellite BeppoSAX[13], de-signed for spectral and variability measurements of X-Ray sources, as well as all-sky mon-itoring for transient phenomena. The monitoring was performed by a so called PhoswitchDetector System (PDS) a crystal scintillator with nearly whole-sky coverage, but poorangular resolution. PDS allowed the first detection of an X-Ray afterglow of a GRB on th of January 1997. But it was the Wide Field Camera instrument with ◦ × ◦ fieldof view and angular resolution of about 3 arcminutes that made a first measurement ofdistance to GRB possible. On th of August 1997 it precisely measured a position of adetected burst and sent it to the on-ground stations, which detected a radio afterglowand, 3 days later, an optical afterglow (Keck telescopes [14]). Precise measurement of theoptical afterglow resulted in determination of the redshift z = 0 . . This meant that thedistance to the burst was about × ly – definitely outside the Galaxy.A very bright afterglow to the burst detected by BeppoSAX was observed on 28thof April 1998. The recorded spectrum and lightcurve were similar to those generated bysupernovae. This lead to a hypothesis that GRBs longer than 10 s are produced during a Soft Gamma Repeaters are sources of repeating gamma radiation outbursts. Such outbursts arebelieved to be produced on magnetastars – neutron stars with huge magnetic fields. BATSE detected nearly 3000 GRBs during its life period. It was discarded due to political reasonsin year 2000. igure 2.2: Galactic coordinates of all GRBs detected by BATSE with their fluence indicatedby colour. collapse of a very massive (more than 20 solar masses) stars – hypernovae.It seems that the final proof was given by GRB030329 , with a rapidly darkeningafterglow, which started to brighten 24 hours after the explosion, having a maximum 1.6days after the trigger. It’s further decay was much more gentle than usual and after 9days it became similar in shape and spectrum to that of supernovae. The interpretationis that the rapidly darkening part of the afterglow was generated by the GRB, while thefollowing raise and gentle decay was created by the hypernova itself.If we assume that the hypernova is responsible for all bursts longer than 10s, whyhaven’t it been seen in all the other GRBs? The answer is that the GRB030329 was closeenough to the Earth (2 Gly), that the hypernova was bright enough to be seen, while forother, more distant bursts, we would only see the GRB and its afterglow.For a long time optical counterparts were believed to be delayed temporally withrespect to the gamma ray radiation. However, development in observational techniquesas well as the increasing number of detections pushed the start of the optical emissionfurther and further towards very beginning of the phenomenon. Important dates arethe nd of January 1999, when small robotic telescope ROTSE discovered an afterglow,with peak brightness m , 20 seconds after the trigger sent from BATSE. For years itremained the fastest observation and the most luminous optical counterpart observed, upto th of March 2008. On this day “Pi of the Sky” experiment was the first one to detectGRB030819B in optical band simultaneously to the start of emission in gamma rays. Theother 2 optical experiments to confirm this discovery were TORTORA[15] and Raptor[16].The burst is so far the brightest observed by man in optical band (reaching . m ), and theobservation confirmed that the optical radiation is produced at the same stage of GRB All gamma-ray bursts are named after the date of their detection. Each pair of numbers after the“GRB” corresponds respectively to year, month and day. If more than one GRB was detected on a singledate, the distinction is made with consecutive letters of alphabet added at the end of the name. Magnitudo m stands for an apparent brightness of an astronomical object in relation to an object ofreference: m = − . II r − m r , where m r is the magnitudo of the object of reference, I r its observedradiation flux and I is the observed radiation flux of the discussed object. This unit has been normalizedaccording to the Ptolemaeus scale, where the most bright stars seen with the naked eye had brightnessof m and the darkest of m . The Sun has a brightness of − m , and Vega m . So far a few thousand GRBs have been detected. Their distribution on the celestial sphereis isotropic and so far no convincing evidence for grouping in any subclass of bursts hasbeen found. They occur outside the Galaxy, mostly on cosmological distances – the mostdistant one occurred at z = 8 . – 13 Gly from the Earth, becoming the second mostdistant object observed by men. The closest GRB had z = 0 . , meaning it was about2 Gly from us.Hardly any regularity is observed in GRBs lightcurves. All exhibit a certain raise anddecay, but these properties can be seen as well in a single or multiple peaks, which varyin shape and duration. The overall duration of this phenomena ranges from millisecondsto hundreds of seconds. Some of the bursts have a precursor – a much smaller burstpreceding the main explosion. So far precursor has only been seen in gamma rays. Oneproperty is common among all GRBs – energy of each one is of the order of − ergs.BATSE observations[23] showed that bursts can be divided into two subgroups – longbursts, with duration of more than 2 seconds and average duration of about 30 seconds,and short bursts with average duration of 0.3 seconds (fig. 2.3). Additionally, short burststend to have a much harder spectra . Long bursts seem to be more luminous and occurat larger distances.There have been numerous attempts to provide other divisions, but so far none has Figure 2.3:
Durations of gamma ray bursts detected by BATSE. T stands for the time inwhich of the burst energy was emitted. Visible is a definite division into two groups – with T < s and T > s. Hardness ratio of a GRB is usually defined as a ratio of two fluences in different energy bands,integrated in time over the duration of the burst. igure 2.4: Fireball shock model of GRB radiation generation (image taken from [24]). Theoutflow from the central engine is accelerated in phase 1 to 2, then in 3 it becomes transparent.The first stage of electromagnetic emission takes places when internal shocks (CI) occur in phase4, then the outflow propagates into external medium, causing forward (6,7) and reverse shocks(5) producing electromagnetic radiation of lower energies. The jet becomes subrelativistic inphase 8. been nearly as convincing as the described one. Perhaps the number of bursts observedwith modern instruments is still too low – in year 2010 SWIFT detected its th GRB,the number still much lower than about 2700 bursts detected by BATSE.
In general, all GRBs are believed to be created in the following way. At the beginning,a so-called central engine starts to emit vast amounts of matter (or energy). While thenature of the central engine is still a subject of speculations and may be different forvarious bursts, the generation of the radiation visible on the Earth is currently describedby a fireball shock model (modified according to recent knowledge)[25], as seen on fig.2.4.The outflow is emitted in relativistic packets (“shells”) with different Lorentz factors Γ (cid:29) . The difference in velocities of the packets cause internal shocks, which in turnlead to emission of radiation in so-called “prompt emission phase”. It has been a recentdiscovery, that not only gamma rays are produced at this stage, but also optical radiation,although via different physical processes[26].The second most important stage of emission is the collision of the outflow with inter-stellar medium (ISM) or wind from the GRB progenitor. The forward shock is responsiblemainly for generation of the X-Ray radiation[27], while the reflected wave propagating The composition of the outflow is still unknown, although weight seems to shift from baryon-dominated to magnetic-dominated recently. to ergs.Data from the SWIFT satellite exhibited another property of GRBs – so-called X-rayflares, which are “bumps” of X-ray radiation appearing in the lightcurve after the start ofgamma emission. Those flares are seen in a significant fraction of GRBs and in most casesare 10 times less fluent then the GRB emission. The important fact that arises from thediscovery of flares is that they are most probably due to the late activity of the centralengine[29].As mentioned before, long GRBs exhibit some properties of supernovae, thus a hy-pothesis of a collapse of a very heavy supernova – hypernova is the leading one. Ingeneral, after the physical mechanism of a large object collapsing, this family of modelsis called “collapsar” models[30]. The nuclei of hypernova collapses forming a black hole surrounded by an accreting disk of matter. Due to the rotation, a stream of relativis-tic barions is ejected which collides with collapsing shell of the star, and the mechanismpredicted by fireball shock model follows.While attributes of supernovae have been seen for long bursts, no such correlationhas been detected for short bursts, thus a central engine driving this type of phenomenais much more of a speculation. However, experimental data such as the hardness ofshort bursts spectra and their correlation with certain regions of the universe suggest a“merger” family of models, in which two compact objects, such as two neutron stars or aneutron star and a black hole revolving around each other, collide and merge, forming acentral black hole with a disk of accreting matter[32]. Following mechanism of producingradiation is in general driven by the fireball shock model.Some of the short bursts observed are probably stronger outbursts of the so-called SoftGamma Repeaters (SGR), which in most cases emit a short pulse of gamma radiationfrom time to time, but in most cases much weaker than in case of GRBs. However, therehave been certain observations of strong pulses coming from SGRs, which, in case ofunknown repeaters may be detected as short GRBs.In most cases, both the collapsar and the merger type models include the spinning ofthe central engine, which is probably the main reason behind the outflow being emittedin jets. When the jest axis is directed at Earth, the electromagnetic emission from ultra-relativistic jets, reaching Lorentz factor of up to Γ (cid:39) , can be seen from cosmologicaldistances, making GRBs one of the most distant objects detected by man. This wouldnot be possible if the emission was isotropic, like we assume in case of supernovae.The nature of jets is not yet well understood. The high brightness of “the nakedeye burst” GRB080319B – the optically most luminous GRB detected so far – reaching . m [26] is attributed to the two-component jet structure. The inner jet with higherLorentz factor and the opening angle of only . ◦ was directed at Earth resulting inhigher apparent brightness of the optical emission. it is possible that more GRBs exhibittwo jet structure, however, in most cases we only see the outer, wider and less energeticjet, thus the higher magnitudo of other GRBs. Hypernova is now generally a name for supernova type Ib and Ic. There is an evidence showing thatsome of the type Ic supernovae are responsible for GRBs, although it is believed that both Ib and Ic canform a GRB, depending on the geometry of the explosion. In most models, the black hole is formed in case of a GRB, while in a case of a common supernovaea neutron star is formed. However, there are some exotic GRB models which introduce a multi-stagecollapse, first to a neutron star, than to a black hole or even a strange star[31]. .2 Other phenomena of interest Novae are perhaps the most famous representatives of so-called cataclysmic variables – abinary systems consisting of a white dwarf and its companion, a donor star, from whichmass is transferred to the dwarf, forming an accretion disk. They exhibit increases ofbrightness from 6 to 18 magnitudo and thus were among the first variable stars to bediscovered. Emission takes place in optical and ultraviolet, as well as X-ray and gammarays[33].The increase has the two-step structure. The first step during which brightness isincreased by m − m takes approximately 1 to 2 days and is often followed by a fewhours to few days plateau and a sudden brightening by about m . After reaching its peakbrightness, star begins to darken – the time of darkening by m indicates the type of nova,from “very fast” (less than 2 days) to very slow (between 151 and 250 days).Spectroscopic observations show that during the explosion the star loses about . of its mass, releasing about − J. Such energies are very likely to be released onlyin nuclear fusion, which is believed to be the mechanism of nova outburst. The fusion isignited by hydrogen from the accretion disk flowing onto the surface of the white dwarf(see fig. 2.5). Due to high gravity, hydrogen is compressed and heated to very hightemperatures. When it reaches about · K, the thermonuclear reaction occurs, andthe hydrogen is burned in the CNO cycle. Remaining gas on the stars surface is blownaway by the explosion, producing light.All novae are believed to be a repetitive phenomenon, for only about . solar massis ejected in the explosion, which is relatively small compared to the white dwarf mass.Then the accretion and flowing of the hydrogen on the stars surface probably repeats .However, the time between outbursts is probably too large to be observed for most classicalnovae. There is a separate class of cataclysmic variables called “recurrent novae”, explodingmuch more frequently (10 to 80 years), but not as luminous, and brightening by m to Figure 2.5:
An artist view of a white dwarf during a nova explosion (by David. A. Hardy). If the white dwarf exceeds its Chandrasekhar limit due to long accretion, the interior may becomea subject to a carbon fusion, destroying the star in a type Ia supernova explosion. m only due to the fact, that the amount of matter participating in the nuclear fusion isabout 10 times smaller than in classical novae. Up to year 1930 there was no distinction between novae and supernovae, for both typesof explosions appeared as new stars to observes and no methods to compare amountsof energy released were available. In that year, a new class of stars – supernovae – wasproposed by Baade and Zwicky [34]. They analysed a nova-like event S Andromedae,which occurred in 1885 in the Andromeda galaxy. The distance to the galaxy, measuredin 1917 by Ritchey, turned out to be much bigger than previously suspected and thusS Andromedae had to release a much larger energy than typical nova.Supernovae are extremely bright explosions, often outshining the entire galaxy theyreside in. They are catastrophic in their nature, completely destroying the star, but re-leasing a shock wave into surrounding interstellar medium, creating a so-called supernovaremnant (SNR). The phenomenon can be triggered by two different mechanisms. The firsthappens for massive stars, which burn out nuclear fuel in their core. They can no longersustain a gravitational pressure and collapse into a neutron star or a black hole. Vastamounts of gravitational potential energy are released in this process, creating a shockwave and heating outer layers of the star. The second scenario is characteristic for whitedwarfs – if the amount of mass accreted from a companion star is large enough for thedwarf to come close to its Chandrasekhar limit, a carbon fusion begins, which becomes aso-called runaway fusion and drives the supernova explosion.At the first glance, supernovae seem likely to be categorized by their triggering mecha-nisms – ignition or collapse of the core. However, they were categorized mainly accordingto their spectroscopic properties (and historical reasons). Type I (consisting mainly ofsubtypes: a, b, c, distinguished by other spectroscopic features) lacks a Balmer series(line of hydrogen) in their spectrum, in contradiction to the type II (consisting mainly ofsubtypes P and L, distinguished by features of a lightcurve). Type Ia is a white dwarfreaching its Chandrasekhar limit and igniting its core, while types Ib,c and type II are
Figure 2.6:
On the left: artist’s impression on supernova SN 1993J, which later occurred to bemixed type II/Ib ( (c) by ESA ); on the right: Hubble telescope reconstruction of supernova SN1994D, type Ia ( (c) by Hubble ). ergs . Most of the supernova explosions are triggered by a collapseof the core (about ). The frequency of supernovae events in the galaxy of size of theMilky Way is estimated to be 1 per 50 years, although the last supernova in Galaxy wasobserved in year 1604. However, about 150 such explosions are yearly discovered in othergalaxies. Flare stars (also called UV Ceti variables) were discovered much later than novae orsupernovae, although they exhibit as violent, or even more violent behaviour than theformer. That is probably due to the fact, that most of them appear as a brightening ofalready visible star, not a new star, and that they are rather a rare phenomenon. Mostof the flare stars observed are closer then 50 ly from Earth – that is due to the fact, thatflaring occurs in red dwarfs, which are very dim . However, some were seen as far as1000 ly away.Their brightness can increase by two magnitudo just in seconds or minutes and dropto normal brightness within minutes to hours. The increase takes place in nearly wholespectrum, ranging from radio to X-ray, and the event is highly unpredictable. The flaringtakes place due to the similar mechanism as solar flares – magnetic reconnection. Fluxpatterns of plasma in the star’s atmosphere are typically aligned with the star’s magneticfield. From time to time the field rearranges, dropping to a lower energy state. Realign-ment of the plasma causes it to collide with itself, heating up even from 3000 K to 10000K, and thus resulting in a star brightening.Such flares are much brighter than solar flares for two connected reasons. First of all,red dwarf mass is only few tenths of the Sun’s mass, thus the material can be ejectedmuch more energetically due to lower gravity. Additionally, red dwarfs are much dimmer Figure 2.7:
An artist view of flare star EV Lacertae during an explosion ( (c) NASA ). Similar as in the Gamma Ray Bursts (sec. 2.1), but released in a much longer period of time. Recent research indicates that flaring may also occur in brown dwarfs. of their surface, is much morenoticeable.The flare star closest to Earth is Proxima Centauri, however, the first to be discov-ered were V1396 Cygni and AT Microscopii in 1924. The most famous one is UV Cetidiscovered in 1948, from which the whole group took name (UV Ceti variables).
Supernovae, novae and flare stars are of special interest to the “Pi of the Sky” project,for they may often appear as a completely new objects on the sky and thus are notsubject to regular observation but rather to detection. However, an important amount ofastronomical knowledge is derived from study of other variable stars. This study can besuccessfully performed by an experiment such as ours.Probably most stars are variable in some way – even Sun has a 11 year solar cycle,during which its energy output varies by about . . There are currently two maincategories of variables, distinguished by the mechanism causing the fluctuations in theirlight flux.Intrinsic variables change their actual luminosity due to some phenomenon happeningin the star. Among this category are flare stars, novae and supernovae mentioned before,belonging to the subtype of eruptive or cataclysmic (or even catastrophic) variables. Thereare also stars which, for example, change their radius due to varying temperature anddensity in their interiors, such as pulsating variables.In case of extrinsic variables external properties cause the variability in the amountof light that can reach an observer. The first subgroup are eclipsing variables, normallybinary systems where one member covers another during their orbiting. The second group– rotating variables – change their apparent brightness due to the different parts of thestar being visible to the observer and releasing different amounts of light. There was a number of objects observed by optical astronomers counted as irregularvariable stars – objects that changed their brightness on a wide range of periods butwithout any apparent pattern. With development in radio astronomy, especially withincreasing the angular resolution of radio telescopes and discovery of quasars in late 1950sit became apparent, that those stars were in fact different objects – blazars.It is believed, that in the centre of most or even all massive galaxies supermassive blackholes exist. Especially in case of giant elliptical galaxies, they form active galactic nuclei(AGNs), which are powered by matter falling onto the black hole. Falling matter formsan accretion disk, which emits vast amounts of energy in form of particles. However, theproperty that decides that an AGN is a blazar are highly relativistic jets (0.95-0.99 c)emitted along the angular momentum of the accretion or the spin of the black hole di-rection. Plasma in the jets interacts with itself and magnetic field emitting radiation.The spectrum ranges from radio to gamma-rays. For some cases TeV photons were alsoobserved. If a jet is pointing in the direction of Earth, the object is called a blazar.The important distinction between blazars and other AGNs comes from the velocitiesof emitted jets. Effects of special relativity, so-called Lorentz boost causes a jet to be afew to few hundred times brighter for an observer than in its rest frame. Additionally,the time scale of the jet variability is much shorter when seen from Earth.18 igure 2.8:
The first discovered quasar and blazar – 3C 273. Image by Chandra X-Rayobservatory.
Blazars are the most luminous persistent sources in the Universe known to man, highlypolarized ( − ) in radio and optical frequencies, showing an irregular variability ontimescales from days to years. Their brightness and variability makes some of them objectswell suited for observation by small, high time-resolution systems such as “Pi of the Sky”. There is a number of visitors to the night’s sky that are not of stellar or quasi-stellarnature. Comets and meteors were observed since the birth of astronomy. Nowadays weknow, that there are many more celestial bodies much closer to Earth than stars, such asasteroids and planetoids. Careful analysis of their lightcurves can reveal their shape andother properties. However, a pure detection of such an object and deriving its trajectorymay prove more important in the future than analysis of its properties – an answer to anunending threat towards Earth.There is, however, a group of objects that became a subject for studies very recently.These are artificial satellites of Earth, both machines still functioning and working forman and cosmic debris – mostly leftovers of past, when man could not imagine, thatpolluting our planet’s orbit may turn out to be a big problem in the future.The study of satellites is important due to an increasing demand for placing new func-tioning objects on the orbit. Finding a free place for them is problematic, for trajectoriesof most cosmic debris (as well as some spy satellites) are unknown. Additionally, debrisis a large threat for objects lifted from the planet’s surface to the orbit.Fortunately, most if not all artificial satellites reflect Sun’s light and thus are visiblefrom the Earth. Careful studies may result in catalogues of cosmic debris with theirparameters. However, this job may only be performed by joint efforts of experimentswith a wide field of view and thus large detection capabilities. One of such experimentsis “Pi of the Sky”. 19 hapter 3The “Pi of the Sky” experiment
Progress in astronomical instrumentation constantly increases our understanding of theUniverse and reveals new, often unsuspected, mysteries. Unraveling some of those ispossible by gathering and careful study of already existing data. However, some newdiscoveries call for a completely new scientific approach. This was the case with gammaray bursts – the main inspiration for the “Pi of the Sky” project.
The instruments that accidentally discovered GRBs had a non-scientific purpose – mon-itoring of large parts of the sky – a purpose strange to astronomers at these times (sec.2.1). This idea became more popular in this field of science along with the understandingof GRB phenomenon. Following 30 years of studies lead to construction of the satelliteinstrument BATSE, then SWIFT and Fermi – satellites capable of constantly observingvery large parts of the sky with high time resolution and thus detecting hundreds of burstsin few years time.Up to now a few thousands of GRBs were detected by dedicated gamma-ray satelliteexperiments[35]. Only a fraction of these bursts were also seen in other wavelengths, andthere were only a few simultaneous observations of optical and gamma-ray emissions. Justone of those measurements had good enough optical time resolution to claim that the firststage of optical and gamma emissions took place at the same time.The main reason for such a difference in the number of gamma-ray and optical de-tections of GRBs is the difference in methodology of experiments. Satellite experimentsare aimed at the detection of GRBs. They constantly monitor a large fraction of the skyand perform a real-time analysis of data, to detect a signal. The idea that stands behindthe vast majority of optical cameras and telescopes is opposite – they are dedicated todetailed observation of specific astronomical sources, thus having very high sensibility(a big “range”), a very small field of view and nearly no automatic detection capabilities.There are also systems dedicated to detection and analysis of variable sources, which havea bigger field of view, such as ASAS[36]. However, unlike satellite experiments, they donot monitor a large fraction of the sky, but perform a survey, taking a photo of a specificpart of the sky every few hours or days. This limits their ability to measure object’svariability on a scale of minutes or shorter and nearly rules out detection of short burstssuch as GRBs.The need for optical monitoring of a large part of the sky is clear, for studying theemission in visible wavelengths from its very beginning can unravel crucial informationabout this phenomenon (sec. 3.5.4). Moreover, if there are other optical transients in the20ky not yet discovered, happening on timescales of seconds, is there a better way to findthem?
The “Pi of the Sky” project was designed to answer the need for optical monitoring of thelarge part of the sky. The name of the project is derived from the concept of watching π steradians of the celestial sphere, almost the whole part available to astronomical obser-vations from a single point on the Earth . Currently we assume that the full system isgoing to cover 1.5 steradians of the sky. This area will hopefully increase in the future.Additionally, a high time resolution (in optical terms) and self-triggering capabilities aregoing to be the main features of the experiment.The field of view of 1.5 steradians is obtained by using 12 cameras, each covering ◦ × ◦ of the sky, using Canon EF lenses with f = 85 mm and f /d = 1 . . Camerasare of a unique construction and design prepared by “Pi of the Sky” project members.Exception is a commercial CCD sensor, with roughly × pixels and a pixelsize × µ m (corresponding to an angular size of 36 arcseconds). This setup, with10 s exposures and 2 s readout time gives estimated range of . m on a single frameand − m on 20 stacked frames. Parameters such as cooling, readout gain, etc. canbe controlled by USB 2.0 or Ethernet interfaces. Cooling is performed with a two-stagePeltier stack and allows reaching K below ambient temperature.Additionally, a custom shutter has been developed. 10 s exposures with 2 s readoutresult in 2000 to 3000 frames per night, 200000 to 300000 exposures after 4 months ofoperation, the maximal lifetime of best commercial shutters. “Pi of the Sky” shutter,sustaining cycles is enough for a few years of uninterrupted operation. Figure 3.1: “Pi of the Sky” final system mount with four cameras installed at the INTA testsite near Huelva, Spain. The whole sky visible from a flat surface would be π steradians, but a large part of the sky closeto the horizon is filled with intense ambient light, which is the main factor making this area not suitablefor scientific observations. Term “pixel” in this work is used in two meanings: pixel as a (in this case) square area on the frameor CCD sensor or as a distance equal to the side of such square.
100 km ←−−−−−→
Figure 3.2:
A scheme of the final “Pi of the Sky” system – two sites, each with 12 camerasinstalled on 3 mechanical mounts. Sites are distant by about 100 km, to allow filtering outnear-earth objects flashes with the parallax measurement.
Four cameras are installed on each of three paralactic mounts (fig. 3.1). Two modesof observation are possible – so called “deep” and “wide”. In the “deep” mode axes of allcameras on a single mount are parallel, so they observe the same field of view. Propersumming of the frames in such a configuration increases the range of the instrument . Incase of detecting an interesting source (i.e. GRB optical counterpart), the flash recognitionalgorithm can automatically send a command to the mount to switch cameras into the“deep” mode. In the “wide” mode, which is the standard mode of operation resulting inthe field of view of 1.5 steradians, cameras on a single mount cover neighbouring parts ofthe sky.Two sites are being prepared for the final “Pi of the Sky” system, aiming to observethe same part of the sky from two locations (see fig. 3.2). With distance of more than30 km (probably about 100 to 200 km) such configuration can be efficiently used toeliminate flashes caused by artificial sources – mostly satellites and cosmic debris reflectingthe light of the Sun using parallax. Depending on their linear and rotational speed,construction and distance to Earth, part of such reflexes are similar in characteristics (bothshape and duration) to optical flashes of astrophysical origin (sec. 3.5.2) and parallax Site separation [km] Orbit [km]
30 92000100 293000200 579000
Table 3.1:
Three exemplary sites separation distances and corresponding maximal orbits ofsatellites that can be rejected by parallax, assuming position difference threshold for automaticevent rejection equal to twice the pixel size (72 arcseconds). The simplest way to increase the apparatus range using multiple images of the same field is to makean average of such images. This way the signal to noise ratio is increased, for the random noise is notadditive, while the light from interesting sources is. Of course, more sophisticated algorithms of suchaveraging can be used, which take into account slightly different position and orientation of a source ondifferent images. igure 3.3: Locations of the existing and a considered sites of the final “Pi of the Sky” system.The existing site A lies near Huelva, while the second, B, near Malaga in Spain. is the only way of undoubtful rejection. Table 3.1 shows maximal orbit radius of anartificial satellite which can be rejected with parallax for specific site separation distances.Separation distance of 30 km should allow for rejection of most near-earth objects such ascommunication satellites, however 100 km separation would be needed to reject very highorbit satellites, such as Vela, which should be sufficient for recognition of all near-Earthobjects in the “Pi of the Sky” range. The 200 km separation distance extends the parallaxusability beyond Moon’s orbit.The whole setup has been designed to work fully autonomously, requiring a directhuman assistance as rarely as possible – hopefully regular maintenance only once or twicea year, or in case of serious emergency. All elements – cameras, mounts, domes, etc. canbe controlled via Internet, thus allowing remote managing of minor disturbances. Theflash detection is done on-line through a multilevel trigger, nearly instantly generatinginformation about detected interesting events. Simultaneously, selected parts of data aregathered for further, offline analysis. The system starts and stops automatically, changingthe modes of observation according to specified rules and scripts. All collected data area subject to data reduction algorithms and the final results are stored in the dedicateddatabase.The installation of the first mount of the final systems with four cameras took place inOctober 2010. The observation site is located near Huelva in Spain in the site of INTA[37].For the second site, CSIC[38] research station near Malaga is considered. The distancebetween the two sites is about 240 km, see fig. 3.3. Nearly all of the algorithms developedfor the “Pi of the Sky” project have been successfully tested with the prototype systemworking for 6 years in Chile. 23 .3 The prototype
The prototype apparatus of the “Pi of the Sky” project (fig. 3.4) was installed in theLas Campanas Observatory[39] on the Atacama desert in Chile in the site of the CarnegieInstitution of Washington, in June 2004. Thanks to the hospitality of the ASAS group[36]it operated in the ASAS dome till December 2009. This site was chosen mainly becauseof conditions, which are among the best on the Earth for observational astronomy – verylow humidity, an unclouded sky most of the year and low sky background – due to theclimate of this place and its high elevation (2380 m above sea). Recently the system wasmoved to the San Pedro de Atacama observatory in Chile, about 1300 km north of LCO,and started gathering data there.Apparatus consists of a paralactic mount with two cameras having the field of view,number of pixels and optics similar to those being installed in the full system. Cameraswatch the same part of the sky and work in a coincidence mode. The distance betweenthem is far to small to make use of the parallax, but this mode allows rejecting falseflashes caused by noise fluctuations and cosmic rays – most significant background sourcefor flash detection algorithms.Since the overall field of view of the prototype is ◦ × ◦ , the standard mode ofoperation is to follow the centre of FoV of the SWIFT. The mode is changed into followingINTEGRAL satellite FoV or observing selected sources, when the SWIFT FoV is belowhorizon. At the beginning and by the end of the night a scan of the whole sky is performed,taking 3 images of each field. Contrary to the normal mode, which is dedicated to detectingoptical transients on scale of seconds to minutes, scan data is suitable for off-line analysisof sources variable on timescales of hours to months. All the modes of operation can beinterrupted by a GCN alert – information about a GRB in progress – which cause themount to point to the burst coordinates and start follow-up observations. Afterwards,the system reverts to the previous task.A very high level of autonomy and reliability has been reached through years of test-ing the prototype. Currently system has very sophisticated self-monitoring capabilities,automatically detecting and solving problems and sending an appropriate e-mail or smsto people involved. No human interference has been required during normal night-shifts Figure 3.4:
The “Pi of the Sky” prototype at the in Las Campanas Observatory on the Atacamadesert in Chile. m − m range, depending on the conditions) and gives a real time informationabout detected transients (sec. 3.4.1). During the day-time, data is being analyzed indetail and catalogued to a database, on which off-line algorithms searching for variabilityon longer time-scales are being run (sec. 3.4.3). The data analysis software is explainedin more detail below. A common approach in astronomy is to gather data during observation time and performa sophisticated data-reduction later. This approach, although appropriate for dedicatedstudies of some objects of special interest, is not suited for a project like “Pi of the Sky”.Single “Pi of the Sky” camera gathers 16-24 GB of data per night. Full system is goingto produce 512-768 GB per night, amount impossible to store for a longer time, notmentioning a detailed analysis. However, only a very small fraction of this is interestingfor us and has to be saved. The question is, how to select this fraction?In its approach to flash recognition, the “Pi of the Sky” project is similar to largeparticle physics experiments, which analyse huge amounts of data in search for very fewinteresting events. In such a case, only interesting events including signatures of infor-mation required for detailed analysis are kept, the rest being discarded as “background”.This task is performed by so-called triggers – multilevel selection algorithms that receivehuge number of events and reduce it by applying cuts on consecutive levels. They startwith very fast and simple rules, eliminating the most obvious and numerous backgroundand finish with sophisticated methods dedicated to analysis of the few events remaining.Same logic is utilized by our project, the only difference is the nature of the data stream– series of images of the sky.At the input to the selection algorithm, the frame, after dark subtraction, is modifiedby a kernel-based transform, similar to a discrete Laplace transform in electromagnetism.In ideal case this would make all real light sources pixel-like, remaining dark pixels abackground. In the real case it sharpens the image and vastly reduces effects caused bybackground light gradients, mainly due to the Moon and clouds. Following steps of thealgorithm are a simple recognition of pixels above estimated background level – T n (afterremoving saturated pixels), and veto on constant objects T v , which eliminates mainlystars – pixels above background that were already present on previous frames. These twosimple cuts gives a list of new candidates for flashes on the frame, still most of them notof astrophysical origin.Following is a number of cuts that eliminate temporal errors of a CCD (like “hot”and “black” pixels). Then a list of pixels, which may be attributed to a single flash istransformed into a list of objects. The last of simple, pixel-based cuts applied removesstacks of close objects, mostly caused by moonlight illuminating dense, moving clouds.A border line between the very fast, “level 1” trigger and the more sophisticated “level2” trigger is the “coincidence”. Lists of objects remaining after the cuts from the twocameras (observing the same field of view) are compared and only those flashes that arepresent on both are left for further analysis. In case of the prototype, this cut eliminatesmostly events caused by cosmic radiation and random noise fluctuations. In the final25 igure 3.5: Number of events ¯B considered in the on-line flash detection algorithm after con-secutive levels of the trigger NC[40]. Different markers represent different nights. The analysiswas performed in 2004, at the beginning of the prototype operation. Later studies allowed toreduce the number of final events by an order of magnitude. system, due to parallax, coincidence will also remove flashes caused by near-earth objects,like satellites and planes. That will vastly reduce or even discard the need for the wholelevel 2 algorithm .Level 2 cuts are designed to remove satellites, planes and fluctuating stars. Thisis mainly performed by discarding events whose position coincides with position of asatellite or a star found in dedicated catalogues. An event is rejected also if flashesobserved on consecutive frames can be grouped into a track or have an elongated shape.The reconstruction of tracks require a number of frames after detecting the flash to bestored an analysed, which causes the time between the flash and a final confirmation ofthe detection to be quite long. Fortunately, as mentioned before, level 2 cuts should notbe needed in the final system.The number of pixels considered on the input of the flash recognition algorithm in asingle frame mode exceeds per night for a pair of cameras (fig. 3.5). Simple, pixel-based cuts reduce number of flash candidates to the order of − and coincidencereduces it to the order of slightly more than events per night. In the time of writingthis work, the number of final events on which an “eye analysis” has to be performed isapproximately few to over a dozen per night.However, more than a dozen events per night for a pair of cameras means even asmuch as 300 hundred events for human examination in the full system. That’s why athird level algorithm is being worked on, which includes even more sophisticated methodssuch as Hugh transform for finding tracks of meteor and plane like events[41].For each of the final events some basic astronomical data has to be calculated, likeits brightness and position in celestial coordinates. These tasks are performed on-line by An alternative to coincidence of flashes on two cameras is a requirement for flash to appear on twoconsecutive frames on a single camera. Although this approach requires an optical transient to be longerin time, it was used when there were problems with one camera in the prototype system. Currently, priorto the installation of the new apparatus of the final system in the second site, it is used in flash detectionalgorithm in INTA. igure 3.6: An aperture used in the fast photometry algorithm[42]. Pixels in the yellow circleare summed as signal, while a median of pixels in the black ring (radius R_in to R_out) is usedas a background estimator. special photometric and astrometric algorithms which are much faster but less accuratethan off-line objects cataloguing procedures.
Astrometry is a term used to describe a way of finding the celestial coordinates of anastronomical object from its image. The first step is determining the position of theobject on the frame, which is not a straightforward procedure. Each optical imagingsystem is characterized by a point spread function (PSF) – a response of the detector(in this case lenses + CCD sensor) to a point source of light . Thus a point-source thatnormally would be contained in a single pixel is spread over several pixels, mainly due todiffraction and optical aberrations. That fact makes the calculation of the source’s positionoften a complicated task and a number of algorithms analysing data stored in multiplepixels were invented. Fortunately, it also allows precision of coordinates determination tobe a small fraction of a pixel, whereas for a single pixel the uncertainty is √ .There is a wide variety of astrometric methods. Starting with simple or more sophis-ticated pixel-based, where the centre of an object is calculated from a set of chosen pixelsas a, for example, centre of mass. This type of procedure is widely used in “Pi of the Sky”.Often slower, but giving more precise results is fitting of a known point spread functionto the object – one of the main topics of this thesis (chap. 6).Astrometry procedures often include also measuring of the object’s brightness, whichis, in general, the aim of a photometric algorithm. These also range from simple, pixelbased algorithms, many of them described as “aperture photometry” , up to those men-tioned above, based on the detailed PSF modelling, called “profile photometry”. Aperture Objects such as stars or most of the artificial satellites of Earth, due to their large distance and anangular resolution of 36 arcseconds are point-sources to “Pi of the Sky” The name “aperture photometry” comes from a defined configuration of pixels attributed to an objectin this algorithm, called an aperture.
Although some very rapid and big changes in brightness of known objects may be detectedby the on-line flash recognition algorithm, most variable stars are recognized through off-line algorithms. The step preceding all off-line variability searches is cataloguing of thestars, i.e. performing astrometry and photometry for each object on single or 20 stackedframes in case of normal measurements and 3 stacked frames in case of images from thesky scan. Results for each object are then put into a database, along with informationabout conditions during the measurement, temperature of the detector, total number ofobjects detected, etc. A database created in this way is then processed by algorithmssearching for novae-like events, flare-like events and those determining variability period.Novae search algorithms look for new stars in the database. For each new objectfound in the night being analysed, a number of cuts is applied, to eliminate false events,such as fluctuations coming from very bright stars, hot pixels, vibrations of the mount,etc. Remaining objects are a subject to an examination by a member of the project.This method is proper for finding novae, supernovae (although no supernovae during the“Pi of the Sky” work period was in our range), blazars, some flare stars or other variablestars with large amplitude, that are normally below “Pi of the Sky” range.Flare star recognition method is similar to the novae recognition. The difference is,that all existing objects in the database are examined in search for a sudden increasein magnitudo compared to previous measurements. Stars satisfying this criteria are asubject to consecutive quality cuts, same or similar to those in the novae search. Resultinglightcurves have to be examined by man.The last one is the algorithm determining the periods of variability for variable starsin the database. Based on this algorithm, periods of already known variables can bedetermined, as well as possible periods for stars suspected to be variable. Also calculated isthe uncertainty of estimated period, amplitude of magnitudo variations and the parameterreflecting quality of the fit – this information is a key to finding new variable stars.28 .5 Scientific results
During 2004-2009 working period, the novae detection algorithm of the “Pi of the Sky”prototype system automatically detected 3 novae stars. The nova V679 Carinae (fig.3.7) was discovered on 1st of December 2008 and was observed for many days later[43].“Pi of the Sky” released an AAVSO[44] notice requesting followup observations whichresulted in precise astrometry and spectroscopic observations by the SMARTS[45] 1.5mtelescope. The star occurred to be a classical nova of the “Fe II” type.On th of September 2007 outburst was detected in the place of 1RXS J023238.8-371812 object. Followup observations showed redshift zero and a blue spectrum, makingthis object possible SU UMa-type or a WZ Sge-type outburst[46]. On th of December2007, “Pi of the Sky” discovered brightening coincident with a . m star GSC2.3 S55U020591,which reached nearly . m . Followup spectroscopic observations with Gunma Astronomi-cal Observatory[47] 1.5-m telescope and a low resolution spectrograph, showed a very bluecontinuum and other features characteristic for a WZ Sge-type dwarf nova[48]. Severalother novae have been automatically detected by “Pi of the Sky”, but we were not thefirst to discover them.The first outburst detected by “Pi of the Sky” was a flare star CN Leo explosion on nd of April 2005[49]. It was so violent that it has been found by an on-line flash recognitionsystem. Off-line flare detection algorithm discovered two additional flare stars. On th ofNovember 2006 a GJ3331A brightened from . m to . m , reaching its peak brightnessin less than 10 minutes[50]. Following was a decay monitored for more than an hour(fig. 3.8).During the night of th July 2008 another outburst was detected, connected to aUSNO1050.00569325 star[51]. It’s characteristics, like significant proper motion ( 0".1/year)and 2MASS J-K color of 0.85 indicate that the object is an M dwarf, suggesting a flareexplosion.Additionally one blazar was observed by the “Pi of the Sky” project (fig. 3.9). It hasbeen automatically detected by the novae detection algorithm in the database, unfortu-nately during its testing period. Thus no regular observations have been performed andits unprecedented bright peak has been missed. However, it would be possibly seen by
Figure 3.7:
V679 Carinae nova explosion, automatically discovered by the “Pi of the Sky”project with the offline novae detection algorithm[43]. “HJD” stands for Heliocentric Julian Day. igure 3.8: Flare star GJ3331A automatically discovered by an off-line flares detectionalgorithm[50]. the full “Pi of the Sky” system and showed “Pi of the Sky” potential at discovering andobserving bright blazars.In addition to observations described above, a catalogue of 725 variable stars withperiods from 0.1 to 10 days, based on 2 years data (2004-2005) has been published[52]. Theclassification was performed with the procedures developed for “Pi of the Sky”, inheritingfrom light-curve shape analysis methods described in the GCVS catalogue. For 15 starswith previously unknown period, the period has been found and some stars were classifiedto a different type than claimed in the GCVS catalogue. Preliminary information onthe variable stars catalogue based on measurements from years 2006-2007 has also beenpublished[53], but the data are still in analysis.
Up to now more than 250 uncorrelated flashes have been detected by the “Pi of the Sky”on-line flash recognition algorithm. For a number of these optical transients their positionshave been very close to some identified sources or archival GRB positions, however this isnot an indicative information. Possibly most of flashes are generated by rotating artificial
Figure 3.9:
Image of the blazar 3C 454.3 taken by the “Pi of the Sky” prototype. Blazar, beingon the edge of the “Pi of the Sky” range, is the dark pixel in the top left corner of the marker. igure 3.10: A false “Pi of the Sky” trigger generated by a rarely flashing satellite (frame 0,central plot). Satellite is also visible on frames -5 (left) and 6 (right), each with 10 s exposureand 2 s readout, thus it blinks approximately once per minute. satellites reflecting the sunlight. Some satellites are observed as a track segment onthe picture, blink often enough to be attributed to a track (fig. 3.10), have movement-distorted shapes or are simply in a catalogue. However, there are some that cannot beattached to any of these groups and those will be distinguishable from an optical transientof astrophysical origin only with the parallax, in the final system.Seven among the uncorrelated flashes were visible on more than one frame (fig. 3.11).There are no known objects orbiting Earth that have small enough angular velocity to beseen in the same place (with “Pi of the Sky” accuracy) and reflect enough light to be inthe range of our experiment. Thus such events are probably of an astrophysical origin.With the full system an instant alert will be released in case of discovery of such an event,requesting follow-up, deeper range observations.
In the time of writing this thesis, nearly 800 gamma ray bursts were detected by satellitesafter the “Pi of the Sky” experiment start-up (tab. 3.2). As expected, about half of themhappened during the day-time in Las Campanas Observatory. About of GRBs were
Figure 3.11:
An uncorrelated bright optical transient visible on 2 frames in their centre, de-tected by “Pi of the Sky” 2005.05.31, 4:33:17. On the left – a frame preceding the flash, followedby frames with the flash. ll det. GRBs apparatus off north hemisph. daytimeEvents
786 138 49 401
Fract. of all below horizon clouds outside FoV inside FoVEvents
82 19 92 5
Fract. of all
Table 3.2:
Statistics of all gamma ray bursts detected from the star of the “Pi of the Sky”prototype on st of July 2004 up to th of October 2010. Only 5 out of almost 800 bursts werein the prototype’s field of view and one of these bursts was detected by “Pi of the Sky”. Up todate statistics can be found on the “Pi of the Sky” project web page: http://grb.fuw.edu.pl/pi. out of the reach of our cameras, being either fully on north hemisphere or just belowhorizon. For nearly of bursts the apparatus was off, but this includes long period ofabout 10 months prior to the installation at the new site.More than a dozen of bursts could not be seen due to dense clouds, while nearly were initially outside of the apparatus field of view – delay caused by a rotation towardsproper coordinates is in most cases too long for us to spot the transient. However, 5bursts were in the “Pi of the Sky” prototype FoV from the very beginning. For 4 of themlimits during and before the burst have been set – still a very rare data in the study ofthis phenomenon. One burst have been autonomously detected by the flash recognitionalgorithm and it proved to be a very significant discovery. The discovery of the GRB080319B optical counterpart by the “Pi of the Sky” project,as well as TORTORA and Raptor experiments was a coincidence of an idea of full skymonitoring and luck. The burst was preceded by an another, which happened roughlyhalf an hour earlier, in coordinates distant by about ◦ . After receiving the GCN alert,“Pi of the Sky” rotated to the GRB080319A and started the follow-up observations, re-maining in the same position for half an hour. Field of view of ◦ × ◦ was enough tosee, in the corner of the frame, the second burst (fig. 3.12). Thus it was a mixture of luck,two explosions having quite close celestial coordinates, and idea of full sky monitoring –without such a big field of view (which is very big in astronomical standards, even thoughit was still just a FoV of one camera), no such detection would be possible.The first, 10 s exposure of the “Pi of the Sky” prototype on which the burst is visible,coincides for the first time in history with the start of emission in gamma rays. There werefew previous simultaneous observations, but the time resolution was far to low to connectobserved optical emission with the start of emission in higher energies. Additionally, alittle bit later the burst became visible to TORTORA experiment, which, with its veryhigh time resolution registered optical peaks similar to these visible in the gamma band(fig. 3.13)[26]. However, the existence of a correlation between lightcurves is still a subjectof discussion.The burst was unusually bright in X-Ray and the most bright ever observed in opticalwavelengths, reaching peak brightness of . m . Measured distance was z = 0 . , i.e. . · ly, meaning it is the brightest object observed so far by man. It’s unusualbrightness and some other lightcurve features can be explained by emission taking placein a two component jet – the narrower, more luminous component directed at Earth,32 igure 3.12: The photograph of GRB080319B in its peak brightness, taken by the“Pi of the Sky” prototype. For the first few frames, before the detector repositioned follow-ing GCN alert, the gamma ray burst was localized in the corner of the apparatus field of view,thus a large deformation of its shape. which may be a very rare occurrence.The discovery of GRB080319B not only showed that there exists a simultaneous emis-sion in optical and gamma bands, which was a hypothesis not accepted by the wholescientific community, but additionally provided scientific data to study that emissionmechanisms. Before this observation, most models predicted that the optical emission, ifobserved from the very beginning of the burst, could be due to a low-energy tail of thehigh energy emission. However, the extrapolation of the measured gamma-ray spectrumto the optical energies (fig. 3.14) showed that in case of GBR080319B the optical fluxexceeded the possible low-energy tail by 4 orders of magnitude. Thus radiation in these
Figure 3.13:
Lightcurve of the initial emission from GRB080319B. Superimposed are pointsfrom “Pi of the Sky” and TORTORA optical experiments, and from gamma-ray detector Konuson the WIND satellite. The first point, taken by the “Pi of the Sky” prototype shows that theemission in the optical band and gamma rays started simultaneously[26]. igure 3.14: GRB080319B spectrum in gamma rays and its extrapolation to optical energiesin time intervals corresponding to the first three images of the burst by “Pi of the Sky”. Threecolourful points in the top left corner of the plot show optical flux for these intervals measuredby “Pi of the Sky”, much higher than the extrapolation. This data was used to prove, that initialoptical emission is not a low-energy tail of the gamma ray emission. two different bands has to be produced by a different physical mechanism – the most nat-ural explanation is synchrotron emission for optical and synchrotron self-Compton (SSC)mechanism for gamma rays.“Pi of the Sky” experiment was fortunate to register GRB080319B from its very be-ginning. However, on the first 3 exposures, before the detector repositioned followingGCN alert, the burst was in the very corner of the frame, and due to optical aberrationsappeared very deformed (fig. 3.12). Standard photometric procedures are optimizedfor symmetric point spread functions, thus such a deformation can introduce additionaluncertainties to measurements of brightness and position.
The quality of brightness measurement for a given object strongly depends on the signalto noise ratio. In an ideal detector with no noise, no nonuniformities, stable and lineargain, photometry uncertainty would be determined by the Poisson statistics of photonsand should be a monotonic function increasing with magnitudo. However, the real case ismuch more complicated. First of all, a saturation threshold makes objects below certainmagnitudo level too bright for the detector. Additionally, there are multiply factors34 M D CCD x [pixel]
CCD y [ p i xe l ] M D Figure 3.15:
Magnitudo uncertainty ∆ M distributions, as a function of the star magnitudo M(left) and as a function of the star position on the frame (right) for stars with < M < . causing signal to fluctuate, like gain variations, tracking instability of the detector orchanges in observation conditions.All these factors result in the dependence of magnitudo uncertainty on magnitudowhich is shown in fig. 3.15 (left). The brightness of stars below ∼ m is very poorlymeasured. This is due to gain nonlinearity at highest charge values and pixel saturation.For stars above ∼ . m a fast increase of uncertainty with magnitudo is observed, whichcontinues up to detector maximal range of around m . However, for stars between ∼ m and ∼ . m the uncertainty in brightness determination does not depend on magnitudo.This behaviour is probably induced by a combination of many factors, such as a non-Poissonian signal fluctuations (for example, fluctuations of the CCD amplifier gain) andthe photometric algorithm, which may be working best in this range.Moreover, the photometry quality strongly depends on the position on the frame (fig.3.15, right). The problem could be caused by a nonuniform illumination of the CCD,but this should be compensated to large extent by proper frame reduction (division by a“flat” frame) and photometric algorithm, which takes into account local background. Thelarge vertical gradient of the magnitudo uncertainty observed in fig. 3.15 (right) indicatesthat part of the fluctuations of the stars brightness could be caused by a vertical chargetransfer procedure. However, slight differences along the horizontal axis are also visible,which may be caused by different shapes of the PSF on different parts of the frame (opticalaxis not exactly perpendicular to the CCD surface). To find celestial coordinates of an object visible on the frame several steps must beperformed. First, a position on the CCD has to be found. Then it has to be transformedinto celestial coordinates, which requires determining pointing direction of the apparatus.The uncertainty of the latter step is rather small, for the calculation involves comparinga great number of stars with a catalogue. In the first step, similarly to the photometricalgorithm, uncertainty depends on signal to noise ratio. Thus it should increase withthe stars magnitudo. However, in case of the “Pi of the Sky” project, the dependence isnot monotonic. As seen in fig. 3.16, the lower edge of the uncertainty band drops with35 x [ p i xe l ] D x [ a r csec ond s ] D M y [ p i xe l ] D y [ a r csec ond s ] D Figure 3.16:
Horizontal (left) and vertical (right) position uncertainty vs measured magnitudoM. The dependence is clearly non monotonic, the uncertainty largest for saturated stars andthose close to the detector’s maximal range. increasing magnitudo up to m − m , and then starts to increase. Additionally, saturationthreshold does not influence astrometry as strong as photometry.The dependence of the position uncertainty on the star position on the frame does notshow an indicative vertical gradient (fig. 3.17). In both horizontal and vertical directionsuncertainties show an approximate rotational symmetry. However, the radial function isnot monotonic, showing a broad minimum at about 800 pixels from the frame centre. Itis also clear that the rotational symmetry is only approximate, as there are also changesalong the azimuthal coordinate. The observed symmetry is a strong indication that theuncertainty in astrometry (and probably in photometry too) is enlarged by dependenceof the PSF on the position on the frame. The PSF should have a significant rotational CCD x [pixel]
CCD y [ p i xe l ] x [ p i xe l ] D CCD x [pixel]
CCD y [ p i xe l ] y [ p i xe l ] D Figure 3.17:
Average position uncertainty vs position on the frame for stars with cataloguemagnitudo m − . m . The dependence shows rotational symmetry, which is consistent withPSF change along radial axis in the axially symmetric lenses. Analysis of the photometry and astrometry uncertainties show, that they depend on pa-rameters such as position on the frame. It seems very likely that this dependence is causedby a change of the PSF on the frame – a fact known, but never carefully studied in anexperiment with such a big field of view as “Pi of the Sky”. Using rotationally symmetricalapertures or PSFs for strongly deformed, definitely not rotationally symmetrical, imagesof stellar objects is a very risky approach. Therefore analysis of shape of the PSF and itsdependence on the position on the frame is well justified, if not for the hope of improvingphotometry and astrometry then for better understanding of the “Pi of the Sky” detectorand future detectors with such FoVs.
PSF parametrisation idea
To analyse the shape of the PSF and its variation with the position on the frame, the shapeitself has to be determined first. Images of the stars observed on frames are convolutionof the PSF and the CCD pixel response. As the pixel size is quite large for “Pi of the Sky”camera, most of single star images even at the very corner of the frame, consist of less than30 pixels, thus giving less than 30 data points for shape analysis (fig. 3.18). While it ispossible to fit a scale and position of a well known profile to this amount of data, it wouldbe very hard to derive a profile’s shape itself. This derivation requires a higher resolution,much better than pixel size. Assuming that a star is placed in a slightly different positionon each frame, relative to pixel centre, a series of such images can result in sub-pixelresolution. Additionally, one can assume that the PSF is invariant under the rotationaround the frame’s centre and take into account all the stars within a certain distancefrom it. This operation performed on a long series of frames should give an average profile
Figure 3.18:
An image of a single star, shown as a “lego” plot. The star is observed close tothe frame corner. The signal S is spread over less than 30 pixels, thus making its shape difficultto extract. hapter 4Cameras laboratory measurements
Astrophysical objects that are of the interest for the “Pi of the Sky” experiment are sofar away, that from the detector’s point of view, can be treated as point sources. Forthe corresponding measurement in the laboratory the apparatus should consist of a CCDcamera imaging a source which can be considered point-like – an object for which animage on the sensor is much smaller than a pixel size. This condition can be fulfilled bya small enough artificial source of light placed far from the camera.Additionally, a high resolution point spread function measurement requires precisecontrol of the relative position of the light source image on a CCD pixel. This couldbe accomplished by moving the camera or the light source. Considering camera shape,dimensions and weight, the latter occurs to be an easier solution. However, cameramovement is the easiest way to perform a big change of the image position on CCD, whenthe PSF is to be measured in a different parts of the sensor. Based on these assumptions,a dedicated laboratory equipment has been prepared.
The experimental setup consisted of a LED diode (red, green, yellow, blue or white) placedbehind a pinhole of . mm diameter at a distance of m from a CCD camera . A pixelsize of the CCD sensor was × µ m , and the camera was equipped with CANON lenseswith focusing length of mm (fig. 4.1) the same as in the “Pi of the Sky” prototype andthe final system. The expected geometrical diameter of the image, neglecting diffractionand assuming perfect optics is 1.5 µ m (0.1 pixel).The diode was placed in a mechanic mount, driven by two step motors, that alloweda precise movement in vertical and horizontal axis (as shown in fig. 4.2). A constantly Figure 4.1:
The schematic layout of the setup used for laboratory measurements. The length of the corridor in the Particle Division building. igure 4.2: The mechanical mount with two step motors used for precise positioning of thelight source (LED diode) in the laboratory measurements shining source was much too bright for PSF measurements, even with the lowest possiblevoltage in the diode working range. Thus a pulse generator was used as a power source,so that the diode image brightness on the CCD could be adjusted by changing the pulselength. Starting exposures, driving the step motors and triggering the generator pulsewere controlled by a computer with self-written, dedicated software and scripts. Thewhole setup was placed in the basement corridor in the Soltan’s Institute for NuclearStudies’ building on Hoża 69 in Warsaw.
The main differences between a real star and the diode acting as one are the finite distanceto the camera, artificial triggering of the light emission and the mechanical mounting.Thus, before performing PSF measurements a proper focusing of the lenses had to befound, as well as mechanical and optical stability had to be tested.
Cameras focusing
The first task was to find the focus of the lenses for a light source distant by 22 m,assuming that the star is best focused when its PSF is smallest. As an estimator of PSFsize we chose root mean square (RMS) of a histogram of diode image pixels in horizontaland vertical axis. Measurements were performed for 5 diode colors and 13 focus settings.As shown in fig. 4.3, there is a non-monotonic chromatic dependence in the focusingfunction, red and green colours having the most distant minima. Perhaps lenses, tocompensate chromatic difference in focusing, try to achieve the minimum in the centre ofthe spectrum. Thus edges of the spectrum – red and blue – have closer best focus valuesthan the centre of the spectrum – green.There are nearly no differences in the horizontal and vertical RMS for the diode ofthe same colour in the central part of the CCD, where the PSF is nearly symmetrical.Such differences are expected and visible for the deformed PSF close to the edge of theframe, but minima of the focusing functions remain independent on the choice of the axis.There is, however, a difference in focusing minimum between the central and the edge40 igure 4.3:
Size of the diode image on the CCD as a function of the focus setting. Left column:horizontal RMS; right column: vertical RMS; top row: for the centre of the CCD sensor; bottomrow: for the edge of the CCD sensor. Line colours represent diode colours, black standing forthe white diode.
PSF, which is smallest for the red diode. Because of the red diode giving the smallestside effect and the fact that the sensor is most sensitive in the red band we decided to usethis colour best focusing f s = 1 . m as the general best focusing setting . Illumination stability
Reconstruction of high resolution PSF from single frames (sec. 4.4.1) requires a diode witha stable light intensity as well as stable measurement of the signal. The latter involvesestimating behaviour of the background coming from ambient light and other (howeverdim) light sources.The background does not need to be constant during the whole measurement serie, andsuch a condition would be difficult to achieve. However, it should not change significantlybetween so called “corridor frames” – background frames taken without the diode pulse,subtracted from the frames with the diode flashing. The background becomes stableabout 300 s after switching off corridor lights . After this time there is no visible trendin environmental light intensity, larger than frame to frame fluctuations (fig. 4.4, left). The best focus setting of lenses for an object distant by 22 m is 1.4 m here due to the fact, that the“Pi of the Sky” cameras are constructed in a different way than standard digital cameras, to which thescale on the lenses refers. The distance between the CCD and the lenses in our case is larger, to allowmore freedom in focusing at real stars. The main factor responsible for the background not being stable immediately after switching offlights is probably deexcitation of some fluorescent signs in the corridor. [s]0 100 200 300 400 500 600 700 B Figure 4.4:
Illumination stability. On the left: corridor background ¯B calculated as an averagepixel value for frames as a function of time t after switching off corridor lights (uncertaintiesare smaller than the marker size); on the right: histogram of the diode signal S calculated in100x100 pixels region around the diode, after subtracting background (signal level without thediode flash). The amount of light emitted during a single diode flash depends upon the shape ofthe current pulse coming from the generator, the diode response to such a pulse andPoisson statistics of photons. According to fig. 4.4 (right), the fluctuations of the sourcelight intensity (after subtracting background) are on the level of . . Nearly ofthose fluctuations can be explained by photon statistics. The results may vary slightlydepending on the pulse length and other conditions, however the observed fluctuationsare in general negligible compared to other uncertainties introduced during further dataanalysis. Position stability
Two series of measurements were performed to analyze the stability of the light sourceposition in time. One serie was taken with mechanical shutter constantly opened (thuswithout vibrations induced by opening and closing of the shutter), second in normal(opening and closing) shutter mode. Fig. 4.5 shows slight change of the position of thelight source image on the CCD, which is in general lower than . pixels per minute (tab.4.1). These results rule out the possibility, that vibrations due to the shutter movement arethe main reason of the position instability. The question is, why vertical shift in positionis bigger for opened shutter? Frames with opened shutter were taken earlier, shortly aftersetting up the camera and the diode mount. Elements of the setup were probably stillin a phase of mechanical stabilization due to thermal changes or gravitation, the latterhypothesis being favored by the observed difference of shift in vertical and horizontal axis.Measurements taken with normal shutter mode were taken after the period of the initialinstability. Increased instability in the time period just after setting up the machinerywas also observed during other data taking series.Although we were unable to remove the position instability during the data taking,the expected shift of position is lower than . pixel on a 100 frames measurement serie(about 6 minutes) and is therefore negligible.42 [s]0 50 100 150 200 250 c m s x [ p i xe l ] c m s y [ p i xe l ] c m s x [ p i xe l ] c m s y [ p i xe l ] Figure 4.5:
Light source position stability (cms x, cms y) over time t. Top row: constantlyopened shutter; bottom row: normal shutter mode. Results of a linear fit to the data points,indicated in the plot, are given in tab. 4.1.
Fig. 4.6 shows pictures taken during tuning of the laboratory equipment. Fig. 4.6(a) is asample “corridor frame” – frame without the diode blink - used for subtracting backgroundfrom the images of the light source. The details of the corridor are quite well visible,although their absolute brightness is small. This shows the capabilities of the camera,especially very high signal to noise ratio. Fig. 4.6(b) shows a saturated diode flash withsubtracted “corridor frame”. Such a bright blink produced many reflexes visible in the fig.4.6(c) – around the actual diode image. Such photographs allowed removing “light leaks”in the setup, which resulted in very clean diode images, such as the one in fig. 4.6(d).After performing laboratory setup tuning and stability tests, the main stage of data-takingcould be started. horizontal [pixel per minute] vertical [pixel per minute]opened shutter (4 . ± . · − (18 . ± . · − closed shutter (8 . ± . · − (8 . ± . · − Table 4.1:
Estimated light source movement ratio – results of a linear fit shown in fig. 4.5. a) “Corridor frame” (b) Saturated diode flash with subtracted“corridor frame”(c) Saturated diode flash magnified (d) Not saturated diode flash Figure 4.6:
Photographs taken using CCD camera during the tuning of the laboratory equip-ment.
The way the CCD sensor is designed causes a single-pixel light sensitivity to be notspatially uniform. That is mainly due to electrodes placed across the pixels and channelstops separating the sensor’s columns. The non-uniformity can be measured with a sourceof light focused in a spot smaller than the pixel size[54]. The geometrical size of the spotin the described apparatus setup is smaller than the pixel size, but the PSF causes thelight to be spread over several pixels. The way to restrain the PSF is to put a circularaperture in the front of the lenses. The smaller the aperture, the smaller the illuminatedlenses area and the smaller the PSF. However, a small aperture causes the spot size to bediffraction limited (tab. 4.2).As a compromise between PSF size and the diffraction size an aperture of mmdiameter had been chosen for the inter-pixel measurements.44perture diameter [mm] red light ∆ CCD [pixel] blue light ∆ CCD [pixel]40 0.1 0.0720 0.2 0.1410 0.4 0.27
Table 4.2:
Theoretical diffraction limited size of the light spot on CCD for different aperturediameters.
The CCD sensor with accompanying electronics is a an analogue-digital converter of lightintensity illuminating the camera into so-called ADU counts. As with most such devicesthere are some issues concerning the conversion process. Perhaps the most visible isthe saturation, which is observed when the potential well is flooded with photoelectronsbeyond its capacity, which causes the charge to escape to neighbouring wells, resulting inthe distortion of the picture of the object and the surrounding area (so-called blooming,see fig. 4.6(c)). Retrieving proper data from such a distorted image is close to impossible,that is why objects saturating the detector are rather not being analysed.However, there are some less obvious concerns like the dependence of the outputsignal in ADU counts on the intensity of light illuminating the CCD pixel. In order todetermine this dependence measurements of single pixel signal for different diode lightintensities were performed for 9 pixels close to the sensor centre. The light intensity wasadjusted by triggering different number of diode pulses of the same length.Fig. 4.7 shows the results of the described measurement. We observe that the depen-dence of signal on the light intensity is clearly not linear. While the nonlinearity broadensthe magnitude range of permissible objects, it requires an application of a nonlinear cor-rection for each pixel in order to determine the real signal. However, we also see that thedependence is very close to linear up to the output signal of about 20000 ADU, and, to
Figure 4.7:
Dependence of the average CCD signal S on the source intensity I for 9 selectedpixels close to the centre of the chip (different markers).
Pixel response function (PRF) describes a single pixel signal value as a function of thespot position relative to the pixel edge. In an ideal case of infinitely small spot size anduniform pixel response, the PRF should have constant value inside the pixel and zerovalue outside. In the real case, it dependents on the pixel sensitivity and on the finitespot size, so that the spot may be only partially contained inside the pixel. The lattercauses a narrow but not-negligible transition of PRF from high to low values close to thepixel border (fig. 4.8).However, the function is also non-zero for spot fully outside the pixel. That may becaused by a finite PSF size or diffraction of the spot, restrained by setup parameters, butstill not negligible. The more interesting possibility is that it is caused by a charge diffusionbetween pixels – illuminating a single pixel causes some charge to be accumulated in aneighbouring pixels as well. In that case, the PRF “tails” contribute also to the observedPSF shape.
Pixel sensitivity function is defined similarly to pixel response function, however an overallCCD signal is studied as a function of the spot position instead of the single pixel signal.Changes in pixel sensitivity are the main factor responsible for signal changes caused bythe image movement across the CCD. With the knowledge of the pixel sensitivity functionand the position of the source’s centre on the pixel one should be able to compensate forthis effect, performing more precise measurement of brightness.The overall signal was estimated by a sum of × pixels around the spot centre.Results of the measurement are shown in fig. 4.9. The maximal observed changes in signaldue to the pixel sensitivity non-uniformity for a red diode are more than , and for a Figure 4.8:
Measurements used to extract pixel response function for the blue (left), red (centre)and white (right) diodes. igure 4.9: Measurements used to extract pixel sensitivity function for the blue (left), red(centre) and white (right) diodes. blue diode more than . However, for a normal PSF image (taken without aperturereducing PSF size), which is spread over more pixels, the amplitude of correspondingfluctuations is lower (around ).A visible vertical structure similar for both red and blue diodes is probably causedby the electrodes. The horizontal structure is clearly colour dependent, probably due todifferent penetration depths of light of different wavelengths. Both structures are muchless visible for the white diode, for which the effect is averaged all over the white lightspectrum, and thus much smaller. A high resolution profile for selected coordinates on the frame was obtained using multiplyimages of a diode. Each exposure was taken for a specific position of the diode’s centre,the full set of images was covering × points inside a single pixel. Additionally, 5images with and without a light pulse from a diode were taken in each position, for noisereduction purposes. All the images were superimposed, taking into account coordinatesof each image.The procedure of superposition requires determining the proper position of the singleimage relative to CCD pixels. The task is not straightforward, for estimators of theposition like centre of mass (cms) tend to group around the brightest pixel centre (fig.4.10(a)). This is caused by a convolution of a CCD pixel structure and a PSF shapeanalysed in a finite region. To simplify the discussion, let’s consider a 1D cross-sectionof the PSF along the X axis and a cms in this direction calculated in three neighbouringpixels, the central one being the brightest one. In this case: cms(x) = S ( − , ∗ ( − .
5) + S (0 , ∗ . S (1 , ∗ . (4.1)where S ( a, b ) = (cid:82) ba P SF ( x ) , a signal stored in a pixel between a and b coordinates,depends on the relative position of the PSF centre to the brightest pixel centre. Shownin fig. 4.10(b) is calculated cms as a function of the diode position for one of the singleruns with the diode moving along the X axis (Y coordinate being fixed close to the pixelcentre). The curve is a fit of eq. 4.1 to these points, assuming Voigt profile[55] shape of the Mass in this case is, of course, signal stored in the pixels of the CCD. ms x [pixel]0 0.2 0.4 0.6 0.8 10246810 (a) Distribution of the centre of mass for positionsuniformly distributed in horizontal direction x [pixel]0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 c m s x [ p i xe l ] (b) Measured (points) and simulated (line) depen-dence of the calculated cms x as a function of thereal spot position x. Voigt PSF is assumed forsimulation. cms x [pixel]0 0.2 0.4 0.6 0.8 1 X ( x ) [ p i xe l ] (c) Corrected position X(x) – the integral of thecentre of mass position cms x
Figure 4.10:
Reconstruction of the spot horizontal position relative to the pixel edge, with thecentre of mass method.
PSF. The χ / NDF of the fit is 102, so apparently Voigt profile does not properly describethe PSF cross section in horizontal axis , but the plot shows that the interpretation ofthe cms grouping effect is correct.To obtain a real position of the diode image centre, the grouping effect of the cms hasto be removed. For an unbiased estimator of the spot position the expected distribution(assuming uniform distribution of the true position) should be flat, unlike the cms xdistribution on the fig. 4.10(a). Therefore, cms has to be transformed in such a way asto result in a flat distribution. This can be obtained by integrating the cms distributionnormalized to 1. The corrected estimator of position, X , is then: X(x) = (cid:90) x P ( x (cid:48) )d x (cid:48) (4.2)where x is the position returned by the cms calculation and P ( x (cid:48) ) is the probabilitydensity of cms( x (cid:48) ) . Calculated values of the corrected position, X(x) , are shown as afunction of the cms x value in fig. 4.10(c). Corrected position for any measurement Theoretically, an ideal Fresnel diffraction shape without any aberrations is described by Besselfunctions. However, in the aberrated case of “Pi of the Sky” cameras, Bessel functions are not expectedto describe PSF any better than the Voigt profile. CD x [pixel]
CCD y [ p i xe l ] Figure 4.11:
Positions of the PSF measurements on the CCD surface in sensor coordinates. Themeasurements were performed for 5 angles (for ◦ , . ◦ , ◦ , . ◦ and ◦ from the horizontalaxis), and 5 or 8 distances from the frame’s centre. For each angle, distances between consideredpositions were about 200 pixels. One quarter of the CCD chip was only considered. can be obtained by interpolating reconstructed X(x) . The process of calculating such anintegral involves some additional, technical steps like smoothing the histogram and fittingthe absolute position of the diode image on the CCD . PSF measurements and reconstruction were performed for a white diode for 5 angles and6 distances from the frame centre, covering 1/4 of the CCD, as shown on fig. 4.11, eachstar representing a place on a CCD in which PSF was reconstructed.A significant deformation develops with the distance from the frame centre, causing notonly the shape of the profile to change, but also the area containing the signal to grow. Thearea covered by PSF profile at half-maximum, as shown in fig. 4.12, increases from slightlymore than 1.5 pixel for the central profile to nearly 9 pixels for the profile 1000 pixels fromthe frame centre. Even ignoring the shape deformation, such a growth may be a non-negligible factor increasing uncertainties of aperture photometry. For profile photometrythe situation is even more difficult, because the PSF shape changes dramatically withradius, as shown on fig. 4.13 (for measurements at angle of ◦ ).PSF reconstruction for the white diode is important from the point of view of study-ing a general image shape for real stars, which are hardly monochromatic. However,parametrization in terms of diffraction model (chap. 5) is much easier when a singlewavelength is considered. Thus a similar PSF reconstruction was performed for the blue(fig. 4.14) and the red diode (fig. 4.15), for 800, 1000, 1200 and 1400 pixels from theframe centre along its diagonal. It can be noticed, that PSFs of the white diode containsuperimposed shape features of PSFs of the blue and the red diode, which is expected.However, polychromatic PSFs tend to be larger due to the fact, that different wavelengthsare focused in a slightly different positions on the CCD.High quality photographic lenses of the same type should have PSFs in the same place The relative position of the diode image on subsequent measurements is well known by the positionof the step motors. However, determining the absolute position is a more complicated task, especially fornon-symmetrical PSFs far from the centre of the frame. [pixel]0 200 400 600 800 1000 1200 1400 A [ p i xe l ] Figure 4.12:
The area A above half of the maximum height for PSFs measured along chipdiagonal, as a function of the distance from the frame centre r. close to identical, assuming the same focus is set. However, as mentioned before, focusingfor the real stars is different than the focusing for the diode, and it is close to impossibleto maintain exactly the same focus for different cameras. To study the possible influenceof the focus setting on the PSF shape, measurements for three different focusing wererepeated for the white diode, for 800, 1000, 1200 and 1400 pixels from the frame centrealong its diagonal (fig. 4.16).The area covered by the PSF changes visibly with the focusing setting, but the generalshape remains similar (blue and green parts). The very centre undergoes bigger changes,especially for 800 and 1000 pixels from the frame centre (red part). PSF seems to be bestfocused for the last setting, which was not the best focus for the central profile. Thisconfirms the observation discussed in section 4.2.1 that the best focus changes with thedistance from the frame centre. 50 a) 0 pixels (b) 200 pixels(c) 400 pixels (d) 600 pixels(e) 800 pixels (f) 1000 pixels(g) 1200 pixels (h) 1400 pixels
Figure 4.13:
PSFs of the white diode measured along the diagonal, for 0, 200, 400, 600, 800,1000, 1200, 1400 pixels from the frame centre, as indicated below the plots. a) 800 pixels (b) 1000 pixels(c) 1200 pixels (d) 1400 pixels Figure 4.14:
PSFs of the blue diode measured along the diagonal, for 800, 1000, 1200, 1400pixels from the frame centre. a) 800 pixels (b) 1000 pixels(c) 1200 pixels (d) 1400 pixels Figure 4.15:
PSFs of the red diode measured along the diagonal, for 800, 1000, 1200, 1400pixels from the frame centre. s = 1 . m f s = 1 . m f s = 1 . m Figure 4.16:
PSFs of the white diode measured for different focus settings. Measurements weredone for f s = 1 . , . and . m (left, central and right columns respectively), and for 800, 1000,1200 and 1400 pixels from the frame centre (from top to bottom). hapter 5Diffraction model of PSF A detailed knowledge about the point spread function for specific coordinates can beobtained directly from measurements, by reconstructing a high resolution profile of thepoint-like light source. However, to model the detector’s response, a description of thePSF all over the frame is needed. For it is close to impossible to measure profiles forall positions on the CCD, the only viable solution is to create a mathematical modelcharacterizing the PSF dependence on coordinates on the frame, based on obtained data.Understanding the detector’s response is inevitably connected with understanding howthe image is formed on the CCD. While the conversion of the photons to electrons, chargeprocessing and digitization are well explored for digital cameras, the PSF formation forsuch a wide-field experiment has not been described in literature, according to the author’sknowledge. In this chapter we try to describe the image projected on the CCD as theresult of a propagation of light through camera’s lenses. Describing results of such apassage, if done properly and with high enough precision, should result in a proper PSFmodel.
In the case of a simple optical system, such as in the “Pi of the Sky” project, the obstacleon the path of the light are lenses, whose size is much bigger than the wave length.This allows an assumption that the field inside their opening is similar to the field ofan undisturbed light. Such an assumption validates use of the scalar diffraction theory,which should be an adequate approximation in the described circumstances.In the basic case of image generation, the main calculation performed is the transforma-tion of the front of the electromagnetic wave coming through an aperture (or lenses) intointensity of light on the screen – the image. In this application, a Rayleigh-Sommerfeldequation can be used as the most general formula in the scalar diffraction theory:
P SF L ( x , y , z ) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) iλ ∞ (cid:90) (cid:90) −∞ U ( x, y, ze ikr r d x d y (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (5.1)where ( x , y ) are the Cartesian coordinates in the image plane, distant by z fromthe aperture, ( x, y ) are corresponding coordinates in the aperture plane, r is the distancebetween chosen ( x, y, and ( x , y , z ) (as shown in fig. 5.1), λ is the wavelength, k is thewavenumber and U ( x, y, is the amplitude of the wavefront on the aperture . The amplitude of the wavefront should be constant on the aperture in the perfect case, but in reality igure 5.1: Scheme of the aperture and screen coordinate systems, as used in the Rayleigh-Sommerfeld formula.
Calculating intensity of light in each point of the image requires integrating the wave-front all over the aperture plane. However, in the scalar theory it can be assumed thatthe wave is nonexistent outside the aperture A in the aperture plane. This leads to theKirchhoff approximation: P SF L ( x , y , z ) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) iλ (cid:90) (cid:90) A U ( x, y, ze ikr r d x d y (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (5.2)In the case of many analysed optical systems, the size of the aperture L and thescreen L are much smaller than the distance between them: | z | (cid:29) L + L , or in a moredetailed formulation: L zλ ≥ . Then z becomes dominant over the x and y and the Fresnelapproximation becomes valid: P SF L ( x , y , z ) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) e ikz iλz (cid:90) (cid:90) A U ( x, y, e ik z (( x − x ) +( y − y ) ) d x d y (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (5.3)The Fresnel formula is in fact a paraxial approximation, dealing only with small anglesbetween optical axis and the image. However, it still takes into account curvature of thewavefront, which can be neglected for bigger distances, where z (cid:29) kL . This is the coreof the Fraunhofer approximation: P SF L ( x , y , z ) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) e ikz iλz e ik z ( x + y ) (cid:90) (cid:90) A U ( x, y, e ikz ( x x + y y ) d x d y (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (5.4)The Fresnel diffraction formula (eq. 5.3) (as well as the Fraunhofer formula) can beeasily transformed into an equation involving a Fourier transform: P SF L ( x , y , z ) = (cid:12)(cid:12)(cid:12)(cid:12) e ikz iλz e i πλz ( x + y ) F { U ( x, y, e i πλz ( x + y ) } (cid:12)(cid:12)(cid:12)(cid:12) which can be calculated using Fast Fourier Transform (FFT) algorithms, greatly re-ducing computing time. However, this is not the case for “Pi of the Sky”, where the can be altered by factors such as lenses transmission. fd = 1 . ( f stands for the focusing length and d for the lensesdiameter). Thus both the aperture and the screen sizes are of the order of the size ofthe focusing length, i.e. z ∼ L , L . This makes use of the Fraunhofer and Fresnelapproximations impossible and forces us to use the Kirchhoff formulation (eq. 5.2).The approximation which can be introduced in the “Pi of the Sky” case is an assump-tion of planar wave reaching the lenses, which is applicable for sources of light with almostinfinite distance from the observer, such as stars. Thus in the perfect case of an emptyopening, U ( x, y, could be treated as having a constant value all over the aperture. In the case of perfect lenses, a point-source of light should be visible as a point on thescreen. However, the passage of light through an opening of a finite size diffracts the waveand causes spread of the image over finite area. The point spread function of a point-source visible through an aperture has been described by Airy and can be easily obtainedfrom the eq. 5.2. However, Kirchhoff formula in such form describes only light passingthrough a perfect optical system in which aperture is a “source” of a perfect, sphericalwave. This is not the case in real optical systems, where introduction of such elementsas mirrors or lenses aberrates the sphericity. Thus, in the case of a monochromatic wave,such deviations are called monochromatic optical aberrations.To obtain an aberrated PSF from diffraction formulas, a deviation W ( x, y ) fromsphericity e ikr of the wavefront has to be introduced into eq. 5.2: P SF L ( x , y , z ) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) iλ (cid:90) (cid:90) A U ( x, y, ze ikr e W ( x,y ) r d x d y (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (5.5) W ( x, y ) may be expressed in multiple forms. Perhaps the most referred to is the sumof Seidel aberrations[56] such as defocus, coma, astigmatism, spherical aberration, etc.This formulation is often used by optical engineers, mainly due to a quite convenient wayof transforming into optical systems parameters.The Zernike concept is generally more popular in the light wavefront analysis appli-cations, mainly due to its mathematical properties. It introduces a set of polynomials,with names corresponding to aberrations in the Seidel formulation, which however areorthogonal on the whole aperture, making this formulation well suited for the wavefrontfitting purposes. In this approach the wavefront can be expressed as: W ( ρ, φ ) = (cid:88) n ≥ m Z mn ( ρ, φ ) + (cid:88) n ≥ m Z − mn ( ρ, φ ) (5.6)where Z mn and Z − mn stand for even and odd Zernike polynomials respectively and n, m are non-negative integers, denoting the order and type of the aberration[57]. Thewavefront distortion function W was expressed in polar coordinates, more natural forradial aperture, where radial position ρ is normalized to the aperture radius, < ρ < .Each polynomial can be written as: Z mn ( ρ, φ ) = R mn ( ρ ) cos( mφ ) , Z − mn ( ρ, φ ) = R mn ( ρ ) sin( mφ ) (5.7)the radial part being: 57 berration n m polynomial wavefront astigmatism 2 2 ρ cos 2 φ coma 3 1 (3 ρ − ρ cos φ spherical aberration 2 0 ρ − ρ + 1 trefoil 3 3 ρ cos 3 φ secondary coma 5 1 (10 ρ − ρ + 3 ρ ) cos φ secondary spherical aberration 6 0 ρ − ρ + 12 ρ − Table 5.1:
Monochromatic optical aberrations and their Zernike polynomials used in the diffrac-tive PSF modelling. R mn = ( n − m ) / (cid:88) k =0 ( − k ( n − k )! k !(( n + m ) / − k )!(( n − m ) / − k )! ρ n − k In fact, Zernike polynomials are compatible with Seidel aberrations if only aberrationsup to rd order are used[58]. A few terms presenting aberrations used later in this workare shown in the tab. 5.1. A starting point for obtaining the model of the point spread function in the “Pi of the Sky”project is to find a proper parametrisation for profiles reconstructed from laboratory mea-surement. Thus the eq. 5.5 was transformed into cylindrical coordinates, with astigma-tism, coma, secondary coma, spherical aberration, secondary spherical aberration andtrefoil aberrations used : P SF L ( ρ , φ , z ) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) iλ (cid:90) (cid:90) A U ( ρ, φ, ze ikr e W ( ρ,φ ) r ρ d ρ d φ (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (5.8) These polynomials were chosen as those capable of providing significantly diverging contributionsinto the PSF, while still forming a rather small set. W ( ρ, φ ) = A ( ρ cos 2 φ ) + C ((3 ρ − ρ cos φ ) + S (6 ρ − ρ + 1)+ T ( ρ cos 3 φ ) + C (cid:48) ((10 ρ − ρ + 3 ρ ) cos φ )+ S (cid:48) (20 ρ − ρ + 12 ρ − (5.9) A, C, S, T, C (cid:48) , S (cid:48) stand for coefficients of polynomials representing respectivelyastigmatism, coma, spherical aberration, trefoil, secondary coma and secondary sphericalaberration. This is the maximal set of aberrations used in the modelling, however it canbe easily extended to a bigger set. Polynomials used commonly for presenting deviationfrom the optical axis (tilt) or defocus are not needed here, since they are naturally presentin the Kirchhoff formulation via coordinates of the image.There are no known methods of analytically performing integration of eq. 5.8, espe-cially in the presence of aberrations. Thus numerical integration has to be performed,which in this case is a complicated task. Due to the large wavenumber k in the term e ikr (for any optical wavelength), the formula consist of a very fast oscillatory integral.Common adaptive integration methods fail at attempt of performing such a calculation.Other methods have to be used, perhaps one of the most obvious and giving good resultsis an integration using Gauss integration points in a two dimensional space.The quality of results strongly depends on the number of intervals used in the integra-tion. Lenses as wide as in the “Pi of the Sky” case (few centimeters) require hundreds tothousands intervals to obtain satisfactory result, the number growing with the distanceof the calculated PSF from the frame centre.Altogether this makes the PSF computing a very time consuming task, reaching sec-onds on modern computers (for example of PSF with 250 points). Finding proper modelparameters with fitting procedures often requires thousands of such calculations, makingthe modelling close to impossible.Fortunately, a work by Miks et al.[59] makes the task of integrating rapidly oscillatingfunctions in the diffraction equations much less time consuming. They obtain an accurateapproximation (based on a Taylor series expansion of the function in the exponent), whichcan be formed as: (cid:90) (cid:90) S f ( x, y ) e ikg ( x,y ) d x d y (cid:39) (cid:88) n f ( x n , y n ) e ikg ( x n ,y n ) sin X n X n sin Y n Y n ∆ x n ∆ y n ,X n = k ∂g ( x n , y n ) ∂x ∆ x n , Y n = k ∂g ( x n , y n ) ∂y ∆ y n (5.10)This formula allows obtaining better results and is a few times faster than previouslymentioned methods – it was used for the purpose of this work.At this point we have to take into account, that the eq. 5.8 describes only the PSFcoming out of the optical system – P SF L , while the image is the convolution of the P SF L and the CCD response P RF : P SF ( ρ (cid:48) , φ (cid:48) , z ) = (cid:90) (cid:90) P RF ( ρ , φ , ρ (cid:48) , φ (cid:48) ) · P SF L ( ρ , φ , z )d ρ d φ (5.11)where ( ρ (cid:48) , φ (cid:48) ) are polar coordinates of the specific point of the final P SF (correspond-ing to CCD pixel centres) in the image plane (the result of the convolution), while ( ρ , φ ) are polar coordinates of the current point of integration in the image plane (same as ineq. 5.8). 59 a) simulated PSF (b) real PSF Figure 5.2:
The PSF description obtained from the diffraction model (left) and the recon-structed PSF 800 pixels from the frame centre (right). The overall shape given by the model isthe same as the measurements, showing that the diffractive approach is proper in general, butit fails to reproduce details.
For modelling purposes, it was assumed that the camera’s CCD-lenses setup can be ap-proximated by an aperture-screen setup, with aperture of 2.5 cm radius (the exitingwindow of the lenses), and the basic distance between screen and aperture of 8.5 cm (thefocusing length of the lenses). At this stage, it was assumed that the amplitude of thewave U ( x , y ,
0) = 1 all over the aperture.The eq. 5.8 approximated by eq. 5.10 and convoluted as in eq. 5.11 was used tomodel PSFs reconstructed in laboratory measurements (fig. 4.14 and 4.15). In generalcase the convolution would require integrating over the infinite space. However, the rapiddrop of the PRF close to the pixel edge allows restricting the integration space. An areaof . × . square pixels around the pixel centre was used, the integration replaced bya simple sum of P RF and
P SF L products on the uniformly distributed net of × point in this space.Model was obtained by fitting, with the least χ method, coordinates of an image ( x, y, z ) (slight deviations from values obtained from the experiment were allowed) andaberration coefficients of formula 5.9 to reconstructed profiles. The fit was performedstarting in the centre of the frame and moving to the expected image position. The firstreason for using such a method was that for large angles the deformation of the PSFcaused by the plain deviation from the optical axis needs to be compensated by comaaberration and the proper coefficient could not be guessed. The second reason is a resultof an observation, that the fitted model may, probably due to some simplifications in themodel, describe the PSF better in a position different than it was actually registered.Figure 5.2 shows the fit result for the red diode PSF measured 800 pixels from theframe centre. Model is very similar to the measured point-source profile, showing, that thediffractive approach is, in general, the proper method of PSF modelling. Unfortunately,60 a) simulated red PSF (b) measured red PSF(c) simulated blue PSF (d) measured blue PSF Figure 5.3:
The PSF description in the diffraction model for red (a) and blue (c) diode, obtainedfrom the simultaneous fit to the red and blue PSF profile measurements results 800 pixels fromthe frame centre, shown in (b) and (d) respectively. the obtained model fails to reproduce exact shape of the PSF centre (shown in red) and,even more, the tails of the profile (shown in blue).The other problem in this approach is that the obtained results do not necessaryrepresent the global minimum of χ . During calculation we noticed, that many local χ minima exist in the parameter space and it is not guaranteed that the chosen “path” forparameters was the correct one.To constrain the fit, an attempt to simultaneously minimize red and blue PSFs wasperformed, based on an assumption, that the aberrations coefficients do not depend on thewavelength. However, results of that procedure were worse for the red PSF, 800 pixels fromthe frame centre than in the single wavelength attempt (fig. 5.3). Fits for larger distanceswere even more diverse from profiles reconstructed from laboratory measurements.61 a) uniform transmission (b) nonuniform transmission Figure 5.4:
Comparison of the PSF description obtained from the model with uniform lensestransmission (left) and with transmission changing with radial coordinate ρ on the aperture as .
73 + 0 . ρ (right). To improve agreement between the model and the data we also considered a fit wherethe U ( x , y , function was modified in a final stage of model tuning to include thetransmission function of the lenses. This however had a small impact on the overall shape,as shown in fig. 5.4. The number of intervals used in the integration also had no impact onthe fitting procedure, for increasing accuracy did not alter parameters significantly. Theinconsistency is definitely not due to approximations used in the calculation procedure.There are a few reasons why the diffractive approach failed to reproduce PSFs inhigh enough detail. The first is the time of computation – although much shorter thanin standard approach, thanks to approximation shown in eq. 5.10, still very long for ageneral search for the shape of the PSF, especially in a situation of multiple local minima.Additionally, fitting procedures used tend to find minima close to the initial parametersof the fit, while the initial parameters here (except for the image coordinates) were a pureguess. The procedure of starting in the frame centre, with aberrations coefficients close tozero and gradually changing them while proceeding to the final image position is also notguaranteed to give proper results. More advanced minimization algorithm, more likely tofind global minimum, should perhaps be used, such as a genetic algorithm or simulatedannealing. However, such algorithms require much longer computation times.It is also possible that the introduced model is too simple to reproduce profiles insufficient detail. A simple aperture-screen setup is based on an approximation of thinlenses, while real photographic lenses definitely do not meet this assumption. In ourmodel aberrations are added only in the aperture plane, while in reality they are formedin multiple planes.It has to be stated, that although the approach failed, it shows promising resultsfor further study. A number of dedicated measurements could be performed, such asdirect measurements of the wavefront, that could aid the modelling. Such a study couldresult in a possibility of modelling PSFs taking into account the spectral type of theobserved source, which in some special cases allowing long computing time, could provideinteresting results. 62 hapter 6Polynomial model of PSF The point spread function describing the star image visible on the frame is a diffractiveimage of the source convoluted with the CCD pixel structure. Thus the most proper, inscientific terms, way to describe the PSF is to find parameters of the lenses and CCD andperform a diffraction calculations, as described in chapter 5. However, such an approachcan only be used for an idealized, simplified version of an experimental setup. Moreover,in this case calculations of a diffractive pattern are a significantly time-consuming task.This makes the main aim of this work – an application of the “diffractive approach” todata analysis and simulation – a troublesome undertaking.A common way to simplify the use of physical models in practice is to introduce aneffective approach. In the case of PSF description a parametrization in the image planecan be sought – directly deriving it from the shape of reconstructed, high resolution point-like source’s profiles. The integration of the wavefront is skipped, this way leaving onlysimple one-time function calculation for each of the PSF points. However, no simplefunction describes aberrated PSF in the image plane (except for some very simple cases).Therefore, to obtain an efficient effective profile a model has to be developed and for thistask a proper mathematical basis has to be chosen.
Deriving the PSF description in the image plane is a task similar to finding a mathemat-ical description of complicated shapes of an aberrated wavefront. Thus similar basis offunctions can be used in both tasks – Zernike polynomials. However, in the case of theimage plane they need to be modified slightly.First modification is due to the fact that the PSF is a function described on an infiniteplane, while the wavefront described by the Zernike polynomials was bound to a finite,circular aperture. Therefore, a transformation of a PSF’s radial coordinate r to argumentof Zernike polynomial, u ( r ) , where we assume that u (0) = 0 and u ( ∞ ) = 1 , has to beperformed. Thus set of modified Zernike polynomials is now defined as: Z mn := Z mn ( u, φ ) , u = 1 − e − rλ (6.1)where u and φ are standard radial and azimuthal coordinates of Zernike polynomialsand λ is a parameter modifying the transformation. Zernike polynomials are defined asin eq. 5.7.Additional modification is needed, for the real PSF has a maximal value around itscentre and asymptotically drops to zero in infinity, while the wavefront has a sharp cut-off63t the border of the aperture. The asymptotic behaviour is introduced to the PSF usinga gaussian-like modification of the Zernike polynomials: P SF L ( r, φ ) = e − . · (cid:80) m,n Z mn ( u,φ ) · r p (6.2)where P SF L stands for a PSF generated by lenses only (no convolution with CCD) and p is a factor modifying the asymptotic behaviour. Additionally, we assume that the PSF issymmetric, thus only symmetric terms are allowed in Zernike polynomials. The obtainedset of functions is not orthogonal, in contrast to the standard Zernike polynomials, butstill is well suited for fitting purposes. However, one has to keep in mind that in thisapproach hardly any physical meaning can be attributed to the parameters of the fit.Similarly to the diffractive approach, PSF visible as an image on the frame is a convo-lution of a point spread function generated by lenses (eq. 6.2) and CCD pixel structure: P SF ( r (cid:48) , φ (cid:48) ) = (cid:90) (cid:90) P RF ( r, φ, r (cid:48) , φ (cid:48) ) · P SF L ( r, φ ) drdφ (6.3)where P RF ( r, φ, r (cid:48) , φ (cid:48) ) stands for pixel response function.In the general case convolution of P SF L and P RF , given in eq. 6.3, should be in-tegrated over the whole CCD. In the real case the integration is performed numericallyon interpolated PRF values (similarly as in the diffractive approach – see sec. 4.3.2) –summed are products of × P SF L points and bilineary interpolated P RF points inthe range of . × . pixels around a pixel centre, which corresponds to the PRF size.The net of points is distributed and weighting for points chosen according to gaussianquadrature integration rules. On of the main aims of this work is to find a universal formula describing point spreadfunction all over the frame. Due to limited time, the function was measured only for afew chosen positions (sec. 4.4.2). Therefore deriving such a universal formula requiresfinding a local model for each of the measured PSFs first and then finding an interpolationmethod between these models.Fitting eq. 6.3 to a reconstructed white diode profile gives very satisfactory results,especially when many (for example first 20) symmetrical Zernike polynomials are used asa basis. However, the bigger the basis the more complicated the mathematical descrip-tion – every extra polynomial means one extra dimension in the parameter space for theinterpolation, not mentioning extra time needed for calculations. Moreover, some poly-nomials of a small significance may lead to distorted results in interpolated PSFs. Thus aproper basis – set of Zernike polynomials that give significant contribution to the shapemodelling needs to be chosen.
The first conclusion of PSF fitting was that p and λ values should be fixed. Theseparameters are highly correlated with other free parameters and cause fits to be unstable.Moreover, the asymptotic behaviour adjusted by them is also altered by the circularZernike Polynomials, i.e. terms 0, 4, 12 and 24. Thus p = 2 and λ = 1 values were chosenwhich are standard for normal distribution, simplify the computation of the eq. 6.2 andproved to work well for all fits. 64o choose optimal terms of Zernike Polynomials, a χ fit of eq. 6.3 to measuredprofiles was performed in a special way. For the first iteration of the fit, the initial valueof parameter 0 (weight for first Zernike polynomial Z ) was set to 1 and all the other to avalue close to 0 . The fit was performed with only the first parameter free, then repeatedwith only the second parameter free and so on until all parameters were tested. Theparameter giving the best χ was chosen for the next iteration, in which the procedureof freeing consecutive parameters was repeated, with the parameter chosen in the firstiteration permanently set as free. In the third iteration parameters chosen in the firstand second were permanently free and so on, until the iteration in which all parameterswere free was performed. This way not only parameters having the largest impact on thefit were chosen, but also their order, important for the future model reconstruction indifferent circumstances. Basis formed of the best parameters for all the diagonal PSFs
The overall fitting procedure could be performed to a single reconstructed PSF as well asto multiple profiles simultaneously, final χ being the sum of χ for single profiles. Thefit to multiple profiles at once was the first attempt to find the optimal basis – performedfor all the diagonal PSFs. While the results of modelling for × pixel area around theprofile maximum were very satisfactory (fig. 6.1), no simple dependence of polynomialterms on the distance from the frame centre could be obtained (fig. 6.2).The attempts to obtain a simpler dependence were numerous. Mostly they focusedon manually restraining the parameter space, such as removing parameters that were lessimportant for some modelled PSFs or fitting separately PSFs at 0 to 600 pixels from theframe centre and 600 to 1400 pixels from the frame centre. Additionally, a constrain thatparameters must follow a parabolic or cubic dependence on the distance from the framecentre was tested. The latter never gave simple dependence on the distance, the formerresulted in unsatisfactory shapes.However, one thing is apparent in most of the described fits – parameters values re-sulting from a fit to PSF at 0 and 200 pixels from the frame centre differ significantly fromvalues for other distances in nearly all cases. Moreover, the fitting procedure attempts tobuild the circularity of strongly circular PSFs – such as these 0 and 200 pixels from theframe centre – with non-circular terms. Those two factors led to a slight modification ofthe procedure – 3 circular terms were always fitted at the beginning.The results of parameter reduction with 3 circular parameters forced into the fit im-proved vastly the circularity of PSF at 0 and 200 pixels from the frame centre and de-scribed well all the profiles with overall number of 17 polynomial terms. The deviationof parameters for PSF at 0 pixels from the frame centre from results at other distanceswas vastly reduced in most cases, but still the parameter vs distance dependency couldnot be described by any simple function. The dependence of two selected fit parameterson the distance from the frame centre is presented in fig. 6.3. Basis formed of the least and the most deformed PSFs
The need to include 3 circular profiles led to an idea, that the basis should include termscharacteristic for the least and the most deformed PSF, i.e. 0 and 1200 pixels fromthe frame centre. Thus choosing of the best parameters was performed for these two For all parameters set to 0 (or 1 in case of Z ) a fit is highly unstable. The fit was more stable,when other parameters values were not initially equal to 0 (or 1), but close to it. Figure 6.1:
Measured PSFs (left column) and their polynomial models (right column), resultingfrom the parameters searching and fitting procedure, for PSF measured 0 to 1400 pixels fromthe frame centre, as indicated in the plot. Continued on the next page. Figure 6.1:
Continued... profiles independently, results joined into a single, 15 parameter basis. Unfortunately,this approach also failed to describe parameters vs distance dependency in any simpleway.
Both basis, the first formed of 3 circular polynomials and those leading to best χ for allPSFs and the second formed by the best terms for the least and the most deformed PSFs,describe measured profiles comparably well and have comparably complicated parameterdependence on the distance from the frame centre. We did not succeed to implementthe dependence in a simple and consistent mathematical model, due to multiple extremanot overlapping in the parameter space (see fig. 6.3) and a large number of parametersrequired to properly model measured profiles.Therefore, for the effective approach, a linear interpolation in distance between pa-rameters was chosen. The linear interpolation is not a perfect choice, for in some cases itleads to a sudden change in a parameter value, but does not cause a suspicious behaviour– false extrema – as in all attempted spline interpolations (example shown in fig. 6.4).67 [pixel]
200 400 600 800 1000 1200 1400 p r [pixel]
200 400 600 800 1000 1200 1400 p Figure 6.2:
Coefficients of Zernike polynomial terms 13 and 14 as a function of the distance fromthe frame centre r – results of the best parameters searching and fitting procedure. No simpledependence of the parameter value on the distance from the frame centre is visible. Results forother parameters are included in appendix A.1.
To test the selected interpolation method, a dedicated fitting procedure was performed.All the measured PSFs were fitted with one of the obtained parametrizations, but withradial coordinate on the frame, determining the shape, set free, although initialized to aproper value. For the comparison purposes, the test was performed for two parametriza-tions – the first with 3 circular and other best terms for all profiles, and the second formedfrom the best parameters of the least and the most deformed profiles. The dependence ofthe PSF position determined from the fit (based on the shape only) on the actual distancefrom the frame centre is shown in fig. 6.5(a) and the difference between the fitted andthe measured position in fig. 6.5(b).The test shows quite good behaviour for diagonal profiles – the deviation from the mea-sured distance from the frame centre does not exceed 40 pixels for the first parametrization r [pixel] p -1-0.8-0.6-0.4-0.200.20.4 r [pixel] p -0.8-0.6-0.4-0.200.2 Figure 6.3:
Coefficients of Zernike polynomial terms 13 and 14 as a function of the distancefrom the frame centre r for the finally chosen basis of 17 polynomial terms – 3 circular and14 giving best χ for all measured profiles. Results for other polynomial terms are included inappendix A.2. [pixel]0 200 400 600 800 1000 1200 1400 p -0.6-0.4-0.200.20.40.60.811.2 Figure 6.4:
Fitted weights for the Zernike polynomial term 30 as a function of the distancefrom the frame centre r (points). Also indicated are results of linear (blue line) and spline (redline) interpolation. Spline interpolation in this case leads to values of parameters far from theseobtained from the fit, thus linear interpolation is safer. and 60 for the second. For other PSFs the result is worse, reaching almost 200 pixels devi-ation for some PSFs in the second parametrization and 100 for the first. This is probablydue to the fact, that the diagonal parametrization does not describe well PSFs for differentazimuthal coordinates on the frame – PSFs are not symmetric in respect to the centre ofthe frame, which was not anticipated at this stage of analysis.In general, the first parametrization gives better results in the distance from the framecentre tests and provides model PSFs more resembling reconstructed profiles. Thereforethis parametrization was chosen for the further use. r [pixel]0 200 400 600 800 1000 1200 1400 [ p i xe l ] r diagonalverticalhorizontal67.522.5 (a) Fitted vs measured distance r [pixel]0 200 400 600 800 1000 1200 1400 - r [ p i xe l ] r -300-200-1000100200 diagonalverticalhorizontal67.522.5 (b) Deviation of the fitted from measured distance Figure 6.5:
Consistency check of the obtained parametrizations – distance ¯r resulting from thefit of model PSF to the PSF profile measured at the distance r from the frame centre. Emptymarkers denote least and most deformed PSFs base, filled the best base for all PSFs with 3circular terms. a) 0 pixels from the frame centre(b) 200 pixels from the frame centre(c) 400 pixels from the frame centre(d) 600 pixels from the frame centre Figure 6.6:
Final results for the point spread function modelling. Shown for different distancesfrom the frame centre are: measured profile (left), polynomial model of the PSF convoluted withthe CCD pixel response (centre) and the obtained shape of the optical PSF before convolutionwith the CCD pixel structure (right). a) 800 pixels from the frame centre(b) 1000 pixels from the frame centre(c) 1200 pixels from the frame centre(d) 1400 pixels from the frame centre Figure 6.6:
Continued... .3 Results of the PSF model As shown in fig. 6.6, chosen parametrization describes the central × pixels of the mea-sured PSFs very well, especially up to 800 pixels from the frame centre. For comparison,results for the optical PSF shape, before convolution with the CCD response function,are also shown. For 1000-1400 pixels from the frame centre proportions of the core (redarea) seems a little bit different then for measured profiles – cores are slightly too thick,especially for PSF 1400. Additionally, the more complicated the shape of the PSF, theworse the reproduction of the deep details – the fact visible most at the distance of 1200pixels from the centre of the frame. However, we have to realize that the differences areonly visible on sub-pixel scales. Therefore we conclude that the obtained model is highlysatisfactory.Deconvoluted profiles smoothly change their shape with the distance from the framecentre, showing that the polynomial base is well suited for interpolation. With such resultsit is very tempting to assume that the deconvoluted fit results are a good approximation ofthe real PSF not convoluted with the CCD pixel structure. Such an approximation wouldbe a great base for testing different parametrizations or even the diffractive approach, forthe whole time needed for calculation of the convolution could be saved. However, this Figure 6.7:
Results of a fit to the PSF measured 400 pixels from the frame centre with differentparametrizations: with 3 circular and other best terms for all profiles (top) and the best terms ofthe least and the most deformed profiles (bottom). Left: resulting model convoluted with
P RF ;right: optical
P SF L only, not convoluted with P RF . Within the available data accuracy goodresults of the final, convoluted fit does not imply similar non convoluted PSFs.
72s not the case – the available accuracy allows to create different, similarly good modelsof the measured profile with very different deconvolution results (an example is shown infig. 6.7).
Parametrization used to describe point spread function measured for a white diode canalso be used to model monochromatic or defocused profiles. It seems however, that theresults of such modelling are not as good as for the white PSFs, as shown in fig. 6.8 andA.3. The size and shape of the profiles (green area) is roughly reproduced, but the detailsof the core (red part) are similar only for the fit of the white diode with f s = 1 . m,although the core in this case seems to be most complicated. That is probably due to thevery short tails (blue part) of this PSF. As for other profiles, the fitting algorithm focusesboth on the core and tails, failing to reproduce any of them in sufficient detail.Worse modelling results despite the shape similarities are another indication that thepolynomial approach does not follow real physical attributes of the PSF. Thus there werenearly no chances that the focusing or wavelength variations would be followed in anyway in the parametrization. Even though, some attempts to check, if this is really thecase have been made.In some cases a general shape of the profile for colour or defocused PSF is better fittedwith the parameter values obtained for the white, focused profile, but with free p and λ parameters in eq. 6.1 and 6.2, as shown in fig. 6.9. However, this is not the case for mostprofiles.Additionally, a simultaneous fit of 3 focus settings for the same distance from the framecentre has been performed for 800, 1000, 1200 and 1400 pixels. No simple dependence ofparameter values on focusing and distance from the frame centre has been observed, eventhough the analysis was performed for just 3 focusing settings. This study shows, that nosimple modelling of chromatic and focusing changes can be obtained with the polynomialapproach. The final test of the obtained PSF model for the white diode was a check of its suitabilityfor photometric purposes. Obtained model was fitted to the single diode images takenfor the PSF reconstruction purposes, for all available distances from the frame centrealong the frame diagonal. Two fit types were performed – based on χ and loglikelihoodminimisation. Additionally, an ASAS photometry and a simple aperture photometryresembling the fast photometric algorithm used in the “Pi of the Sky” were also tested onthe same data set.Results presented in fig. 6.10 show that the obtained model of the PSF can servethe photometric purposes very well. There is no clear dependence of the brightnessmeasurement spread on the distance from the frame centre for the model fits as wellas ASAS photometry. Some dependence is visible for the “Pi photometry”, probably dueto the fact, that this method uses a single aperture size, while the diode PSF size changessignificantly with its position on the frame (see fig. 4.12).73 a) Blue diode with nominal focus f s = 1 . m(b) Red diode with nominal focus f s = 1 . m(c) White diode with focus changed to f s = 1 . m(d) White diode with focus changed to f s = 1 . m Figure 6.8:
Measured PSFs for colour diodes and for the white diode with changed focusing(left column) compared with their models obtained with the developed parametrization (rightcolumn), for distance of 800 pixels from the frame centre (results for other distances are shownin appendix A.3). a) Measured PSF (b) Refit of the chosen basis (c) Fit improve attempt Figure 6.9:
Comparison between measured PSF for the red diode, 800 pixels from the framecentre (a) and its models resulting from two different fitting approaches (b and c). Model (b) isa fit of the standard polynomial basis, as obtained for the white diode, while model (c) was thestandard basis with parameter values taken from the white diode fit, but with refitted p and λ parameters of eq. 6.1 and 6.2. r [pixel]200 400 600 800 1000 1200 1400 I/I D Standard fitLoglikelihood fitAsas photometryPi photometry
Figure 6.10:
Photometry of the diagonal white diode images (used for diode PSF reconstruc-tion). Black and red points stand for the fit of the obtained PSF model (respectively χ andloglikelihood fits), green for the ASAS photometry and blue for the single aperture photometrybased on the fast “Pi of the Sky” photometric algorithm. Error bars for the fitted model aresmaller than markers. hapter 7Modelling real sky data Laboratory measurements of the point spread function were a necessary step for obtaininga model properly reproducing images of the stars seen by “Pi of the Sky” cameras. Themodel described in the previous chapter successfully reproduces the PSF of an artificialpoint source and gives very promising photometry results. However, it has to be comparedto real sky data to find out if it achieves its goal.A quick glance at data from real sky images taken with the prototype “Pi of the Sky”detector in Chile shows that the star images are not identical in shape to the laboratorymeasurements, although have similar general properties. This difference can result fromdifferent focusing of lenses used in laboratory and sky measurements. Moreover, starimages show that the real PSF does not change with distance from the frame centre only,but also strongly depends on the azimuthal coordinate. This has to be due to mechanicaldifferences appearing in the assemble procedure – the dependence on azimuthal coordinateis most likely a result of a non-negligible tilt of the optical axis with respect to the theCCD (sensor board) plane.Different focusing and slight, at a first glance, but non-negligible after deeper studydeviations in cameras assembly enforce a recalculation of model parameters for each set ofequipment. However, this is an iterative procedure and requires use of the model obtainedin the laboratory measurements at the first step. Recalculation of the model parameters requires obtaining high resolution profiles of realstars. This task proved difficult if not impossible for highly deformed PSFs withoutknowledge of their general shape (as discussed in sec. 3.6.3). However, intermediatemodel obtained in laboratory measurements gives a general shape of stars for differentcoordinates on the frame and can be used for obtaining the high resolution profile.A set of 285 bright stars up to m , without bright neighbours, scattered all over theframe was chosen. For every star scale (signal), position and background level were fittedwith the intermediate model on each of the consecutive 172 frames of a single field. Usingobtained positions, a high resolution profile of each star was prepared, similarly as inthe chap. 4.4.1. The results of this reconstruction were successful and, after subtracting Camera used in laboratory experiments and those taking real sky images in Chile were identical intechnical specification and were equipped with the same lenses. igure 7.1: Dependence of the selected model parameters on the star position on the frame.On the left: in Cartesian coordinates, on the right: after transformation to polar coordinates.The azimuthal angle φ is calculated from (0,1) versor, counterclockwise. neighbouring stars were more than sufficient for model parameters refit.Similarly as in the intermediate model, PSF has to be described for every position onthe CCD, while fits were performed only for a finite set of stars with specific positions onthe sensor. Thus each parameter value has been plotted in Cartesian coordinates and aDelaunay triangulation has been used for interpolation between measured points. Resultsfor most of the parameters show general rotational symmetry and it seemed that the polarcoordinate system would be the most natural for the interpolation purposes, described inmore detail in sec. 7.1.2.However, for most parameters the symmetry axis of the plot was shifted from the centreof the CCD, see fig. 7.1. To estimate the actual position of the optical axis on the CCDsensor we assumed that the radial dependence of the parameter should be described as a th order polynomial in distance from the axis and looked for axis position which resulted For the reconstruction purposes, stars without bright neighbours were chosen. However, dimmerneighbours, often not visible during the procedure of star selection, but visible on the high-resolutionprofile made from 172 frames remained.
77n the best description of given parameter (neglecting azimuthal dependence). The newcentre of the frame was obtained from weighted average of centres for all parameters andwas used for transformation into polar coordinates.The model recalculated in described way and interpolated after transformation to polarcoordinates could already be used as a general representation of PSF for the discusseddata. However, to ensure high quality of the parametrisation an iterative approach shouldbe used – the procedure of profiles reconstruction, refitting of the model parameters andinterpolation should be repeated several times. In our case the final model was obtainedafter 6 iterations. However, the number of iterations was induced more by the fact, thatthe procedure was being improved and polished than by the better quality obtained ateach next run. A well known practice for such PSF iterative modelling (although for caseswith much simpler PSF) is to perform 3 repetitions of the procedure.
Figure 7.2:
Two sample reconstructed star profiles (on the left), used to obtain the finalPSF model (on the right), for about 1300 (top) and 100 (bottom) pixels from the frame centrerespectively. The model reproduces profiles in high detail, although far tails and asymmetricparts are lost. .1.2 Quality of the model The final model recalculated for the real sky data seems to reproduce collected imagesquite well, for slightly asymmetric profiles close to the frame centre and for very deformedones close to the frame edge, as shown in fig. 7.2. However, as in the laboratory mea-surements, far “tails” of deformed PSFs are not properly described. Moreover, one cansee that the reconstructed profile is slightly asymmetric – an effect that was not visibleduring laboratory measurements and is probably induced by the slight tilt of CCD surfacerelative to the lenses axes.As mentioned before, the general shape of the PSF and its development on the CCDis similar for the artificial point source (laboratory measurements) and the real sky data.The similarity of the PSF shape evolution is visible on the plots showing area above of the maximal signal for the final model (fig. 7.3) and the laboratory model (fig. 4.12).The plots cannot be directly compared, for in the former case a mathematical function isintegrated, while in the latter real measurements are analysed. Still, the behaviour is thesame: signal is contained in small number of pixels close to the frame centre and near theframe corner and spread over area even few times larger for intermediate positions. Themaximum position here, however, is oscillating around 1200 pixels from the frame centre,while it was about 1000 pixels for the laboratory measurements. The difference betweenmaxima for different quarters of the CCD is probably due to the sensor tilt describedpreviously.In general, the number of considered star measurements (number of frames) was suffi-cient for profile reconstruction and model fit. However, the number of bright stars in thevery corners of the frame seem to be somewhat scarce compared to the rapid changes ofthe PSF shape in this area. This obstacle cannot be easily solved, for there are simplyvery few stars in the corner of the frame compared to the other parts of the frame dueto the sensor geometry. Stars from more fields should be gathered but this requires extratime, needed for preparing each field, and introduces some additional uncertainties. Ingeneral, the lower density of stars per area in corners of the frame, the lower quality ofthe model for these coordinates.Additional systematic uncertainties in the model are due to the selected method of r [pixel] A [ p i xe l ] Figure 7.3:
Sensor area A covered by the signal larger than of the maximal signal for thefinal model of PSF, applied along the frame diagonals (in distance from the frame centre r). Fourcurves represent fits to the results obtained at different quarters of the CCD.
The final model of the point spread function (introducing 2D Delaunay triangles interpola-tion in polar coordinates) has been used as a basis for simple photometry and astrometry.Each star on a frame has been fitted with a PSF specific for its coordinates, only scale(signal level), background and position on the pixel were treated as free parameters, notgiven by the model. Position on the pixel in x and y coordinates were the only numeri-cally fitted parameters, while the scale was being calculated analytically, according to theformula: S = (cid:88) x,y f ( x, y ) s ( x, y ) + bσ ( x, y ) − b ∗ (cid:88) x,y f ( x, y ) σ ( x, y ) (cid:88) x,y f ( x, y ) σ ( x, y ) (7.1)where x, y stand for pixels coordinates on the image, f ( x, y ) denote the PSF modelfunction value for centres of image pixels (for a given PSF position), b is the backgroundand σ ( x, y ) stands for the uncertainty of the signal in pixel ( x, y ) . We assume here that σ (cid:39) (cid:112) s ( x, y ) + b , where s ( x, y ) is the signal in pixel ( x, y ) . The eq. 7.1 is obtainedfrom an analytic minimization of a χ equation for the model function to real signal andbackground at given PSF position. Although background can be fitted analytically aswell, results turned up to be most stable if we set background as to a trimmed median of all pixels within the area × pixels around the star’s centre.The position on the pixel in the fitting procedure was being initialized with a positionestimated from the transformation of celestial coordinates of the star to frame coordinates,using general astrometry of the frame. However, this transformation procedure does nottake into account different shapes of the PSF and thus different positions of a geometricalcentre of the star’s image on the frame . Thus a correction had to be determined – the The trimmed median in this case was median of values of all pixels, after removing pixel valuesoutside mean ± · RMS . The term “geometrical centre” of the star is probably a little bit misused in this case. In truth, it isthe centre of the modified radial Zernike polynomials set used for obtaining the model PSF. This pointwas the best centre of the star’s shape according to the numerical procedure fitting the model parameters,but its physical meaning had not been determined and is rather not relevant. In the general case it isdifferent from the “centre of mass” position. igure 7.4: Fit area A F giving the smallest brightness measurement spread vs. magnitudo Mfor a test set of stars. systematic shift between position from transformation and fitted position on the pixel allover the sensor.The other important step for photometry and astrometry preparation was a determi-nation of an optimal area around the star centre for which the fitting procedure should beperformed. Photometry have been made on 172 frames for 78 chosen stars and 8 differentfitting areas: 5 pixels (a cross), 9 pixels (a × square), 13 pixels (as for 9, but with4 additional), 21 pixels (a × square without corners), 25 pixels (a × square), 49pixels (a × square), 81 pixels (a × square) and 121 pixels (a × square). Foreach star the area with the smallest spread of the fitted signal as calculated from thephotometry results on the 172 considered frames was chosen. Results of the best fit areavs. magnitudo are visible on fig. 7.4. No straightforward dependency of the fit area onthe star’s brightness is visible on this plot, but it forms an input for an educated guess:for stars brighter than m the 49 pixel area and for dimmer the 25 pixel area have beenchosen as resulting in best results for most stars.It has to be emphasized at this point that the discussed photometry and astrometryprocedure is a very simple one, performed rather to test the PSF model than to competewith professional profile algorithms. The main difference is the fact that all the starswere fitted in a one-stage, straightforward manner, without, for example, subtractingneighbours of each fitted star – a common procedure in many profile algorithms. The detailed photometry studies have been performed on two fields of one night and onefield of a different night. All results were very similar, so only the results for 172 framesof one field are shown here. Results of the fitted magnitudo spread vs magnitudo, asobtained from the polynomial profile photometry and compared to ASAS photometry areshown in fig. 7.5. The smallest obtained brightness spread is about . m for nearly thewhole available brightness range (fig. 7.5, left). Up to m the spread rarely exceeds . m , The magnitudo spread is obtained from the fitted signal, without applying any corrections based onreference stars, due to the correction calculation time and its small relevance for comparison of photometrymethods. M D polynomial photometryASAS photometrypolynomial photometry profileASAS photometry profilepolynomial photometryASAS photometrypolynomial photometry profileASAS photometry profilepolynomial photometryASAS photometrypolynomial photometry profileASAS photometry profilepolynomial photometryASAS photometrypolynomial photometry profileASAS photometry profile M A M D / P M D Figure 7.5:
Measured magnitudo spread ∆M for polynomial (blue points) and ASAS photom-etry (red points) (left) and the ratio of the polynomial to ASAS photometry spread ∆M P ∆M A (right)as a function of the magnitudo M of the measured star. but for dimmer stars it steeply ascends to reach close to . m for ∼ m . The best rangefor polynomial photometry all over the frame is around m − m . For brighter stars fitsare less stable probably due to more evident “tails” of PSF not properly described by themodel. For dimmer stars the reason is simply the descending signal to noise ratio andneighbouring brighter stars disturbing the measurement – the number of these increasingwith magnitudo.The polynomial photometry is comparable to ASAS aperture photometry, similar orslightly worse for most stars up to ∼ . m , and slightly improving with magnitudo (fig. 7.5,right). However, for the most stars the ratio of brightness measurements uncertaintiesremain close to 1 all over the magnitudo range. There are probably several reasons ofthe model giving in general results not better then the aperture photometry. The firstis the possible instability of the fit in the profile photometry, especially when there isanother star in proximity, while the aperture photometry is always stable. The second isthe fact, that ASAS algorithm calculates a mean value of 4 apertures of different sizes.The result of such calculation is shown here, while results of single aperture are in mostcases worse than the polynomial photometry. Furthermore, brightness spread for eachstar was calculated only for those frames accepted by the ASAS photometry, the numberof such frames significantly decreasing with magnitudo, thus results shown favour ASASalgorithm slightly.Figure 7.6 shows comparison of the two photometry methods as a function of theposition on the chip. The discussed ratio of the profile photometry spread to the ASASphotometry spread is smallest for an annulus about 500-1000 pixels around the framecentre (fig. 7.6, left). This covers the range of an extended PSF size, according to fig.7.3. Additionally, the ratio is smallest for the top-right corner of the frame. The reasonof this may be related to the change of the PSF shape with the azimuthal coordinate onthe frame, but higher statistics in corners is needed to check this hypothesis.The polynomial photometry is in most cases better than ASAS photometry for starsabove m (fig. 7.6, right). Based on this results, it can be stated that the polynomialphotometry is better for dim stars far from the frame centre. However, according to theresults shown in fig. 7.5 (right), even for bright stars there is a number of objects for which82 CD x [pixel]
CCD y [ p i xe l ] A M D / P M D r [pixel]0 200 400 600 800 1000 1200 1400 M A M D / P M D Figure 7.6:
Ratio of the polynomial to ASAS photometry magnitudo spread ∆M P ∆M A as a functionof the star position on the chip (left), and as a function of star magnitudo and distance from theframe centre (right). the model performs better. This shows that it is definitely worth trying to fit “objectsof special interest” with the polynomial method – it is possible that the frame-by-framefit in human controlled condition would eventually give a better result than the aperturephotometry. The polynomial astrometry is performed in the same manner and in the same time asphotometry, but the subject of analysis is not the analytically fitted scale (signal) of astar, but its numerically fitted position on the CCD. The position of a star relative tothe pixel edge cannot be accurately determined from a simple centre of mass calculation(as explained in sec. 4.4.1), but the accuracy for bright stars is much improved for
CCD x [pixel]
CCD y [ p i xe l ] polynomial astrometryASAS astrometry (a) Position of a star (b) Distance between two stars Figure 7.7:
Changes in the position of the selected bright star on the CCD during 172 exposures(left) and the distribution of the measured distance between two stars D (right). Compared areresults from polynomial and ASAS astrometry. [pixel] r [ p i xe l ] D polynomial astrometryASAS astrometry r [pixel]0 200 400 600 800 1000 1200 1400 M A r D / P r D Figure 7.8:
Comparison between polynomial and ASAS astrometry results. Left: positionspread ∆r for stars with magnitudo M < m as a function of the distance from the frame centrer. Right: ratio ∆r P ∆r A of position spread of polynomial and ASAS astrometry as a function of starbrightness and distance from the frame centre. more sophisticated algorithms, like ASAS aperture astrometry and discussed polynomialastrometry performed by fitting the modelled PSF to the star. The star position on theCCD of the “Pi of the Sky” camera can be a subject to significant frame-to-frame changes,as seen on fig. 7.7(a), thus the accurate astrometry on a single frame is important. Themovement of a bright star ( ∼ m ) during 172 exposures is clearly visible in results fromboth polynomial and ASAS astrometries and is similar for all the bright stars on theframe. However, slight differences between results of the two algorithms remain and aredirectly related to the differences in their accuracies, becoming more visible for dimmerstars. The question is, how to compare astrometry accuracies?Perhaps the simplest idea would be to analyse the difference between a star positionobtained from the astrometric algorithm and coordinates from a stars catalogue. However,this method involves transformation of the position on the frame to the celestial coordi-nates. This introduces significant uncertainty related to the estimation of the frame’scentre and orientation in celestial coordinates.Nevertheless, all the stars on CCD shift from frame to frame in approximately the samemanner , so the distance between two stars should remain constant on all frames. Thusany spread of measured distance between stars should be a direct result of inaccuraciesin position determination and can be used to quantify astrometry quality.For a given star up to 100 stars in × pixel proximity and ± . m range arefound. For each neighbouring star a histogram of its distance from the given star iscalculated for all the frames out of 172 on which ASAS astrometry returned results forboth stars. The RMS of such a histogram (example presented in fig. 7.7(b)) is usedas the spread of the distance between two stars. The spread is calculated for each pairand the average of all results gives the final estimate of position uncertainty for a givenstar. The number of frames analysed decreases with the magnitudo of a star while thenumber of neighbouring star in ± . m range increases with magnitudo. As a result the The general movement direction and range on the considered frames is similar for all stars, howeverdetails in tracks are reproduced only for stars in proximity. This fact is probably due to the rotation ofthe sensor while following the sky movement. better, being as much as 2.5 times better for bright stars − pixelsfrom the frame centre, the distance uncertainty being ∼ . pixel for the polynomialand ∼ . pixel for the ASAS astrometry (see fig. 7.8, left). The ratio of the estimateduncertainties is the smallest for the part of the frame containing the highest number ofstars, due to geometrical reasons. The quality of the algorithm drops for the distance veryclose and far from the frame centre, where the PSF of a star becomes smaller – probablydue to the fit performing better on larger shapes.Comparison of the two astrometry methods as a function of the distance to framecentre and magnitudo is shown in fig. 7.8 (right). In general, for most of the magnitudo-position space polynomial astrometry performs better than ASAS. The ratio best values,reaching 0.3 for bright stars slightly increases with magnitudo, but the range in the dis-tance from the frame centre where the profile astrometry performs much better becomeswider. The brightness dependency shows, that ASAS astrometry gives better results forstars very close and far from the frame centre, but only for magnitudo ≤ . . More-over, closer examination of these regions reveals that the higher ratio is rather due toinstabilities in the fit than due to systematic differences.Results of the polynomial astrometry show that this algorithm should be used in nearlyall cases for which a longer processing time is acceptable. Such cases include positiondetermination for objects of special interest, for which such a precise measurement isimportant, for example, unidentified artificial satellites of Earth, for which determiningorbital parameters would help in creation of catalogues of “cosmic debris”. Applications of the obtained PSF model are not limited to photometry and astrometry.One of other possible utilizations is a dedicated search for signal at specific coordinateson the frame. An example presented here is the search for optical precursor[62] to “thenaked eye” burst GRB080319B, which, due to its unprecedented brightness perfectly suitsthe task. “Pi of the Sky” started observing the position of the burst more than 19 minutesbefore the trigger, providing required data. The position of the GRB prior to the explosionand on the first 3 exposures was close to the edge of the frame, thus the burst image andits possible precursor were a subject to a large PSF deformation. In this circumstancesmeasurement of a possible precursor signal with the dedicated PSF model should givea larger signal to noise ratio than a standard profile or aperture photometry based oncircular/elliptical PSFs.
The analysis was performed by fitting a model PSF to GRB coordinates on the frame onall the frames covering 19 minutes prior to the explosion, on two cameras (with internalnames: k2a and k2d) of the “Pi of the Sky” prototype. Figure 7.9 (left) shows thefitted signal value II for k2a camera, for all considered frames. To suppress systematicuncertainty fitted signal I is referred to the scale I of a nearby . m star. No signal85 [minutes] t-t-20 -15 -10 -5 0 I/I -0.02-0.0100.010.020.030.040.050.06
95% CL s M I/I-0.04 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.040123456789
GRB positionControl position
Figure 7.9:
Left: the measured signal II in the expected place of the optical precursor to theGRB080319B as a function of time before the GRB t − t and the corresponding limits on itsbrightness for k2a camera. Right: distribution of the measured signal for 84 frames before theburst with fit of normal distribution indicated. exceeding σ limit has been found on k2a camera. Standard approach in case of “nosignal”, at least in astronomical observations, is to quote σ limit assuming measuredvalue of signal to be zero. Estimation of all the limits is based on an assumption thatthe measurement error for a small signal that would be emitted by the optical precursoris similar to the fit error for the sky background at this position. The uncertainty σ wascalculated from a fit of the normal distribution to the histogram of signal ratio II (shownin fig. 7.9, right). The result is consistent with distributions of signal in 3 empty controlareas in the burst proximity which were considered as a cross-check. In this approach theresulting σ limiting magnitudo for k2a camera in the polynomial photometry is . m .The published limit given by the standard “Pi of the Sky” photometry was ∼ . m , soa ∼ . m increase in limit value was obtained by using the new method. On the confidence level ( = 1 . σ ) the limiting magnitudo for k2a camera is . m .The advantage of the polynomial photometry is that it also allows setting limits onparticular frames based on the signal level in the GRB coordinates (or rather, in thiscase, background fluctuation in the place where we seek for the signal). The limit wasnumerically calculated according to the formula given in [63], which allows extractinglimits with well defined confidence level (CL) all over the measured signal range, includingnegative signal fluctuations. The σ limits ( . CL) calculated on single frames of k2acamera range from . m to . m , fluctuating in most cases between . m and . m .The confidence level limits (which are commonly used in high energy physics andseem more appropriate for single frame analysis) range from . m to . m , fluctuatingin most cases between . m and m . The σ limit on the last frame before the GRBis . m , and this value can be compared to the limit set by standard “Pi of the Sky”photometry of . m (which is actually the limit based on the same frame, but withassumption of no signal). In this case we obtain . m limit improvement when usingpolynomial photometry. In this case not the ASAS photometry, but so called “fast Pi photometry”, a slightly different aperturebased algorithm, was used. [minutes] t-t-20 -15 -10 -5 0 I/I -0.02-0.0100.010.020.030.040.050.06
95% CL s M I/I-0.03 -0.02 -0.01 0 0.01 0.02 0.03012345678
Figure 7.10:
Left: the measured signal II in the expected place of the optical precursor tothe GRB080319B as a function of time before the GRB t − t and the corresponding limits onits brightness combined for two – k2a and k2d – cameras. Right: distribution of the combinedsignal for 84 frames before the burst with fit of normal distribution indicated. The polynomial photometry, measuring the signal in the expected position of the precur-sor, allows also combining measurements from different cameras, increasing the signal tonoise ratio and thus giving a possibility for deeper search or setting better limits. Priorto the comparison, the signal from k2a had to be normalized to the level of the signal onk2d camera, due to differences in amplification and other hardware related factors. Thenormalization constant was taken from the comparison of the signal between the star ofreference on two cameras. After the normalization, noise levels (readout errors) for eachcamera have been calculated as the standard deviation of normal distribution fitted tosingle-camera II histograms. The combined scale for precursor I was then computed asan weighted average of scales fitted on both cameras. Similar procedure (but taking intoaccount the star signal error) was used for the scale I of the star of reference.No signal above σ level is visible in the combined signal distribution, as shown infig. 7.10, (left). However, most limits greatly improve. Standard σ limits calculatedassuming zero signal is m , same as limiting magnitude given by standard “Pi of the Sky”photometry for 20 coadded frames and higher by . m than the single camera limit. The CL limit also increases by . m reaching . m . The σ limit based on measuredsignal for single frames is between . m and . m for most frames. For CL limit iscontained between . m and . m , and for most frames between . m and . m . The σ limiting magnitudo on the last frame before the GRB is m , this time backgroundfluctuation into low values in k2a camera is compensated by the opposite fluctuation ink2d camera, resulting in signal consistent with zero.In general, the combined limits based on signal measured in two cameras (fig. 7.10,right) are more stable than limits on the single, k2a camera. The measurement error wasreduced by nearly , from . for k2a to . for combined signal. This is the mainreason for limit improvement with two cameras.The example of the precursor search shows that the polynomial model of PSF givesbetter chances in dedicated searches for a small signal in specific coordinates. In case ofno signal found, it allows us to use more advanced limit setting procedures, giving betterresults. Moreover, it allows combining the signal from multiple cameras of “Pi of the Sky”87bserving the same place, which provided important result in this case and may lead toeven better outcomes in the final “Pi of the Sky” system, where 4 cameras will observethe same part of the sky (in “deep” mode of work).88 hapter 8Simulator The polynomial model of point spread function described in this work can be used inmultiple applications, as described in chap. 7. However, the main idea behind its de-velopment was to gain a better understanding of the “Pi of the Sky” detector. Duringthe process of the model development we found out how the star image (PSF) in ourcameras is formed (in general) and described it in a parametric way along with spacialand energetic pixel response functions. The obtained knowledge allows reconstruction ofa “Pi of the Sky” frame in high detail. Good description of real data opens the possibilityto extend the study and to test frame processing as well as dedicated analysis algorithmsin controllable circumstances using simulated data samples. This approach, often referredas Monte Carlo method, is widely used in high energy physics, to test and develop newdetectors, predict future results and analyse those already obtained. Inspired by meth-ods from experimental particle physics, dedicated simulator code has been created, itsfeatures, tests and possible applications are presented in this chapter.
The main aim of the simulator is to generate a frame, which would correspond to a framemeasured by “Pi of the Sky” camera for selected sky coordinates. This task is performedin following steps:
Sky background generation
All pixels on the frame are initialized to a value of sky background, which is one ofthe simulator parameters. Star distribution reproduction
The frame range on the sky is computed from position of its centre in celestialcoordinates. All stars positions and brightnesses in this range are read from TYCHOcatalogue.
Frame coordinates fluctuations
If requested, star positions on the CCD surface are shifted by a value randomlyselected from an assumed distribution of pointing fluctuations (flat distributionin a single pixel range is used as default); one shift vector is generated for the Lenses transmission is taken into account by multiplying the background by a transmission functionreconstructed from the real flat-field frames. If requested, flat-divided frame (i. e. without correction)can be generated. igure 8.1: Comparison between simulated (left) and real (right) images of the sky. Somenon-star objects are visible on the real frame but not included in the simulation, such as a verybright, saturated satellite flash visible near the lower edge of the frame. whole frame. This reproduces slight changes in frames coordinates due to mountmovement. Obtained star positions form an input for each star image creation (fig.8.1).
Star brightness fluctuations
If requested, each star’s brightness is multiplied by a value randomly selected froman assumed distribution (normal distribution centered in 1 and with standarddeviation was used in results presented further in this chapter). This step reproducessource brightness fluctuations due to, for example, turbulence of the atmosphere. Star image reproduction
An image of each star is computed from the model PSF for given coordinates on theframe. After normalizing to the star magnitudo, taking into account fluctuations,if included, it is projected on the frame.
Poisson statistics application
Each pixel value is changed according to Poisson distribution, which describes fluc-tuations in number of electrons collected by CCD in given period of time. Thistakes into account fluctuations of both incoming photon flux and CCD quantumefficiency.
Readout noise
CCD readout noise is added to each pixel. The noise is modelled by a Gaussdistribution. The default noise level (sigma of the distribution) is 16e – a valuecharacterizing simulated camera used in laboratory measurements.
Gain fluctuations
Collected charge is converted to the output signal units (Analog-Digital converterUnits: ADU) based on the assumed CCD gain (default is 2.7e/ADU). In additionfluctuations of the CCD gain, approximated by a Gauss distribution, can be appliedto each pixel. 90 igure 8.2:
Comparison between zoomed parts of simulated (left) and real (right) frames. Thecore of PSF is similar, while the real data show a second order “halo” effect for the brighteststars, not included in the model and thus the simulation. The bright flash visible near the middleof the right border of the real frame part is due to an artificial satellite.
Energetic pixel response function application
Signal measured in all pixels is recalculated according to the energetic pixel responsefunction, which introduces non-linearity of generated signal for bright light sources.After the requested number of frames is generated, the result are stored in FITS files –a standard format for astronomical data. This output is suitable for all software workingon FITS frames including most of the “Pi of the Sky” analysis algorithms. Additionally,results can be stored in the ROOT framework[64] TTree classes, which were the basicdata format for analysis performed in this thesis.From the point of view of external programs, generated frames are identical to real skydata. Thus the simulation of a real part of the sky as seen by the “Pi of the Sky” camerais a good tool for testing photometric, astrometric, sources detection and cataloguingalgorithms as well as other analysis programs in well defined conditions. Moreover, thePSF model can be easily changed to describe different apparatus.
Fig. 8.1 shows the simulated and the real frame of the same part of the sky. The spatialdistribution of stars is the same, although the real frame reveals some objects not includedin the catalogue of stars, such as satellites or planets. Clearly visible in 8.1 is a flash ofone of the Earth’s artificial satellites: the brightest point near the lower edge of the framecreating a bright belt of charge going up to the upper edge . The background on the realframe is not as uniform as on the simulated frame due to true background variations oran imperfect flat frame used in the reduction process. Nevertheless, function describingspatial distribution of the background can be easily added to the simulator. The belt o charge is a result of so-called “blooming”. Charge from saturated pixels floods neighbouringpixels and leaves a trail in the most brightest pixels’ columns during the frame readout process. I/I D real starssimulated starsreal stars profilesimulated stars profile M I/I D -1 real starssimulated starsreal stars profilesimulated stars profile Figure 8.3:
Polynomial photometry results for simulated and real sky images. Dependence ofthe relative brightness measurement error ∆II on the star magnitude M is shown in linear (left)and logarithmic (right) scale.
The differences between apparent brightnesses of bright stars on the real frame arebigger than on the simulated frame. This is an illusion caused not by star’s image realbrightness but its size. As mentioned in chap. 7, far “tails” of the PSF are not reproducedby the model. Additionally, a kind of “halo” is visible around the brightest stars (fig.8.2), which is a second-order effect, not taken into account in the model. The “halo” isthe main cause behind the brightest stars appearing larger and thus brighter on the realframe compared to the simulated frame. However, the magnification of the frame revealsthat dimmer stars are reproduced in high detail.It has to be noted, that the results of the visual comparison of the real and thesimulated frame, although good, are rather of a lower importance. What matters mostis the degree of similarity between results of analysis obtained on real and simulateddata. Therefore a polynomial photometry has been performed on a series of 25 simulatedcatalogue frames.The dependence of the stars signal spread on the brightness is very similar for realand simulated data, as shown in fig. 8.3. The slightly better match of results could beobtained by fine-tuning of the simulation parameters describing measurement fluctuations.It is visible that the number of poorly-measured stars is similar all over the magnitudorange. Additionally, other details of the dependence are well reproduced, such as the highinstability of the brightest stars, caused probably by the CCD energetic response non-linearity, and the subgroup of dim stars with low brightness measurement uncertaintybetween ∼ . m and ∼ . m .The spatial distribution of the photometry uncertainty on the frame is also very similarfor real and simulated data, as shown in fig. 8.4. The similarity proves that the PSF isthe main factor responsible for varying photometry quality on the frame.The key to the good agreement between results of photometry performed on generatedimages and real data is not only the PSF, but also proper simulation of the frame variabil-ity – the fluctuations of stars brightness caused by miscellaneous factors. The describedsimulation takes into account 5 factors contributing to the frame-to-frame variability. Theoverall signal for each star (mainly due to atmospheric turbulence) is modified by gaussianfluctuations with σ = 1% , and the star position is shifted by a random value (uniform92 [pixel] M I/I D (a) simulated data r [pixel] M I/I D (b) real data Figure 8.4:
Dependence of the relative signal measurement uncertainty ∆II on the star mag-nitudo M and the distance from the frame centre r, for polynomial photometry performed onsimulated sky images (left) and on real data (right). inside a pixel). Charge deposited in each pixel of the frame is smeared according to thePoisson distribution, readout noise of 32 electrons and gaussian gain fluctuations of . are also applied .The impact of each variability factor on the photometry quality is shown in fig. 8.5.The biggest uncertainty is caused by gain fluctuations, which dominate in nearly the wholeregion except for the bright stars below . m , where energetic pixel response function non-linearity becomes significant. In this region the most important factor is the uncertaintycaused by position fluctuations. For stars dimmer than ∼ m Poisson fluctuations becomemore significant than position fluctuations. Star fluctuations (due to turbulence) add
M4 5 6 7 8 9 10 11 12
I/I D photon statisticsgain fluctuationsreadout noiseposition fluctuationsstar turbulationscombined fluctuationsdata M4 5 6 7 8 9 10 11 12
I/I D -2 -1 photon statisticsgain fluctuationsreadout noiseposition fluctuationsstar turbulationscombined fluctuationsdata Figure 8.5:
Relative error of the photometry performed on the simulated data, as a function ofthe star brightness (points) and the estimated contribution from different sources of variability(colour lines). The error is shown in linear (left plot) and in logarithmic scale (right plot). The readout noise for the described camera as measured in laboratory conditions was 16 electrons,but the value of 32 electrons seems to describe real data better. This can be attributed to the contributionfrom the dark frame subtraction process, as well as to aging of the electronics.
93 factor constant (compared to the other factors) nearly all over the magnitudo rangeand were introduced here mainly for the simulator functionality presentation purposes.The readout noise is of a small importance over the whole magnitudo range. Combinedcontributions from the separate factors follow the shape of the data coming from thesimulation with all the factors included.All the variability parameters, except the photon statistics, were an educated guessbased on the photometry of real data. However, the proportion between factors should notchange significantly after tuning. The important conclusion coming from the simulationis that the biggest improvement in photometry stability could be obtained by reducingCCD gain fluctuations. Still, the uncertainty caused by position fluctuations is significantand should be a subject of further studies. Reducing readout noise of the detector is of arather minor importance.
Generation of series of frames consisting of catalogue stars is a good tool for generalsoftware and system tests, and for planning future measurements. It can be used forpredicting cameras range in different circumstances, estimating photometry/astrometryefficiency in different parts of the frame and for different star’s magnitudo, etc. Electronicsnoises and mechanical tracking precision can be easily adjusted and thus those calculationscan be easily repeated to study the all possible impacts of system electronics and/ormechanics upgrade.Although this mode of operation allows general testing of algorithms, their deeperanalysis and tuning can be performed in a much more efficient way. The key is a generationof a frame dedicated for such tasks, consisting of a small number of stars with specificparameters, like coordinates, brightness, etc., and their changes in time. These featurescreate two possible groups of additional simulator uses.
Performing analysis on a frame of “randomly” scattered stars has a particular drawback– the number of stars changes with the distance from the frame centre and is especiallysmall in the corners of the frame. Thus corners of the frame always suffer from a smallnumber of measurements. The simulator capacity of generating sets of stars for requestedpositions allows for a detailed study of algorithms performance on any position on theframe.However, generating stars in requested positions aids studying more sophisticatedtasks. The first among them is the study of the possibility of using parallax phenomenonin the “Pi of the Sky” project. Frames reproducing parts of the sky as seen from differentsites on the Earth can be generated, with an image of a near-Earth object appearing indifferent places. Images shown in fig. 8.6 were generated assuming altitude of 23000 kmfor a satellite directly above “Pi of the Sky” sites distant by 100 km, resulting in ∼ (cid:48) angular separation of the images.In this case simulated data can help determining the maximal distance of the objectfrom Earth for which it could be stated, that it appears in two cameras at differentcoordinates. Assuming that the smallest measurable position difference of 0.1 pixel ( . (cid:48) )allows identifying near-Earth objects maximally ∼ km from Earth surface forsites distant by 100 km and ∼ km for sites distant by 30 km, in both cases behind94 a) site 1 frame (b) site 2 frame Figure 8.6:
Simulation of the parallax phenomenon: satellite 23000 km above Earth as seenfrom two “Pi of the Sky” sites distant by 100 km. the Moon orbit. However, the real separation distance changes (mainly due to the PSFchanges) with coordinates on the frame and can be calculated using the simulated data.A similar task is studying efficiency of analysis of crowded fields. Simulated images ofoverlapping stars (fig. 8.7) of different magnitudos can be used for determining when thestars become indistinguishable, how their position/brightness measurements depends ontheir distance, magnitudo, position on the frame, etc. Additionally, the simulator allowsspecifying a function describing a star movement through the frame during specified time.This capacity mimics mount tracking errors and can help in compensating stars variabilitycaused by it.
In the project such as “Pi of the Sky”, dedicated to observing rapidly variable phenom-ena on the sky, a special stress is put on the algorithms analysing stars variability anddetecting new, explosive sources, such as gamma ray bursts. So far optimal parame-ters and efficiency of the flash detection algorithm have been estimated in a very simplesimulation[40]. Similar estimations have been performed for other detection and variabil-ity analysis software used in the project.Also in this case, the simulator aids much more sophisticated analysis of such al-gorithms. It allows creation of star-like objects with specified lightcurve, reproducingrequested variability, as shown in fig. 8.8. Thus for periodic variable stars a period andamplitude determination algorithms can be tested. For flare, novae stars or optical tran-sients thresholds for detection can be estimated, etc. Additionally, all other features of thesimulator can be used simultaneously, like generation of a neighbouring star with specificparameters, mount tracking instabilities, electronic noises, etc., giving testing possibilitiesunprecedented in previous simulations. 95 igure 8.7:
Example image from the simulation of a dense field – three overlapping stars of m , m and m . The design of the simulator allows easy creation of extensions, such as introducing newcatalogues of objects like planets or stars, which would help replicating the real sky viewin more detail. Simulating near-earth objects could help develop algorithms for orbitsdetermination. Additionally, more features related to the electronics performance can beadded, such as blooming (signal saturation) effect, hot pixels, nongaussian distribution ofnoises, etc., which in some conditions affect the data analysis.Presented features of the simulator make it a very powerful tool for data analysis al-
Frame number0 20 40 60 80 100 II - -0.1-0.08-0.06-0.04-0.0200.020.040.06 aperture 0aperture 1aperture 2aperture 3aperture 4 Figure 8.8:
ASAS photometry results (signal I shifted by mean signal ¯I to cover the same range)for a simulation of a sinusoidally-variable star of ∼ . m . Differences between amplitudes givenby apertures of different sizes are probably caused by a non-linearity of energetic pixel responsefunction. hapter 9Summary “Pi of the Sky” is a set of robotic telescopes aiming for continuous observation of a largepart of the sky with high temporal and optical resolution using wide field-of-view instru-ments. Its primary goal is to look for optical afterglows associated with the gamma raybursts (GRB), but it is also well suited to study any kind of short time scale astrophysicalphenomena. The very-wide sky coverage of the final “Pi of the Sky” system, about 1.5 sr,requires large field of view of a single camera, ◦ × ◦ . Lenses allowing such cover-age introduce large deformations of the point spread function. This results in increaseduncertainties of measured brightness and positions of the observed objects.The main goal of the presented study was to develop a detailed model of the cam-era response including a parametrization of the star image shape (point spread function– PSF). For this purpose the author of this thesis prepared and configured a uniquelaboratory setup, which allowed precise measurements of the detector response in fullycontrolled conditions. The point spread function, along with the pixel sensitivity functionas well as spatial and energetic pixel response functions, for multiple wavelengths were re-constructed. These results were a starting point for development of the diffraction-basedmodel of the PSF.The diffraction-based model describes the observed PSF shape as a result of the lightpropagation through camera’s lenses. In the considered approximation the wavefrontaberrations were parametrized as a sum of the so-called Zernike polynomials. In thisapproach PSF parametrization required integration of the wavefront, which is a rapidlyoscillating function, over the aperture. Therefore calculations could not be made effi-ciently with standard numerical algorithms. Non-standard method of integration wassuccessfully implemented along with numerous numerical simplifications. The diffractivemodel reproduced measured profiles quantitatively, proving that the approach is proper ingeneral. However, numerical calculations involved are still very time consuming and notprecise enough to allow for practical application of this method for modelling or simulationpurposes.As an alternative, an effective model parametrizing the PSF in the image plane hasbeen developed based on laboratory measurements. This required constructing a dedi-cated mathematical base, derived from Zernike polynomials, then selecting proper modelparameters and parameter interpolation methods. After adjusting the model to describethe real sky data from the prototype, it was incorporated into photometric and astromet-ric algorithms and tested on collection of real sky images. No improvement in photometrycompared to standard aperture based algorithm was obtained. However, significant im-provement of more than a factor of two (for bright stars) has been reached in astrometryaccuracy. Additionally, new limits on the optical precursor to “the naked-eye burst” have98een set, based on the signal level reconstructed on each frame in the place of GRBand combined from two cameras. This approach was not possible prior to the modeldevelopment.The improvement in astrometry convinced us that the effective model describes dataproperly. Therefore, the model was used as an input for a general purpose package pre-pared for Monte-Carlo simulation of sky images as seen by the “Pi of the Sky” camera.The simulator can be used to generate sky frame for given celestial coordinates properlyreproducing the deformed stars shapes. It is also capable of creating frames of specialinterest, consisting of objects with specified parameters (eg. GRB flashes or satellites),which can then be used for testing analysis algorithms. The description of measurementfluctuations introduced in the frame generation was successful at reproducing the quanti-tative behaviour of photometry errors in the real data, thus proving the simulator a wellsuited tool for development of data analysis algorithms.The development of the CCD camera response model, which was the subject of thisthesis, turned out to be a very complicated task. At many stages of the analysis dif-ferent approaches were possible. Because of time limitation, not all the possibilitieswere explored in details. Nevertheless some very important results have been reached,mainly a more than significant improvement in astrometry precision and new limits on theGRB080319B optical precursor, improving previous limits by . m . The results of thesimulation indicated the most significant factor responsible for stars brightness measure-ment errors – the gain fluctuations of the camera’s electronics. Up to now, this factor hasbeen underestimated and this new information may help to improve future “Pi of the Sky”cameras significantly.It is clear that further improvements are possible both in the PSF modelling and inthe frame generator. A valuable achievement would be to complete the diffractive modelof the light passage through lenses. Although this approach would be, in any case, to slowto be incorporated into the data analysis algorithms, it could be used as a reference fortuning of the effective model. Additionally, it could allow for predicting the PSF changeswith spectrum of light and different focusing of cameras. However, the complexity of thistask calls for a separate, dedicated research.We hope that the detailed PSF description will be useful in analysis of the datafrom the final “Pi of the Sky” system. The first mount of the full system shows bettertracking precision and reduced mechanical vibrations. This should result in reduced shapefluctuations of real star images and allow for much more precise measurements. Reducedreadout noise and gain fluctuations of new cameras should also increase the impact ofproperly parametrized star shape on objects brightness and position determination. Wedo hope that the brightness measurement accuracy, much better than in the prototypesetup can be obtained. First results from the new detector in Spain should be availablesoon.We expect that the developed model, after possible improvements, will be used for anumerous new studies not considered in this work. Precise astrometry may result in moreaccurate summing of star images from consecutive frames, extending the “Pi of the Sky”range. Simultaneously, profile photometry of dense fields could be now much more precise.In both cases, tests and optimization of the analysis algorithms can be performed withthe developed simulator.The full “Pi of the Sky” system, with its large sky coverage, will have the unprece-dented capability of detecting short optical transients on the sky. Moreover, undoubtfulidentification of flashes of an astrophysical origin, thanks to parallax. This is yet an un-touched field in optical astronomy and as such it may lead to many surprises, such as not99et discovered phenomena.Together with still the biggest optical sky coverage, “Pi of the Sky” gives hope for largediscoveries. However, these detections will require studies of the parallax determination,as well as, assuming “Pi of the Sky” is going to be the source of optical transients alertsdistributed to other observatories, most precise measurements of brightness and position.In all these, the developed simulator package should be of great help. Therefore, weare sure that the results presented in this thesis will be helpful for future studies. Full“Pi of the Sky” system with its large optical sky coverage and parallax gives us hope forimportant discoveries in the coming years. 100 ppendix A Number of plots included in the main text of this thesis was reduced for clarity. Shownin this appendix are additional plots complementing results presented in chapter 6: • dependence of all coefficients of Zernike polynomial on the star position, for the firstparameter set, is shown in fig. A.1 (refer fig. 6.2) • same dependence for the second parameter set is shown in fig. A.2 (refer fig. 6.3) • comparison of measured PSFs for colour and white diodes with their model descrip-tion, for different focus settings and distances from the frame centre are shown infig. A.3 (refer fig. 6.8) 101 [pixel]
200 400 600 800 1000 1200 1400 p r [pixel]
200 400 600 800 1000 1200 1400 p r [pixel]
200 400 600 800 1000 1200 1400 p r [pixel]
200 400 600 800 1000 1200 1400 p -50510152025 r [pixel]
200 400 600 800 1000 1200 1400 p r [pixel]
200 400 600 800 1000 1200 1400 p Figure A.1:
Coefficients of Zernike polynomial as a function of the distance from the framecentre r – results of the best parameters searching and fitting procedure. [pixel]
200 400 600 800 1000 1200 1400 p r [pixel]
200 400 600 800 1000 1200 1400 p r [pixel]
200 400 600 800 1000 1200 1400 p r [pixel]
200 400 600 800 1000 1200 1400 p r [pixel]
200 400 600 800 1000 1200 1400 p -2024681012 r [pixel]
200 400 600 800 1000 1200 1400 p -1.6-1.4-1.2-1-0.8-0.6-0.4-0.200.20.4 Figure A.1:
Continuation... [pixel]
200 400 600 800 1000 1200 1400 p -1012345 Figure A.1:
Continuation... [pixel] p r [pixel] p -12-10-8-6-4-202 r [pixel] p -0.100.10.20.30.40.5 r [pixel] p -0.500.511.52 r [pixel] p r [pixel] p -202468 Figure A.2:
Coefficients of Zernike polynomial as a function of the distance from the framecentre r for the finally chosen basis of 17 polynomial terms – 3 circular and 14 giving best χ forall measured profiles. [pixel] p -1-0.8-0.6-0.4-0.200.20.4 r [pixel] p -0.8-0.6-0.4-0.200.2 r [pixel] p -0.6-0.4-0.200.20.40.60.81 r [pixel] p -1.5-1-0.500.511.52 r [pixel] p -3-2-10123 r [pixel] p Figure A.2:
Continuation... [pixel] p -0.8-0.6-0.4-0.200.2 r [pixel] p -0.3-0.2-0.100.10.20.30.4 r [pixel] p -0.100.10.20.30.40.50.6 r [pixel] p -0.6-0.4-0.200.20.40.60.811.2 r [pixel] p -0.6-0.4-0.200.20.40.60.81 Figure A.2:
Continuation... a) Blue diode with nominal focus f s = 1 . m(b) Red diode with nominal focus f s = 1 . m Figure A.3:
Measured PSFs for colour diodes and for the white diode with changed focusing(first rows) compared with their models obtained with the developed parametrization (secondrows), for distance of 1000 (left column), 1200 (centre), 1400 (right column) pixels from theframe centre. Plots for 800 pixels from the frame centre were shown on fig. 6.8. a) White diode with focus changed to f s = 1 . m(b) White diode with focus changed to f s = 1 . m Figure A.3:
Continued... ibliography [1] Z. Fu-Yuan, R. G. Strom, and J. Shi-Yang, “The guest star of AD185 must havebeen a supernova,”
Chinese Journal of Astronomy and Astrophysics , vol. 6, no. 5,pp. 635–640, 2006.[2] A. Hellemans and B. H. Bunch,
The timetables of science: A chronology of the mostimportant people and events in the history of science . Simon and Schuster, 1988.[3] B. R. Goldstein, “Evidence for a supernova, of A.D. 1006,”
Astronomical Journal ,vol. 70, pp. 105–+, Feb. 1965.[4] G. W. Collins, II, W. P. Claspy, and J. C. Martin, “A Reinterpretation of HistoricalReferences to the Supernova of A.D. 1054,”
The Publications of the AstronomicalSociety of the Pacific , vol. 111, pp. 871–880, July 1999.[5] O. Krause, M. Tanaka, T. Usuda, T. Hattori, M. Goto, S. Birkmann, and K. Nomoto,“Tycho brahe’s 1572 supernova as a standard type ia as revealed by its light-echospectrum.,”
Nature , vol. 456, no. 7222, pp. 617–9, 2008.[6] W. P. Blair, K. S. Long, and O. Vancura, “A detailed optical study of Kepler’ssupernova remnant,”
Astrophysical Journal , vol. 366, pp. 484–494, Jan. 1991.[7]
The Cambridge History of Science, Volume 3: Early Modern Science . CambridgeUniversity Press, 2006.[8] B. Zhang and P. Meszaros, “Gamma-Ray Bursts: Progress, Problems & Prospects,”
Int. J. Mod. Phys. , vol. A19, pp. 2385–2472, 2004.[9] R. W. Klebesadel, I. B. Strong, and R. A. Olson, “Observations of Gamma-RayBursts of Cosmic Origin,”
Astrophysical Journal , vol. 182, p. L85, 06 1973.[10] “The third interplanetery network.” .[11] G. Schilling,
Flash!: The Hunt for the Biggest Explosions in the Universe . CambridgeUniversity Press, 2002.[12] “BATSE.” .[13] “BeppoSAX.” .[14] “W. M. Keck Observatory.” http://keckobservatory.org .[15] G. Greco, G. Beskin, S. Karpov, S. Bondar, C. Bartolini, A. Guarnieri, and A. Pic-cioni, “High-Speed and Wide-Field Photometry with TORTORA,”
Advances in As-tronomy , vol. 2010, 2010. 11316] “Raptor.” .[17] “SWIFT.” http://heasarc.nasa.gov/docs/swift/swiftsc.html .[18] J. E. Carson, “GLAST: physics goals and instrument status,”
J. Phys. Conf. Ser. ,vol. 60, pp. 115–118, 2007.[19] V. Cocco et al. , “The science of AGILE: part I,”
Nuclear Physics B - ProceedingsSupplements , vol. 113, no. 1-3, pp. 231 – 238, 2002.[20] K. Mitsuda et al. , “The X-Ray Observatory Suzaku,”
Publications of the AstronomicalSociety of Japan , vol. 59, pp. 1–7, Jan. 2007.[21] R. L. Aptekar et al. , “Konus-W Gamma-Ray Burst Experiment for the GGS WindSpacecraft,”
Space Science Reviews , vol. 71, pp. 265–272, Feb. 1995.[22] B. J. Teegarden and S. J. Sturner, “INTEGRAL Observations of Gamma-RayBursts,” in
AAS/High Energy Astrophysics Division , vol. 31 of
Bulletin of theAmerican Astronomical Society , pp. 717–+, Apr. 1999.[23] W. S. Paciesas et al. , “The Fourth BATSE Gamma-Ray Burst Catalog (Revised),”
Astrophys. J. Suppl. , vol. 122, pp. 465–495, 1999.[24] G. Vedrenne and J. L. Atteia,
Gamma-Ray Bursts The Brightest Explosions in theUniverse . Springer Praxis Books, 2009.[25] P. Meszaros, “The Fireball Shock Model of Gamma Ray Bursts,”
AIP Conf. Proc. ,vol. 526, pp. 514–518, 2000.[26] J. L. Racusin et al. , “Broadband observations of the naked-eye γ -ray burstGRB 080319B,” Nature , vol. 455, no. 7210, pp. 183–188.[27] A. Panaitescu, “Decay phases of Swift X-ray afterglows and the forward-shockmodel,”
Philos Transact A Math Phys Eng Sci. , vol. 365, no. 1854, pp. 1197–1205,2007.[28] P. Meszaros and M. J. Rees, “Optical and Long-Wavelength Afterglow from Gamma-Ray Bursts,”
Astrophysical Journal , vol. 476, p. 232, Feb. 1997.[29] A. D. Falcone et al. , “Observations of X-ray Flares from Gamma Ray Bursts,”
AIPConf. Proc. , vol. 1000, pp. 91–96, 2008.[30] A. I. MacFadyen and S. E. Woosley, “Collapsars: Gamma-Ray Bursts and Explosionsin “Failed Supernovae”,”
Astrophysical Journal , vol. 524, pp. 262–289, Oct. 1999.[31] B. Paczyński and P. Haensel, “Gamma-ray bursts from quark stars,”
Mon. Not. Roy.Astron. Soc. Lett. , vol. 362, pp. L4–L7, 2005.[32] E. Nakar, “Short-hard gamma-ray bursts,”
Phys. Rept. , vol. 442, pp. 166–236, Apr.2007.[33] “Fermi detects ’shocking’ surprise from supernova’s little cousin.” . 11434] W. Baade and F. Zwicky, “On Super-Novae,”
Proceedings of the National Academyof Sciences of the United States of America , vol. 20, pp. 254–259, 1934.[35] “GRB Coordinates Network.” http://gcn.gsfc.nasa.gov .[36] G. Pojmański, “The All Sky Automated Survey. Catalog of Variable Stars. I. h − h Quarter of the Southern Hemisphere,”
Acta Astronomica , vol. 52, no. 397, 2002.[37] “Instituto Nacional De Tecnica Areoespacial.” .[38] “Consejo Superior de Investigaciones Cientificas.” .[39] “Las Campanas Observatory.” .[40] L. W. Piotrowski, “Poszukiwanie błyskow światła widzialnego towarzyszących poza-galaktycznym błyskom gamma,” Master’s thesis, University of Warsaw, May 2005.[41] M. Sokołowski et al. , “Automated detection of short optical transients of astrophysicalorigin in real time,”
Advances in Astronomy , vol. 2010, no. 2010.[42] M. Sokołowski,
Investigation of astrophysical phenomena in shor timescale with“Pi of the sky” apparatus . PhD thesis, Soltan’s Institute for Nuclear Studies, 2007.[43] “Nova V679 Carinae discovery.” http://grb.fuw.edu.pl/pi/var/catac/nova081125.htm .[44] “American Association of Variable Stars Observers (AAVSO).” .[45] “Smarts – small and moderate aperture research telescope system.” .[46] “1RXS J023238.8-371812 outburst detection.” http://grb.fuw.edu.pl/pi/var/catac/11922153.htm .[47] “Gunma astronomical observatory.” .[48] “GSC2.3 S55U020591 outburst detection.” http://grb.fuw.edu.pl/pi/var/catac/14839438.htm .[49] “Flare star CN Leo discovery.” http://grb.fuw.edu.pl/pi/events/200504/01/Frame00445/flare.html .[50] “Flare star GJ3331A discovery.” http://grb.fuw.edu.pl/pi/events/flare/GJ3331A .[51] “Flare star USNO1050.00569325 explosion.” http://grb.fuw.edu.pl/pi/events/200807/18/event0.php .[52] A. Majczyna et al. , “The catalog of short periods stars from the “Pi of the sky” data,”
New Astronomy , vol. 13, pp. 414–417, August 2008.[53] A. Majczyna et al. , “ “Pi of the sky” catalogue of the variable stars from 2006 - 2007data,”
Proceedings of SPIE , vol. 7745.11554] H. Toyozumi and M. C. B. Ashley, “Intra-pixel sensitivity variation and charge trans-fer inefficiency – results of CCD scans,”
Publications of the Astronomical Society ofAustralia , vol. 22, pp. 257–266, 2005.[55] J. J. Olivero and R. L. Longbothum, “Empirical fits to the Voigt line width: A briefreview,”
Journal of Quantitative Spectroscopy and Radiative Transfer , vol. 17, no. 2,pp. 233 – 236, 1977.[56] P. L. von Seidel
Astronomische Nachrichten , p. 289, 1856.[57] G. Conforti, “Zernike aberration coefficients from seidel and higher-order power-seriescoefficients,”
Optics Letters , vol. 8, pp. 407–408.[58] “Seidel vs. zernike.” .[59] A. Miks, J. Novak, and P. Novak, “Method for numerical integration of rapidly oscil-lating functions in diffraction theory,”
International Journal for Numerical Methodsin Engineering , vol. 82, pp. 525–536.[60] B. Delaunay, “Sur la sphere vide,”
Izvestia Akademii Nauk SSSR , pp. 793–800, 1934.[61] D. Shepard, “A two-dimensional interpolation function for irregularly-spaced data,”
Proceedings of the 1968 ACM National Conference , pp. 517–524, 1968.[62] B. Paczyński, “Optical Flashes Preceding GRBs,” 2001.[63] G. Feldman and R. Cousins, “Unified approach to the classical statistical analysis ofsmall signals,”
Physical Review Letters D , vol. 57, pp. 3873–3889.[64] “Root framework.” http://root.cern.chhttp://root.cern.ch