UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

White dwarf populations in globular clusters Goldsbury, Ryan 2016

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata

Download

Media
24-ubc_2016_may_goldsbury_ryan.pdf [ 13.51MB ]
Metadata
JSON: 24-1.0228563.json
JSON-LD: 24-1.0228563-ld.json
RDF/XML (Pretty): 24-1.0228563-rdf.xml
RDF/JSON: 24-1.0228563-rdf.json
Turtle: 24-1.0228563-turtle.txt
N-Triples: 24-1.0228563-rdf-ntriples.txt
Original Record: 24-1.0228563-source.json
Full Text
24-1.0228563-fulltext.txt
Citation
24-1.0228563.ris

Full Text

White Dwarf Populations in Globular ClustersbyRyan GoldsburyB. Sc., Clemson University, 2009M. Sc., The University of British Columbia, 2011A THESIS SUBMITTED IN PARTIAL FULFILLMENTOF THE REQUIREMENTS FOR THE DEGREE OFDoctor of PhilosophyinTHE FACULTY OF GRADUATE AND POSTDOCTORAL STUDIES(Astronomy)The University of British Columbia(Vancouver)March 2016c© Ryan Goldsbury, 2016AbstractThis thesis details three distinct projects that explore stellar populations in Milky Way globularclusters. In the first, a method of modelling mass segregation in clusters is presented. The modelis fit to 54 clusters and the best fit parameters are presented in tabular form. The newly derivedparameter that indicates the amount of mass segregation correlates strongly with other dynamicalcluster parameters. In the second study, white dwarf data in the cluster 47 Tucanae are used toconstruct an empirical relation between temperature and time for these stars. The modified data arecompared to theoretical cooling models from four different research groups. We find disagreementbetween all of the models and the data. The models are also inconsistent with each other. In thethird investigation, new UV white dwarf data in 47 Tuc is used to constrain the hydrogen massfraction and neutrino production rates in cooling white dwarfs. A much different approach fromthe second project is used. The data are left untouched and the model is transformed to the spacein which the data exist. Using the unbinned maximum likelihood statistic, the model’s parameterspace is explored with MCMC sampling. A constraint on the rate of neutrino production in whitedwarfs comes from this analysis.iiPrefaceAll of the work presented in this dissertation was conducted between October 2011 and October2015 at the University of British Columbia in Vancouver.Chapter 1 is original work not published elsewhere.Chapter 2 is published in Goldsbury et al. [2013]. The data were acquired in 2007 from HSTprogram 10775 with Principle Investigator Ata Sarajedini. The first stages of image processingwere performed by Jay Anderson. I am responsible for the concept, data analysis, and manuscriptcomposition.Chapter 3 is published in Goldsbury et al. [2012]. The data were acquired in 2010 from HSTprogram 11677 with PI Harvey Richer. The first stages of image processing were done by JasonKalirai. The original concept for the project came from Jeremy Heyl. I am responsible for all of thesubsequent analysis including the composition of the manuscript.Chapter 4 has been submitted for publication. The data were acquired in 2013 from HST pro-gram 12971 with PI Harvey Richer. I designed the observing program as well as the artificial startests. The first stages of image processing were done by Jason Kalirai. Jeremy Heyl introduced meto the statistical methods used in the analysis. I performed all of the analysis and composed themanuscript.iiiTable of ContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiiTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ivList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vii1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.1 Globular Clusters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.2 47 Tuc . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.3 The Colour Magnitude Diagram (CMD) . . . . . . . . . . . . . . . . . . . . . . . 11.4 Three Projects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21.5 Chapter 2: Quantifying Mass Segregation and New Core Radii for 54 Milky WayGlobular Clusters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41.6 Mestel’s Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51.7 White Dwarf Cooling Models . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71.8 Chapter 3: An Empirical Measure of the Rate of White Dwarf Cooling in 47 Tucanae 81.9 Chapter 4: White Dwarf Cooling and 47 Tuc Revisited . . . . . . . . . . . . . . . 92 Quantifying Mass Segregation and New Core Radii for 54 Milky Way Globular Clus-ters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112.2 Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122.3 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122.3.1 Solving for Projected Density Distributions of King-Michie Models . . . . 132.3.2 Luminosity vs. Star Counts . . . . . . . . . . . . . . . . . . . . . . . . . 142.3.3 Selecting Groups Along the Main-Sequence . . . . . . . . . . . . . . . . . 162.3.4 Fitting Models to Subgroups . . . . . . . . . . . . . . . . . . . . . . . . . 182.3.5 Assigning Masses to Bins and Fitting Power Laws . . . . . . . . . . . . . 192.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 202.5 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22iv3 An Empirical Measure of the Rate of White Dwarf Cooling in 47 Tucanae . . . . . . 253.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 253.2 Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 263.3 Selecting our White Dwarf Sample . . . . . . . . . . . . . . . . . . . . . . . . . . 263.4 Fitting White Dwarf Temperatures . . . . . . . . . . . . . . . . . . . . . . . . . . 283.5 Determining the Timescale . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 313.5.1 A Sorted List of Temperatures . . . . . . . . . . . . . . . . . . . . . . . . 313.5.2 The Rate at Which Stars are Leaving the Main Sequence . . . . . . . . . . 323.6 The Cooling Curve . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 393.6.1 The Models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 393.6.2 Model Comparison . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 393.7 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 444 White Dwarf Cooling at High Temperatures: The Effect of Neutrino Production onCooling Rates . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 474.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 474.2 Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 484.3 Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 484.3.1 Artificial Star Tests . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 484.3.2 MESA Models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 534.3.3 Moving the Model to Data-Space . . . . . . . . . . . . . . . . . . . . . . 554.3.4 The Unbinned Likelihood . . . . . . . . . . . . . . . . . . . . . . . . . . 574.3.5 Degenerate Parameters and Prior Knowledge . . . . . . . . . . . . . . . . 594.3.6 Restricting the Fitting Region . . . . . . . . . . . . . . . . . . . . . . . . 594.3.7 Markov Chain Monte Carlo Sampling . . . . . . . . . . . . . . . . . . . . 604.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 614.4.1 Hydrogen Layer Thickness (qH) and Neutrino Scaling (An ) . . . . . . . . . 614.4.2 Neutrino Species . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 624.4.3 Fit Quality . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 644.4.4 Cooling Models From Other Groups . . . . . . . . . . . . . . . . . . . . . 664.5 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 675 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 705.1 Cluster Mass Segregation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 705.2 White Dwarf Cooling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72vList of Tables2.1 Fit Correlations (the blue dashed lines from Figure 2.5) . . . . . . . . . . . . . . . 222.2 Fit Powerlaw Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 244.1 Exposure times in seconds, for each of the nine visits . . . . . . . . . . . . . . . . 484.2 Values for N˙ are determined as described in the paragraph following this table. Theyare provided in units of objects per Myr. . . . . . . . . . . . . . . . . . . . . . . . 594.3 The “fit probability” is defined as the probability that data drawn from the givenmodel would produce a likelihood value less than the real data do when comparedback to that model. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67viList of Figures1.1 A colour magnitude diagram labeled with the most prominent features. . . . . . . . 32.1 Normalized surface density profiles. We have recalculated Figure 4-9b from BT87,showing the projected density distribution for various boundary conditions. Thesedistributions correspond to (from left to right)Y(0)=s2= 12;9;6;3 or c= 2:74;2:12;1:26;0:67,since there is a one-to-one relation between Y(0)=s2 and c. . . . . . . . . . . . . 152.2 The left hand panel shows core radius values taken from Harris [1996] plottedagainst values determined using our approach of fitting the cumulative star countdistribution. The center panel shows values determined by fitting the cumulativedistribution in which objects are weighted by their luminosities. The right panelshows fits to the mass weighted cumulative distributions. The red line highlightsa non-biased relation in each case. The horizontal error bars represent 1s and theHarris values have no error plotted. . . . . . . . . . . . . . . . . . . . . . . . . . . 162.3 The algorithmically determined fiducial sequence and the regions for the ten se-lected bins are shown for the cluster 47 Tucanae (NGC 104). . . . . . . . . . . . . 172.4 The best-fitting power-law model and likelihood of fit parameters for 47 Tuc. Simi-lar figures can be found in the online supplement for Goldsbury et al. [2013], show-ing power laws fit to both r0 and Rc for all 54 clusters in the catalogue describe in2.2. A and B are defined in equation 2.12. . . . . . . . . . . . . . . . . . . . . . . 202.5 Correlations between our fit parameters and the core relaxation time of the cluster.Here A has been scaled by the distance to each cluster so that it is in units of parsecs.Red points indicate core-collapsed clusters not used in calculating the correlationstatistic. The insets show the histogram of measured Pearson correlation coefficientsafter bootstrapping and taking into account uncertainties in the fit parameters. . . . 213.1 The footprint of the ACS (central star shaped region) and WFC3 (surroundingswath) observations. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 273.2 The CMD after making selection cuts in error, sharpness, and c2 of the PSF fit. Thefinal WD selection region is shown in the enclosed area (red). 887 objects are foundwithin this enclosed region. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28vii3.3 A model WD spectrum with an effective surface temperature of 10,000K and asurface gravity of log(g)=7.5. Throughput curves for our two filters are shown over-plotted on a separate axis scale on the right. The spectrum shown is the resultof 3D hydrodynamic simulations of white dwarf atmospheres. Qualitatively, it is ablackbody spectrum with hydrogen absorption features corresponding to the Balmerseries transitions (n=2). The continuum suppression to the short wavelength side of364nm is referred to as the Balmer Jump in stellar astrophysics. The underlyingblackbody spectrum would follow the continuum seen to the right of this jump. Thesuppression to the left is due to continuum absorption. Photons of higher energythan 364nm can continue to ionize hydrogen atoms in the n=2 state, and so thecontinuum is suppressed beyond this limit. . . . . . . . . . . . . . . . . . . . . . . 303.4 The standard deviation of input – output magnitude derived from artificial star testsdescribed in Kalirai et al. [2012] (for all fields combined). . . . . . . . . . . . . . . 303.5 Completeness fraction of our WD sample as a function of F606W magnitude (aver-aged over all fields). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 323.6 The fitted Teff as a function of object number for our WDs. The object number isproportional to the time on the white dwarf sequence, and is calculated as shownin equation 3.4. The y-axis is just the temperature fit to each object’s SED, as dis-cussed in section 3.4. The results from both data sets agree well within uncertainties.Residuals are shown as Data−ModelData (the fractional residual). . . . . . . . . . . . . . 333.7 A comparison of the cumulative radial distributions of the entire WD sample and thelower RGB stars. A KS test indicates that there is a 73% chance that two samplesof this size drawn from the same underlying distribution would differ by at least thismuch. We therefore find no evidence of significantly different radial distributionsfor the WD and RGB stars. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 343.8 The CMD from Stetson [2000], including only objects that fall within the region ofour WFC3 data. The region of the CMD between points P1 and P2 is fit with modelsfrom Dotter et al. [2008a] to determine the time it takes an object to evolve betweenthe points. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 353.9 Likelihood of ∆A fit to our selected magnitude range, after marginalizing out allother model parameters. ∆A corresponds to the time it takes an object to evolvefrom point P1 to point P2 in Figure 3.8. . . . . . . . . . . . . . . . . . . . . . . . . 363.10 A cumulative distribution of Stetson’s RGBmagnitudes between P1 and P2 in Figure3.8 compared to the expected distribution assuming a constant evolution rate. AKS test suggests that there is a 40% chance that a sample drawn from the modeldistribution would differ by at least as much as this data. . . . . . . . . . . . . . . 373.11 The assumed constant formation rate is shown in red. The formation rate shown inblack is the same as derived in Equation 3.11. . . . . . . . . . . . . . . . . . . . . 383.12 The cumulative number of white dwarfs formed in our field over the given time isshown. The red line results from assuming a constant rate of formation. The blackline results from the assumptions of a Salpeter IMF outlined above. . . . . . . . . . 38viii3.13 The measured 1 and 2 s uncertainties of our cooling curve (the shaded blue regions).The shape in this region is well constrained as shown in Figure 3.6, but the agescaling leads to a large systematic uncertainty in the x-direction. This is the regionrepresented by these contours. Models of white dwarf cooling rates from Fontaineet al. [2001a], Lawlor and MacDonald [2006b], Hansen et al. [2007a], and Renedoet al. [2010] are shown for comparison. . . . . . . . . . . . . . . . . . . . . . . . 403.14 This plot shows measured cumulative distribution of our objects found to have tem-peratures between 7,000 K and 30,000 K (813 objects WFC3, 188 objects ACS), aswell as those predicted by the four models shown in Figure 3.13. . . . . . . . . . . 413.15 Cumulative distribution of our fit temperatures compared to the models when as-suming values at the extreme ends of literature distance measurements for 47 Tuc. . 413.16 Results from our iterative error analysis using data from WFC3. Each distributionis a histogram of D values found when repeatedly comparing samples from eachmodel back to the theoretical model itself. The blue line indicates the value of Dcalculated when comparing our real data to each model. Integrating from the blueline to infinity gives the probability that data drawn from each model would differ byat least as much as our real data does. All of these cooling models are inconsistentwith our data at the 3s level. This is true even after considering a conservativelylarge systematic uncertainty in distance and a random uncertainty in temperature. . 423.17 Results from our iterative error analysis using data from the ACS field. These showthe same thing as figure 3.16, but using the ACS data rather than the WFC3 data. . 433.18 Cumulative distributions of our object magnitudes in each filter compared with thosepredicted by the models. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 453.19 The luminosity function of our WFC3 data compared to the models over the F390Wrange [ 21 , 26.5 ]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 464.1 The pattern of the observed fields in the core of 47 Tuc. Fields imaged by WFC3are shown in blue, while those observed by ACS are in pink. In the final data setonly nine of the ten orientations for both WFC3 and ACS were used, due to loss ofguide stars during one exposure. . . . . . . . . . . . . . . . . . . . . . . . . . . . 494.2 The colour magnitude diagrams of the two fields. The central UV data are shownon the left. The visible CMD from the surrounding annulus is shown on the right.In each plot there are three distinct sequences, although in the left plot the centralsequence is quite faint and sparse. From left to right these are the white dwarfsequence of 47 Tuc, the main sequence of the Small Magellanic Cloud, and themain sequence of 47 Tuc. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 504.3 The real data are shown in black and the artificial star input grid is shown in red.The magnitude plot is centred on the white dwarf sequence. . . . . . . . . . . . . . 514.4 Completeness fraction in magnitude space for a fixed radial distance is shown on thetop. Completeness fraction along the two-dimensional plane of R and F225Winputwith a fixed value of F336Winput is shown on the bottom. Radial distance is in unitsof WFC3 pixels, which correspond to 0.04 arc seconds each. These surfaces are thevalues from C described in Equation 4.1. C is indexed by three parameters. In eachplot above, one of these three is fixed, andC is shown varying across the remainingtwo. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52ix4.5 1, 2, and 3s contours (enclosing 68%, 95%, and 99:7%) of the photometric errordistribution are shown at two points in the data-space. The top panel shows the errordistribution for faint stars close to the outside of the inner WFC3 field. The lowerpanel shows the error distribution for brighter stars closer to the centre. R is againin units of pixels. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 544.6 The colour magnitude diagrams of the two fields. Overplotted are the region inwhich we fit our model in blue, and the white dwarf cooling sequence in red for themean prior values given in Table 4.2, before convolution with the error model. . . . 604.7 The auto-correlation of each chain parameter as a function of the lag between pointsin the chain. No memory of the past position in parameter space is retained beyond40 steps. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 614.8 A comparison of the likelihood distribution for each of our model parameters shownfor the entire chain in red. The thinner black lines are independent segments 1000steps long. Each of these separate chains is much longer than the length scale onwhich correlations persist in the chains. All of them average to the same broaddistribution that we find with the entire chain, but with more statistical noise. . . . 624.9 Contours of the joint fit distribution of An and qH. 1, 2, and 3s contours (enclosing68%, 95%, and 99:7%) are shown. . . . . . . . . . . . . . . . . . . . . . . . . . . 634.10 A histogram of the posterior distribution of qH values from the Markov chain. Thisis the likelihood distribution of qH after marginalizing over all other parameters inthe model with prior distributions on each, except for An . The range of 2:0×10−5to 8:8×10−5 contains 95% of the distribution (shown in red). . . . . . . . . . . . 644.11 A histogram of the posterior distribution of An factor values from the Markov chain.This is the likelihood distribution of An after marginalizing over all other parametersin the model with prior distributions on each, except qH. The range of 0:80 to 1:21contains 95% of the An fit distribution (shown in red). . . . . . . . . . . . . . . . . 654.12 The relation between luminosity and time shown for varying An values while qH isheld fixed. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 664.13 The relation between luminosity and time shown for varying qH values while An isheld fixed. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 674.14 The likelihood of sin2qW and n. For the accepted value of sin2qW (0.23155 fromOlive et al. [2014]), C′V is very small. This means that the dependence of An on nis very weak. As a result, our constraint on the neutrino production (An ) does notreasonably constrain n. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 684.15 The distribution of likelihood values expected when sample data are drawn from thebest fitting cooling model is shown as the black distribution. The red line indicatesthe log-likelihood value found when comparing the real data to the model. Wewould expect to draw data from our best fitting model that differs at least as muchas the real data do from that model 4.9% of the time. . . . . . . . . . . . . . . . . 69xChapter 1Introduction1.1 Globular ClustersGlobular clusters are gravitationally bound, spherical populations of stars found around virtually alllarge galaxies [Harris, 1991]. Globular clusters are distinct from galaxies which are much larger,and also from similar sized objects such as dwarf spheroidal galaxies in that globulars do not showany evidence of significant dark matter components. In our own Galaxy there are approximately150 known globular clusters with stellar populations from 104 to 106 stars. Milky Way globularclusters are found in a spherically symmetric distribution centred on the Galactic core. Their threedimensional density in the Galaxy is roughly r(R) µ R−3 where R is the distance from the centreof the galaxy [Harris, 1976]. A correlation between galactocentric distance and globular clustermetallicity is observed, with more metal rich clusters found closer to the centre of the Galaxy[Harris, 1991]. The globular clusters in the Milky Way appear to have formed over the same epochas the Galaxy itself [Krauss and Chaboyer, 2003], and their stellar populations have not gained anysignificant amount of stars from the rest of Milky Way population since their formation.1.2 47 TucChapters 3 and 4 both deal exclusively with the globular cluster 47 Tucanae. 47 Tuc is located inthe Southern hemisphere in the constellation Tucanae. This cluster is of interest for a number ofreasons. 47 Tuc is the seventh closest globular cluster to us at around 5 kpc, but it is the secondbrightest since it has many more stars than any of the other nearby clusters. The interstellar spacebetween us and 47 Tuc happens to be relatively free of dust and gas, and so the scattering andabsorption of light between us and 47 Tuc is minimal. These factors make 47 Tuc a better target forstudying stellar populations than any other globular cluster [Harris, 1996].1.3 The Colour Magnitude Diagram (CMD)Throughout this thesis, the Colour Magnitude Diagram (CMD) is referenced as a way to graphicallydescribe the population of stars in a cluster. To generate a CMD of a population of stars, measure-ments of magnitude or intensity in two filters are required. The difference of these two magnitudes,the bluer filter minus the redder filter, is then used as the “colour” while either filter individually can1be used as the “magnitude”. The colour is plotted on the x-axis and the magnitude is plotted on they-axis. For a population of stars that is at a fixed distance from the observer, the x-axis is related tothe surface temperature of those stars, while the y-axis is related to their intrinsic luminosity. Themajor features of the 47 Tuc CMD used in Chapter 4 can be seen in Figure 1.1.On the right side of the CMD there is the main sequence. At the top of this sequence starsbend off to redder colours (toward the right of the diagram) as they enter the giant branch. At thispoint, they expand significantly and their surface temperature drops drastically. All of the stars inthis cluster have formed at approximately the same time, and more massive stars will have a shorterlifespan on the main sequence. As time goes on, this bend off of the main sequence will move downto lower and lower mass stars.After shedding much of their mass (around 40% for the stars currently leaving the main sequencein this figure) the stars end up as very hot and very dense white dwarfs. On the left side of the CMDthe sequence of these objects can be seen. Newly formed white dwarfs first arrive at the top of thissequence, and then travel down as they cool over billions of years. The time spanned by the red linealong the white dwarf sequence is approximately two billion years. White dwarfs at the bottom ofthis sequence today would be found at the top if we could observe the cluster from two billion yearsago.1.4 Three ProjectsEach chapter of this dissertation, other than this one, is a discussion of a published article that wascompleted during my graduate work. Because the work summarized in this dissertation has alreadybeen published elsewhere as self contained documents, each chapter can be read on its own andunderstood without reference to or understanding of the dissertation as a whole. The three mainchapters describe three separate but interconnected projects.All of the following projects involve the study of various stellar populations in Milky Wayglobular clusters. They touch on stellar dynamics, stellar evolution, and statistical analysis. Theproject in Chapter 2 discusses mass segregation in globular clusters and how to model it in a simpleway. The remaining two chapters focus on how white dwarfs cool, from a data driven perspective.I have used the data in a way that puts the tightest possible constraints on what current theoreticalmodels can say about cooling white dwarfs. This has largely been made possible by two recentadvances. The first is the acquisition of Hubble Space Telescope (HST) data in the core of nearbyglobular clusters in which white dwarfs are easily resolvable and distinguishable from the rest ofthe stars in the cluster. The second is the advancement of computing power and software. All of thestatistical analysis in this thesis was performed on a single desktop computer. Twenty years ago thiswould not have been possible with the methods we have used.Our recent HST data sets have globular cluster white dwarf samples that number in the thou-sands, rather than the tens or hundreds, as had previously been the case. By looking at a globularcluster, rather than the sky as a whole, we are able to eliminate a number of corrections that wouldneed to be made to a more diverse sample, such as a field in the disk or halo of the Galaxy.Rather than having to determine distances to individual white dwarfs, we know that all of thewhite dwarfs in our sample are at essentially the same distance from us. It is still necessary to deter-mine the distance to the cluster, but this shift affects the data set as a whole, rather than introducingvariance within the photometry of individual white dwarfs. The same is true of the reddening ofthe white dwarfs’ spectra, since all of the light from the cluster is passing through the same inter-2Figure 1.1: A colour magnitude diagram labeled with the most prominent features.3stellar material on its way to Earth. Globular clusters also form over a very short time comparedto their current age, and form from a uniform chemical composition of gas. The time for a star toburn through all of the hydrogen in its core depends on the mass of the star. Observing the clusterbillions of years after its formation, as we do, means that all of the very massive stars have alreadyleft the main sequence and all of the very low mass stars have not yet. All of the stars that we cur-rently see leaving the main sequence and evolving up the giant branch are roughly the same mass.Consequently these will form white dwarfs of roughtly the same mass. This means that all of theyoung white dwarfs observable in a cluster at a given time will have more or less the same massand age (age since the formation of their progenitors, not since their transition to white dwarfs).Carrying out a study like this on the white dwarfs in the Galaxy as a whole necessarily brings inthe history of star formation in the Galaxy, and the range of chemical compositions in the materialfrom which these stars formed. When looking at white dwarfs in an individual globular cluster, itis approximately correct to think of the formation history as being one instantaneous burst from aprogenitor gas cloud of uniform chemical composition.Strictly speaking, none of these assumptions are true. The stars are all roughly at the samedistance, but the cluster does have some physical size, and so the distance to individual white dwarfswill vary across the extent of the cluster. Almost all globular clusters exhibit multiple distinct stellarpopulations, formed in separate bursts during the first billion or so years of the clusters’ lifetimes[Gratton et al., 2012]. Clusters also have a population of massive blue straggler stars, which arethought to be formed through mergers of less massive stars. These stars are continuously forming,evolving, and creating their own populations of more massive white dwarfs [Heyl et al., 2015]. Someclusters also have significant differential reddening, in which one part of the cluster is obscured bymore interstellar material than another part. So every single assumption listed in the precedingparagraph has its caveats.However, the globular cluster in which we study white dwarf cooling is 47 Tuc. This cluster isreasonably close by, experiences very little interstellar reddening, and has no measurable differentialreddening. There is some evidence that 47 Tuc exhibits multiple stellar populations [Milone et al.,2012], but these populations differ minimally. The white dwarfs that recently evolved from thesepopulations will not be measurably separate given even the best currently available data (which wepossess).The uniformity of distance, extinction, chemical composition, and age of white dwarfs in 47Tuc allows for a much closer look at the physics of their cooling than would be possible with a morediverse sample population.1.5 Chapter 2: Quantifying Mass Segregation and New Core Radiifor 54 Milky Way Globular ClustersThe first globular cluster structure models were introduced contemporaneously by Ivan King [King,1966] and Richard Michie [Michie, 1963], and these models are commonly referred to as King-Michie models. Such models describe the distribution of massive particles in a smooth gravitationalpotential. They rely on the assumptions that all of the objects in the system are of the same mass,and that the system is in equilibrium.Another commonly used assumption is that globular clusters have a constant mass-to-light ratio.Combining these models with this assumption leads to the conclusion that the surface brightnessdistribution of a cluster is the same as the mass distribution of the cluster, given the appropriate unit4conversion. All of these assumptions turn out to be false when probed in detail.Chapter 2 focuses on how these assumptions break down, and how to account for this as simplyas possible. Firstly, globular clusters have stars with a range of masses from sub-stellar browndwarfs (0:08 M⊙) up to massive blue stragglers (1:6 M⊙), so the assumption that all of the objectsin a globular cluster have the same mass is incorrect. This has been known even since King andMichie first introduced their models. These assumptions were necessary to make the calculationstractable without modern computing power.Removing the assumption that all objects in the cluster are the same mass does not have to breakthe proportionality between surface brightness and mass density. If all of the objects move aboutthe cluster in the same distribution regardless of their mass, then the preferred model would still bea linear combination of solutions to the King-Michie models. But, it turns out that this assumptionis also untrue. The distribution of stars within the cluster depends strongly on their mass.Over time periods that are long compared to the time it takes a single object to cross the cluster,individual stars will have many interactions with other objects in which kinetic energy exchangeoccurs. The distribution of stars will tend toward equipartition of kinetic energy. Stars that are moremassive will end up with smaller velocities and will “sink” to tighter orbits near the centre of thecluster. Less massive stars will have larger velocities and will move to larger orbits and a morediffuse distribution.Even with this further complication, we could still salvage the analogy between the surfacebrightness distribution and the mass distribution of the cluster if the brightness of individual starsis proportional to their mass. But, like all of the assumptions so far, this one breaks down too.Massive stars put out much more energy per unit mass than light ones (L µM3:5 for main sequencestars). So this means that the measured surface brightness profile of a cluster is not the same as themass distribution of the cluster, nor is it the same as the distribution of stars in the cluster. Becausebrighter stars are more centrally concentrated, both the number of stars and their total mass will bemore diffuse than the distribution of surface brightness in the cluster.Despite the fact that all of these basic assumptions are understood to be false in detail, for thesake of simplicity most published papers still use a single core radius value to describe the overalldistribution of a globular cluster. The work in this chapter expands on this assumption by presentinga way to describe how the core radius changes with the mass of objects. The analysis indicates asystematic bias in core radius values currently found in the literature, which follows directly fromthe breakdown of the assumptions required by the King-Michie models. To account for this, weprovide a simple expanded model. This model introduces one additional parameter (the slope of apower law) that accounts for how the central concentration of stellar populations changes with mass.With this parameter, the distribution of objects of any mass can be approximated within a globularcluster.1.6 Mestel’s TheoryMost historical summaries of white dwarf cooling naturally start with Mestel [1952]. Mestel wasthe first to publish a predictive model of white dwarf cooling. This model is quite simple, butremarkably accurate in describing how white dwarfs decrease in luminosity in the middle (6,000K< T < 20,000K) of the sequence. The qualitative picture of the white dwarf in Mestel’s modelis as follows. There is an isothermal degenerate core that makes up 99:99% of the white dwarf’smass. This core is surrounded by a thin (relative to the core) non-degenerate hydrogen atmosphere5in thermal equilibrium. The opacity of the atmosphere is determined only by bound-free and free-free absorption. Bound-free absorption occurs in this case when an electron bound to a hydrogenatom is ejected due to the absorption of a photon. Free-free absorption occurs as free electronsinteract with photons traveling through the atmosphere and they exchange energy. The white dwarfhas no nuclear energy sources. It is no longer contracting and so there is no energy radiated dueto gravitational collapse. The sole energy source in the white dwarf is the loss of heat from ionsin the core. Mestel’s model predicts that the Luminosity of a white dwarf should be proportionalto time−75 , and this is roughly what is observed for white dwarfs in this middle temperature range.With a power law such as this, the white dwarfs cool more and more slowly the cooler they get. Ifone were to look at a very large sample of white dwarfs, one would expect these to be filled withmany more cool white dwarfs than hot ones. This also means that Mestel’s model does a fairlygood job of predicting the overall distribution of observed white dwarf temperatures for very oldpopulations.A summary of Mestel’s model following Van Horn [1971] is presented here. Mestel’s modelassumes that the only source of energy in the white dwarf is the cooling of the ions. This means thatthe heat capacity per ion isCv =32kAH; (1.1)where k is Boltzmann’s constant, A is the atomic mass, and H = 1:66044×10−24g per atomicmass unit. The luminosity comes entirely from the core temperature (Tc) changing with time:L=−32kMAHdTcdt: (1.2)The second major assumption Mestel makes is that the opacity of the atmosphere is describedby Kramer’s opacity law:k = k0rT−3:5: (1.3)k0 is some collection of constants that we aren’t concerned with. The general form and depen-dence on T−3:5 are what is important. Next we consider the conditions of radiative equilibrium andhydrostatic equilibrium, respectivelyLr =−4pr2 4sc3T 3krdTdr; (1.4)where s is the Stefan-Boltzmann constant and c is the speed of light, anddPdr=−rGMrr2: (1.5)Rearranging these and inserting in Equation 1.3 for k to find dPdT givesdPdT=4sc34pGMk0LrT 6:5: (1.6)6At this point we substitute in the equation of state for an ideal gas,P=1mkHrT; (1.7)to eliminate r and then solve Equation 1.6 for P, which givesP2 =18:54sc34pGMk0LkmHT 8:5: (1.8)We also consider the relation between density and temperature at the boundary of the degenerateand non-degenerate parts of the star (from [Schwarzschild, 1958]):1mr = 2:40×10−8T 1:5: (1.9)After combining Equations 1.7 and 1.9 and substituting into Equation 1.8 to eliminate r and Pwe find that:L µ T 3:5c : (1.10)Taking the result of Equation 1.10 with 1.2 leads directly to:L µ t−75 : (1.11)1.7 White Dwarf Cooling ModelsBut, the assumptions inherent in Mestel’s model do not hold in the case of very young and hot whitedwarfs or very old and cool white dwarfs. Above temperatures of 30,000K the dominant source ofenergy loss from white dwarfs is the production of neutrinos in the core from the interaction of pho-tons and plasma [Itoh et al., 1996]. White dwarfs also undergo nuclear burning in their atmosphereand gravitational contraction at these hot temperatures early in their cooling. As the white dwarfsbecome very cool, Mestel’s assumptions break down in a different way. First, convection begins tooccur at the interface between the degenerate interior and the hydrogen/helium layers. Eventuallythis becomes the dominant mode of energy transport into the upper atmosphere. Then, at cool tem-peratures, the white dwarf interior undergoes phase changes from a degenerate gas, to a liquid, andeventually to a solid crystal. Mestel’s simple model does not include any of the above effects, andso there has been plenty of room for improvement.In 1971, Van Horn systematically laid out how to improve on each of the assumptions in theMestel model [Van Horn, 1971]. Van Horn argued that nuclear burning and gravitational contractionare unlikely to contribute significantly to the cooling of white dwarfs. However, he found that:plasma neutrino losses are important at high temperatures; the electronic heat capacity can be upto 50% of the heat capacity of the ions for carbon white dwarfs; Coulomb interactions lead tophase changes in the interior of cool white dwarfs; and convection zones that extend from theenvelope into the degenerate interior are also important. A few years later Lamb and Van Horn[1975] presented the results of simulations of cooling and crystallizing pure C white dwarfs. Thephase changes in white dwarfs release latent heat and have a very dramatic effect on the predictedluminosity function at the cool end. Winget et al. [1987] used the models of Lamb and Van Horn7[1975] along with a number of additional assumptions about the formation of our Galaxy to predictthe luminosity function of white dwarfs in the Galactic disk. In what was probably the first exampleof white dwarf cosmochronology in practice, these models were used to constrain the age of theGalaxy. Hansen and Phinney [1998] introduced hydrodynamic atmosphere simulations that have aconsiderable effect on the convection that occurs in cooler white dwarfs.All of these developments were incorporated by Fontaine et al. [2001a], who produced manypublicly available cooling models. These models have varying core compositions from pure carbonto pure oxygen. They looked at the effects of the distribution of carbon and oxygen within the core,and also the distribution and amounts of helium and hydrogen within the atmosphere. The resultsare available, and useful to anyone wishing to model white dwarf populations.The work in Chapter 3 and Chapter 4 looks at these models and compares them to very recentdata sets. Models from other groups generated in the last decade are also included [Lawlor andMacDonald, 2006a, Hansen et al., 2007b, Renedo et al., 2010]. Chapter 4 focuses mostly on modelsmade available through an even more modern stellar evolution code known as MESA [Paxton et al.,2011], which stands for Modules for Experiments in Stellar Astrophysics. The benefit of MESA isthat one can run stellar evolution simulations all the way from the pre-main sequence to very cooland low mass white dwarfs. This allows one to connect the physics of stellar evolution and massloss along the giant branch with the physics of white dwarf cooling, and to constrain one consistentmodel by fitting simultaneously to both types of stars in a cluster.1.8 Chapter 3: An Empirical Measure of the Rate of White DwarfCooling in 47 TucanaeThe goal of the project Chapter 3 summarizes is to generate an empirical relation between thetemperature of a white dwarf and the time since the mass loss that brought it to the white dwarfsequence.The first step in this process is to identify which objects in our data are white dwarfs, and to fiteach of these objects with a white dwarf spectral model to determine their temperatures. After eachobject is assigned a temperature, the list of temperatures is sorted from hottest to coolest. Becausewe expect white dwarfs to decrease in temperature monotonically with time, the ordering in this listis a proxy for the time since the formation of the white dwarf. The positions in the sorted list areconverted to times by assuming a constant formation rate of the white dwarfs. This assumption issupported by an analysis of the cluster’s giant branch. The cumulative distribution of stars along thegiant branch is compared to the cumulative distribution predicted using a main sequence isochronefrom Dotter et al. [2008b] and an assumed constant formation rate. These distributions are foundnot to be measurably different.After this relation between temperature and time is generated from the data, it is comparedto predictions from various white dwarf cooling models. Models from four independent groupsare compared to the relation generated from the data. None of the models are able to accuratelyreproduce the structure seen in the data. The data used in the comparison includes many more hightemperature (and therefore young) white dwarfs than any previous analysis. This fact means thatthe quality of the model fit hinges on a region of the cooling curve that has not been well studied,and it is therefore not surprising that the models do not do a good job of fitting these data.81.9 Chapter 4: White Dwarf Cooling and 47 Tuc RevisitedThe project Chapter 4 summarizes expands the white dwarf cooling model testing done in Chapter3 to a larger data set and to even hotter stars. When white dwarfs leave the giant branch and shedtheir outer layers, they initially have surface temperatures of over 100,000K. Within only ten millionyears, they have cooled to 30,000K. It takes another hundred million years to reach 20,000K andalmost a billion years to reach 10,000K [Fontaine et al., 2001a]. This means that very hot whitedwarfs (> 30,000K) are extremely short-lived objects. White dwarf luminosity functions studied inthe past have been completely dominated by old cool white dwarfs.To reliably find young hot white dwarfs, we look in the core of one of the nearest clusters,47 Tuc. For the reasons described in Chapter 2, the core is where the highest density of verymassive stars are found. Young white dwarfs were just recently giant stars. Shedding some oftheir mass and transitioning from giants to white dwarfs, has happened relatively recently comparedto the dynamical timescale of the cluster. These objects have lost around 40% of their mass, andequipartition should lead to them receiving higher velocities and larger cluster orbits as a result.But, not enough time has yet passed for that to happen, and so their radial distribution should stillbe the same as those of the giants. This means that one should look for young white dwarfs in thesame place that one would look for giants; namely the core of the cluster. The problem with such asearch is that white dwarfs are typically faint compared to the giants that crowd the cluster centre.By looking in the near-UV (225 and 336 nm) we are able to easily locate thousands of white dwarfseven in the dense centre of the cluster. In fact, the brightest objects in our entire field (at 225 nm)are white dwarfs with temperatures of around 100,000K. This is despite the fact that they are morethan 10,000 times smaller than the giants.The motivation to look at hot young white dwarfs is two-fold. At very high core temperatures,the energy loss from white dwarfs is actually dominated by neutrinos. Our data only containsinformation about the photons leaving the star, but the photometric data reveal the temperature ofthe white dwarf, and how this temperature changes with time is heavily dependent on the amountof neutrinos produced. So even though we cannot observe these neutrinos directly, we can stillconstrain their effects by measuring how the spectrum of the white dwarf changes with temperature,and how that changes with time. The second benefit of looking at very young white dwarfs is thesensitivity to hydrogen layer thickness. As white dwarfs leave the giant branch they immediatelybegin to contract. How large they are initially and how quickly they contract are also dependent onthe amount of hydrogen left in their atmospheres.The data set analysed in this chapter is unprecedented in its ability to constrain these two key pa-rameters of white dwarf cooling models. Unlike previous chapters of this thesis, these observationswere designed specifically for this purpose.The approach used to constrain models in this analysis is quite different from what is describedin Chapter 2. The philosophy behind this analysis is to leave the data in as raw a state as possible.Rather than using artificial star tests to correct for systematics and completeness of stars, the resultsof these tests are built in to the model. The data are compared to the model using the unbinnedmaximum likelihood which requires no binning, as the name suggests. This approach works verywell for our data. We measure a magnitude in two filters and a position in the field for each starin our data set, which contains a few thousand stars. If we were to break this data space into onlyten bins along each axis (filter 1, filter 2, and radial position), that would be 1000 total bins. Onaverage there would only be 3–4 stars per bin. Reducing the data like this produces something that9does not meet the requirements of a statistic such as c2. The variance in the counts of each binwould be described by the Poisson distribution, and not the Gaussian distribution. In fact, this isalways true for counts in bins, but the for large N the difference between these two distributionsis not significant. Binning multi-dimensional data necessarily leads to small counts in individualbins, which then fail to satisfy the assumptions of commonly used statistics. It is possible to builda likelihood statistic assuming that the counts in each bin are described by the Poisson distributionrather than the Gaussian distribution. However, the binning in this case does not provide any benefit.If one is going to use Poisson statistics anyway, one might as well go all the way to an unbinnedapproach, which does not require bins to be defined and is computationally simpler. The unbinnedapproach is derived in in Chapter 4.The results of this analysis provide new constraints on both the thickness of the hydrogen layerin white dwarfs, and the amount of neutrino production in their cores.10Chapter 2Quantifying Mass Segregation and NewCore Radii for 54 Milky Way GlobularClusters2.1 IntroductionComprehensive catalogues containing Milky Way globular cluster parameters have been compiledby many different groups over the last two decades. The first large scale effort was put together byTrager et al. [1995], in which they presented surface brightness profiles for 125 Galactic globularclusters, and included parameters determined from model fitting for 63 of them. This work wasexpanded upon by McLaughlin and van der Marel [2005] who used data from the previous paper aswell as new data, and fit multiple classes of models to the surface brightness profiles. The resultsfrom both of these papers, as well as many others, have been compiled in Harris [1996] (also 2010edition).As shown in this chapter, the limitation of the approaches used by these previous studies isthat they do not distinguish between the surface brightness distribution and the underlying massdistribution of the cluster. In both Trager et al. [1995] and McLaughlin and van der Marel [2005] thecluster is assumed to have one mass-to-light ratio that is constant everywhere. This leads directly tothe conclusion that the surface brightness distribution is equivalent to the mass density distribution.In fact, all of the clusters in this study exhibit some degree of mass segregation. Mass segregationis a result of the dynamical interactions of stars within clusters. Over long timescales, in whichstars undergo many individual interactions with each other, they tend toward equipartition of theirkinetic energies. When two stars interact, there is a statistical tendency for the higher energy starto decrease its velocity and for the lower energy star to increase its velocity. Over a large numberof interactions this process leads to slower tighter orbits for more massive stars and faster moreextended orbits for less massive stars. Additionally, the more massive stars in the cluster generatemore luminosity per unit mass than the less massive stars.Taken together, these two effects generate the following picture. Looking at samples of starsthroughout the cluster, one would find a higher average mass the closer one looked toward thecentre of the cluster. Because more massive stars have a higher luminosity per unit mass, the mass-to-light ratio will decrease toward the centre of the cluster. Because the mass-to-light ratio is not11constant throughout the cluster, and in fact changes in a systematic and quantifiable way, the surfacebrightness distribution will not be an appropriate analog for the underlying mass distribution of thecluster if these effects are ignored.The approach presented in this chapter differs from these previous studies in that we do notassume a single mass-to-light ratio for the cluster. We instead work directly from star counts andconsider bins 0.8 magnitudes wide to break each cluster down further, fitting the distribution of starsin each bin independently. Using stellar evolution models from Dotter et al. [2008a] we can assignmasses to each bin and analyze how the stellar concentration changes with mass. This is equivalentto analyzing how the mass-to-light ratio changes with cluster radius. This approach avoids assuminga constant mass-to-light ratio for a cluster. Given the presence of any amount of mass segregation theassumption of constant mass-to-light ratio will not be true since mass-to-light ratio is a function ofmass, and the average mass of stars changes as a function of distance from the cluster core. Becauseof this, this surface brightness profile does not necessarily reflect the underlying mass density in acluster.2.2 DataAll of the data used in this study are from the ACS Survey of Galactic Globular Clusters [Sarajediniet al., 2007]. Each cluster was observed with the Advanced Camera for Surveys (ACS) / WideField Camera (WFC) on Hubble Space Telescope (HST) for one orbit through each of F606W andF814W . These are filters centred at 606nm and 814nm. Each orbit contained one short exposure (tofill in the saturated stars) and four to five deep exposures, which were stepped over the ACS chipgap for more uniform spatial coverage. The catalogue was constructed by analyzing the field patchby patch in all the individual exposures simultaneously. In each patch, the brightest stars were foundfirst and fitted with a Point Spread Function (PSF) that was tailored to each particular exposure, thensubtracted to allow fainter neighbors to be found. Only stars that clearly stood out above the knownsubtraction errors and known PSF artifacts in multiple exposures in both filters were included inthe catalog. The very few stars that were saturated even in the short exposures were measured byfitting the PSF to the surrounding unsaturated pixels. This entire data catalogue is now public onDr. Sarajedinis personal webpage [Sarajedini, 2011]. A very thorough discussion of the reductioncan be found in Anderson et al. [2008].2.3 MethodsThere is some ambiguity in the literature regarding a few of the parameters commonly used tocharacterize the concentration of star clusters. In an attempt to avoid this, we will use the followingsymbols and definitions for the remainder of the chapter:• core radius; Rc ;S(Rc) =12S(0); (2.1)12• King radius; r0 ;r0 =√9s24pGr0; (2.2)• tidal radius; rt ;• concentration; c ;c= log10(rtr0): (2.3)The quantity S here refers to the projected density distribution of the cluster. The core radius(Rc) is the projected radial distance from the centre of the cluster at which the projected density ofthe cluster drops to half of the central value. The quantity s is the velocity-dispersion parameter,and r0 is the central density of the King model. These two parameters will be explained in moredetail in the next section. The King radius (r0) can be calculated from the previous two quanitiesand describes the scale of the model similar to the core radius. The concentration (c) is definedfollowing the convention in Binney and Tremaine [1987] (hereafter BT87). This disagrees with thedefinition used by Harris [1996], in which c= log10(rt=Rc). For a given r0 and rt a King model hasonly one possible value of Rc, but it is not equivalent to r0. It is also important to note here that weuse the capital R to denote a projected radius, while the lower-case r refers to a three-dimensionalradius. Equation 2.1 is intentionally ambiguous as it does not indicate whether the surface densityis in units of luminosity, mass, or number of stars. This will be expanded upon in Section 2.3.2.2.3.1 Solving for Projected Density Distributions of King-Michie ModelsOur method involves generating a number of King-Michie [King, 1966, Michie, 1963] models,which can be derived as follows. First, consider The Poisson equation for a spherically symmetricgroup of stars with distribution function over energy and momentum f (E;L):1r2ddr(r2dFdr)= 4pGr = 4pG∫f(12v2+F; |r¯× v¯|)d3v¯: (2.4)This can also be written in terms of a relative energy e and relative potential Y defined asY≡−F+F0 and e ≡−E+F0 =Y− 12v2: (2.5)The simplest type of spherically symmetric models are those that have an isotropic velocitydistribution. This means that f = f (e). The distribution function King used is one such distribution:fK(e) =r1(2ps2)− 32 (ee=s2−1) e > 00 e ≤ 0:(2.6)The second case ensures that stars with velocities greater than the escape velocity at a givenradius (vescape =√2Y) will not be a part of the distribution. If we substitute the relative energy13from Equation 2.5 and integrate over velocity from 0 to√2Y, the density of the King model isfound to berK(Y) = r1[eY=s2erf(√Ys)−√4Yps2(1+2Y3s2)]: (2.7)Poisson’s equation for Y can be writtenddr(r2dYdr)=−4pGr1r2[eY=s2erf(√Ys)−√4Yps2(1+2Y3s2)]: (2.8)Y(0)=s2 is referred to as the central dimensionless potential and is sometimes written Y0. Wecan also see from Equation 2.7 that there will be some radius at whichY= 0 and so rK = 0. This isthe tidal radius rt. There is a one to one relation betweenY(0)=s2 and c= log10(rtr0). These relativeparameters control the shape of the resulting density distribution. The scale of the distribution isdetermined by any one of Y(0), s2, r0, or rt. In practice, we do our fitting in the space of c and r0.We follow the prescription in Section 4.4 of BT87, ending with a series of models giving nor-malized surface density (S=(r0r0)) as a function of radius (r=r0). To calculate density profiles fromEquation 2.8, we use two boundary conditions. The first is always the same: dY=dr = 0 whenr = 0. The second is the value of Y(0), which determines the concentration of the resulting densitymodel for a given s . Since we are only concerned with Y(0)=s2, which will control the shape ofthe density distribution, and notY(0) or s independently, which determine the absolute scale of thesystem, we only consider Y(0)=s2 as a single combined parameter.We use the fourth-order Runge-Kutta method (RK4) [Kutta, 1901] to solve this equation nu-merically. The step size varies from 0:01r0 for the most concentrated (small c) to 10−4r0 for themost extended (large c) models in our grid. The truncation error is less than 0:1% for all modelscalculated. A solution to this equation for a given set of boundary conditions gives Y(r). This canthen be fed into Equation 2.7 to go from r(Y) to r(r).Finally, to transform from a three dimensional density distribution to a projected density distri-bution, one must perform an Abel Transformation:S(R) =∫ ¥Rr(r)rdr√r2−R2 : (2.9)A range of the resulting projected density distributions is shown in Figure 2.1. The brief sum-mary in this section is not intended to be comprehensive, only to allow the reader to replicate ourprocedure for calculating these models. A more thorough description, including motivation for thedistribution function from which these models are derived, can be found in King [1966] or BT87.2.3.2 Luminosity vs. Star CountsIn the interest of presenting a straightforward comparison to previously determined values of Rc wefirst consider all of the stars in the catalogue for each cluster. Our fitting method will be discussedin Section 2.3.4. The model grid described in the section above is fit to the stars of each clusterin three different ways. In the first case, each star is weighted equally. The King radius r0 andconcentration c are determined for each cluster from the star count distribution. The value of the14Figure 2.1: Normalized surface density profiles. We have recalculated Figure 4-9b fromBT87, showing the projected density distribution for various boundary conditions.These distributions correspond to (from left to right) Y(0)=s2 = 12;9;6;3 or c =2:74;2:12;1:26;0:67, since there is a one-to-one relation between Y(0)=s2 and c.core radius Rc can easily be calculated from each best fitting model as defined in Equation 2.1. Acomparison between the values currently in the literature and our resulting best fit values and theiruncertainties is shown in the left hand panel of Figure 2.2. Our uncertainties are estimated withan iterative fitting procedure that considers uncertainty in the location of the cluster center as wellas uncertainty due to sample size through bootstrapping. We find that by fitting to the cumulativestar count distribution and assuming a conservative 20% uncertainty in Harris’ values the meandifference between our fit values and Harris’ values is inconsistent with zero at slightly more than4s , where the error used is the standard error. We tend to measure statistically significant largercore radii. This is not surprising, given that the Rc values reported in Harris [1996] are determinedfrom fitting the luminosity profiles of clusters.If we repeat our fitting procedure, but now weigh each star’s contribution to the cumulativedistribution by its luminosity, we then produce the result shown in the centre panel of Figure 2.2.These core radii now share the same bias as those in the Harris catalogue and the mean differencebetween our values and those of Harris is within 1:2s of zero.The explanation for this discrepancy is as follows. The mass-to-light ratio is a function of radiusin most globular clusters due to mass segregation; more massive (and therefore more luminous)stars will be more centrally concentrated than less massive (less luminous) stars. This implies that15Figure 2.2: The left hand panel shows core radius values taken from Harris [1996] plottedagainst values determined using our approach of fitting the cumulative star count distri-bution. The center panel shows values determined by fitting the cumulative distributionin which objects are weighted by their luminosities. The right panel shows fits to themass weighted cumulative distributions. The red line highlights a non-biased relation ineach case. The horizontal error bars represent 1s and the Harris values have no errorplotted.the mass-to-light ratio decreases toward the centre of the star cluster and so converting betweena distribution measured from surface brightness to the underlying distribution of stars does notinvolve a constant factor. In fact, how this mass-to-light ratio changes with radius also differsbetween clusters. So, even if we consider mass-to-light ratio as a function of R, there is no universalcorrection to be applied.If we wish to determine the mass distribution of the cluster, then fitting the star count distributionstill has an inherent bias for the same reason discussed above. We make a rough approximation ofthe core radius of the true underlying mass distribution by fitting cumulative distributions in whichstars are weighted by mass. The masses are assigned by fitting isochrone models from Dotter et al.[2008a] to the cluster CMDs. Isochrones are models that predict the distribution of stars in colour-magnitude space at a single instant in time. Modeling the globular cluster CMD with an isochroneassumes that all of the stars have the same age, which they do. After weighting the stars by mass,the results are shown in the right panel of Figure 2.2. These values fall somewhere between theprevious two, but the key point is that they are inconsistent with those determined from fitting thesurface brightness profile. This indicates that the surface brightness profile is not an appropriateproxy for the mass surface density of a cluster. For the remainder of the chapter we will attempt toquantify these effects by analyzing how subgroups of varying mass are distributed in each cluster.2.3.3 Selecting Groups Along the Main-SequenceFor each cluster, we create ten independent groups along the main sequence. We begin by fittinga fiducial sequence to each cluster color-magnitude diagram (CMD) in F606W vs. (F606W −F814W ). The fiducial sequence is created using an iterative algorithm. A manual starting pointis selected at the base of the CMD. The algorithm then chooses a direction to move by using theAMOEBA algorithm from Press et al. [2007] to minimize the median of the absolute separation ofthe nearby points perpendicular to each proposed step. There is also an additional cost added based16on the angle between sequential steps to prevent the sequence from bending back on itself. Theresult is a fiducial that traces out the centre of the sequence. Stated another way, the algorithm stepsalong the density distribution of points in Figure 2.3 while trying to maintain the steepest possiblegradient perpendicular to the path it takes. Once this fiducial is established for each cluster, theturn-off is then defined as the magnitude that produces the most negative derivative (dcol=dmag) ofthe fiducial. Ten magnitude bins are created starting from 0.5 magnitudes below this turn-off andcovering eight magnitudes in total. The widths of the bins are determined by calculating the standarddeviation of the CMD in the colour direction from the fiducial. Only objects that fall within ±3sare included. While this approach does require an initial point to be selected manually on the CMD,it was found that various starting points would converge to the same sequence within approximatelyfive iterations. These initial burn-in points were removed and the fiducial was linearly interpolatedbackward from the later points. This effectively removes all human input, leading to a bin selectionmethod that is entirely automated and consistent across all of the clusters used in this study. Figure2.3 (right panel) illustrates the location of these bins on a CMD. These bins could also have beenarbitrarily assigned by a human without any serious change to the results. Using the algorithmdescribed above avoids having to do this for all 54 clusters.Figure 2.3: The algorithmically determined fiducial sequence and the regions for the ten se-lected bins are shown for the cluster 47 Tucanae (NGC 104).172.3.4 Fitting Models to SubgroupsFor each cluster, we fit our grid of King-Michie models to the cumulative radial distribution ofeach of our 10 sub-groups. We fit only out to projected radii of 100 arc-seconds, due to the sizelimitations of a single ACS field. Out to 100 arc-seconds from the centre (Rmax in Equation 2.11),each annulus is fully within the field. Beyond this, more and more of each successive annulus willfall outside of the field. It is possible to take this effect into account when fitting the models, butincluding these regions does not significantly improve our ability to constrain the best-fitting modelparameters, and so they are ignored in the interest of simplicity.To generate a model to fit to our data we need to specify two of the three parameters r0, rt, andc, since the third can always be determined from the other two. The tidal radius of the cluster in allcases will be many times the maximum radius that can be observed with our field, and we expectthis to be a constant for all subgroups within the cluster since escape velocity is independent ofobject mass. Therefore, instead of fitting for two parameters, we use the values of c and Rc fromHarris [1996] to calculate rt for each cluster and hold this value fixed.Before fitting the models to the data, we also need to take into account the completeness ofthe stars in each bin. To do this, we use the results of artificial star tests. These tests are done byinputting many thousands of artificial stars into the science images and cataloging them using thesame data reduction procedure as was used for the real stars. This allows us to determine how wellstars of a known magnitude can be found across the image. From the results of these test we cancalculate a function called I. The incompleteness I(R) is the probability that a star in one of ourdefined magnitude bins would be found at a given position in the field. These corrections are thenapplied to the model density distributions before fitting to the data.Given the projected density distribution S(R) (stars per solid angle) and data incompletenessI(R) (fraction of input stars recovered), we calculate the cumulative distribution asC(R) =∫ 2p0∫ R0 S(R′)I(R′)R′dR′dq∫ 2p0∫ Rmax0 S(R′)I(R′)R′dR′dq: (2.10)The angular component is trivial since neither of these components have any q dependence. Itcancels on the top and bottom so we haveC(R) =∫ R0 S(R′)I(R′)R′dR′∫ Rmax0 S(R′)I(R′)R′dR′: (2.11)Here I is only dependent on R since we are considering a group of objects that fall within anarrow magnitude range. The cumulative distribution is normalized by the integral over the entirerange of R. This means that C(0) = 0 andC(Rmax) = 1.In determining the best fit to the cumulative distribution of each bin, there are two major sourcesof uncertainty. The first is the location of the centre of the cluster. We use the centre coordinatesand uncertainties from Goldsbury et al. [2010], which were determined using this same data set.The second source of uncertainty is the sample size in each bin. Both of these are taken intoaccount by using a bootstrap fitting method as follows. In each iteration, a sample of the samesize as the full bin sample is chosen with replacement. The location of the centre to use is drawnfrom a two-dimensional Gaussian distribution centred at the value from Goldsbury et al. [2010]with the appropriate s . The cumulative radial distribution is calculated and the best-fitting modelis determined by minimizing the maximum separation between the real distribution and the model,18which is the Kolmogorov-Smirnov statisticD. In each iteration, we record a single value for the bestfitting r0. After many iterations we can then generate a histogram of fit r0 values. This distribution isroughly Gaussian and is broadened by both sources of error discussed above. We take the mean andstandard deviation of these values as our best fitting r0 and 1s uncertainty for each sub-group in eachcluster. The dominant source of error in most cases is the√N random uncertainty. However, in verywell populated bins the systematic uncertainty from the location of the cluster centre contributes asmuch to the total error as the Poisson counting uncertainty.2.3.5 Assigning Masses to Bins and Fitting Power LawsAfter iterating the fitting procedure described in Section 2.3.4 for all sub-groups in all clusters,we can show how the concentration of stars along the main-sequence changes with magnitude foreach cluster. However, this is not particularly useful for making comparisons among clusters, sincethe apparent magnitude of an object is a result of a number of factors: distance, extinction, mass,metallicity, and cluster age. In terms of dynamics, we would like to look at how the concentrationchanges with object mass. So, to assign a mass to each magnitude bin, we need to consider theeffects of the other four physical properties.To do this, we again use an iterative approach. In each iteration we draw values of the four pa-rameters listed above from Gaussian distributions. The parameters of these distributions are drawnfrom a number of sources. All of the distances, metallicities, and extinctions are taken from Harris[1996]. The vast majority of these have no reported uncertainties in the original references, and sowe assume conservative estimates: for distance we assume a standard deviation of 0.2 magnitudes;and for extinction and metallicity we assume 20% uncertainties. The ages and uncertainties aretaken from Dotter et al. [2010], with the exception of five clusters that were not listed in that paper.These are NGC 1851 and NGC 2808 [Koleva et al., 2008], NGC 5139 [Forbes and Bridges, 2010],NGC 6388 [Catelan et al., 2006], and NGC 6715 [Geisler et al., 2007].The distance, metallicity, and age are then used to define a stellar isochrone model from Dotteret al. [2008a]. For each magnitude bin, we assign masses to all of the objects using this model andthen take the mean mass of all objects in each bin. Along the main sequence, there is a one-to-one relation between mass and magnitude (complicated slightly by the uncertainties that we discusshere). This one-to-one relation means that the variation in the mass of the objects in each bin will besmaller than the differences in mass between bins. For each iteration, we record one mass value foreach bin. We then use these mass values along with our best fitting r0 value and uncertainties to fita two parameter power-law model to the relation between King radius and mass using a maximumlikelihood approach. The isochrone model grid does not agree well with the structure of most CMDsat the very low mass end, so we ignore all bins for which the mean mass is less than 0.2 M⊙. Wealso ignore points for which the best fitting r0 falls outside of our field, since in these cases, r0 isnot well constrained by our data. For each iteration, we save the likelihood results over the twopower-law parameters shown in Equation 2.12. This power-law model is not physically, but ratherempirically motivated, as it fits a wide range of clusters well with a small number of parameters:r0 = A(MM⊙)B: (2.12)In each iteration we find the likelihood across the two-dimensional space of A and B. After manyiterations as described above, we can average this stack of two-dimensional likelihood surfaces. This19amounts to a Monte Carlo integration over the likelihood of the four other parameters (distance,extinction, metallicity, and cluster age) with prior distributions on each parameter (the Gaussiansdescribed in the second paragraph of this section). This leaves the likelihood of only the two power-law parameters, but with the uncertainties in the other four propagated through. An example fit for47 Tuc is shown in Figure 2.4.Figure 2.4: The best-fitting power-law model and likelihood of fit parameters for 47 Tuc. Sim-ilar figures can be found in the online supplement for Goldsbury et al. [2013], showingpower laws fit to both r0 and Rc for all 54 clusters in the catalogue describe in 2.2. A andB are defined in equation 2.12.2.4 ResultsThe parameter A is a normalization that corresponds to the King radius of cluster members at 1M⊙. This value is given in projected coordinates, which are related to the true three-dimensionalKing radius by the distance to the cluster. The parameter B describes how much the concentrationof stars changes with their mass, and is independent of the total scale of the cluster. Although wechoose to fit our models for r0, we could just as easily parametrize them by Rc and fit power-lawsto this parameter in the same way as in Equation 2.4. For all 54 of the clusters included in thisanalysis, these two power-law parameters, as well as their uncertainties and the c2 values for the20best fit (the value given is the total c2, not per degree of freedom) are included in Table 2 for bothr0 and Rc as a function of mass. For any of these clusters, this empirical relation can then be usedto estimate the distribution of stars of a given mass. The values of A in this table are in units ofarc-seconds, and the parameters correspond to a power-law function in the form of Equation 2.12.The two right-most columns of the table are the parameters for fits to Rc(M) rather than r0(M).Although r0 is determined from the three dimensional density distribution, our values come fromfitting the projected density distribution and so we fit in units of angle rather than distance. Sincewe are fitting in projected space, we leave our results in these projected units. These results can beconverted to physical units by using the distance given in Harris [1996] or another source.Parameter A has units of angle and our measurement of it depends on our distance from thecluster. Parameter B, on the other hand, is unitless and would be measured the same regardless ofour distance from the cluster. This parameter is a quantitative descriptor of the mass segregationpresent in a cluster. We can use the apparent distance modulus and extinction from Harris [1996]to convert our fit values for A in arc-seconds (or pixels) to parsecs, which gives us a parameter todescribe the intrinsic size of the cluster. After doing this, we find that A and B correlate stronglywith each other such that larger clusters tend to have less mass segregation than small clusters. Bothof these quantities also correlate with the core relaxation time from Harris [1996]. Figure 2.5 showsthese correlations.Figure 2.5: Correlations between our fit parameters and the core relaxation time of the cluster.Here A has been scaled by the distance to each cluster so that it is in units of parsecs.Red points indicate core-collapsed clusters not used in calculating the correlation statis-tic. The insets show the histogram of measured Pearson correlation coefficients afterbootstrapping and taking into account uncertainties in the fit parameters.The distribution inset in each diagram indicates the Pearson correlation coefficient (R) distribu-tion determined through bootstrapping. A principal component analysis of this three-dimensionalspace indicates only a single significant component that links all three dimensions. This can be un-derstood as follows: large clusters have longer crossing times and therefore longer relaxation times.The amount of mass segregation present in the cluster should approach some limit near equipartitionas many relaxation times pass. Therefore, these larger and slower to relax clusters will exhibit lessmass segregation than those that are small and quick to relax.Lines of the following forms are fit to each correlation:21log10(A) = b1+S1 ∗B; (2.13)log10(tc) = b2+S2 ∗B; (2.14)log10(tc) = b3+S3 ∗ log10(A): (2.15)These fits are shown as blue dashed lines in Figure 2.5, and the values of the best fit parametersas well as their variances and covariances are listed in Table 1. As with the calculation of the corre-lation coefficients, core-collapsed clusters are excluded from these fits. Core-collapsed clusters areclusters for which the luminosity profile continues to increase all the way to the centre of the clus-ter. These profiles do not level off as shown in Figure 2.1. Classification of core-collapsed clustersis taken from Harris [1996]. These are derived from luminosity profiles and so this classificationpotentially suffers from the biases introduced by mass segregation that are discussed above. Of theeight clusters that are classified as core collapsed, we find that the projected density distributionmeasured from star counts is well modeled by just a power law for six of them. Two of the clus-ters (NGC 6541 and NGC 6723) show clear evidence of levelling off toward the centre and do notappear to have collapsed cores when the density is measured through star counts rather than lumi-nosity. Similar discrepancies have been found by Miocchi et al. [2013], where clusters are foundto show a central cusp in their surface brightness profile but no such feature in their stellar densityprofile. The authors attribute this to the presence of a few anomalously bright stars, but we believethis discrepancy could also be caused by the systematic effects of mass segregation. Despite the factthat clusters that appear to be core collapsed in surface brightness may not actually be so, we havestill excluded all of these clusters from the correlation analysis shown in Figure 2.5, since this clas-sification still carries a potential bias through the calculation of the relaxation time. Core collapsedclusters are defined as those that do not appear to show the surface brightness distribution level offtoward the centre of the cluster as seen in Figure 2.1. In these clusters, the power law continues allthe way to R= 0.Table 2.1: Fit Correlations (the blue dashed lines from Figure 2.5)Relation b S Var(b) Var(S) Cov(b,S)log10(A) vs. B 0.35 0.29 0.007 0.003 0.005log10(tc) vs. B 9.06 0.68 0.019 0.012 0.013log10(tc) vs. log10(A) 8.25 1.84 0.003 0.043 −0.0012.5 ConclusionsWe have shown that measuring the core radius of a cluster through the surface brightness introducesa bias that makes clusters appear more concentrated than the actual distribution of their main se-quence stars or total mass on the main sequence. This bias is explained by the presence of masssegregation within clusters. Mass segregation creates a mass-to-light ratio that changes with radius,and so there is no constant factor that can be used to convert between the density distribution of22a cluster as measured through the surface brightness profile and the true underlying stellar den-sity. Cluster model parameters measured in this way will be biased to the most massive stars in thecluster, and will not accurately describe the distribution of lower mass stars.We have attemped to quantify the amount of mass segregation in 54 clusters using public cata-logues [Sarajedini et al., 2007]. This was done by fitting King-Michie models to a series of narrowmass bins in each cluster, allowing us to measure how the King radius varies with mass. We foundthat a simple two-parameter power law fits these data well and we have presented the parametersof these fits for all of the clusters in Table 1. We also present the best-fitting power law and thelikelihood contours of fit parameters for each cluster in the supplementary material. We have shownthat a single concentration parameter does not accurately reflect the distribution of most stars withina cluster, and the models we have presented could be used to estimate the distribution of stars as afunction of mass within the clusters.We find that the exponent of our power-law model correlates strongly with both the core relax-ation time and the absolute size of the cluster. This suggests that these models could be used as aproxy for the underyling dynamical state of the cluster. However, the scatter in these correlationsis larger than can be attributed to the errors alone. In the case of the second and third correlations,this is partly attributable to relaxation time values having no associated uncertainty. However, thescatter still appears to be larger than expected in the first correlation, for which uncertainties havebeen properly considered. This potentially indicates the effects of other factors that can influencethe dynamical evolution of clusters, such as differences in the individual cluster’s Galactic orbits, orthe presence of intermediate mass black holes in some clusters.This work was performed in 2013. At which point in time I was unfamiliar with the unbinnedlikelihood. In retrospect, this technique would also work extremely well for this sort of analysis.Redoing this analysis with these new statistics would likely yield much better constraints.23Table 2.2: Fit Powerlaw Parametersr0 RcCluster A B c2 A B c2NGC0104 32.3± 0.7 −0.96±0.05 27.3 24.4± 0.6 −0.95±0.05 27.1NGC0288 64.9± 7.6 −0.67±0.19 2.5 48.5± 5.7 −0.58±0.19 2.0NGC1261 16.9± 0.4 −1.55±0.03 13.1 14.2± 0.5 −1.15±0.03 32.8NGC1851 15.6± 0.5 −1.74±0.07 50.9 12.0± 0.4 −1.61±0.07 59.1NGC2298 11.7± 0.5 −1.54±0.08 9.4 8.9± 0.4 −1.48±0.08 9.9NGC2808 17.0± 0.9 −1.48±0.11 37.7 11.9± 0.3 −1.60±0.02 38.8NGC3201 44.8± 2.7 −0.96±0.11 12.6 33.8± 2.0 −0.93±0.11 12.3NGC4147 8.3± 0.6 −1.11±0.10 35.3 6.3± 0.4 −1.09±0.07 34.3NGC4590 23.0± 0.7 −0.89±0.05 6.8 17.3± 0.6 −0.87±0.05 6.8NGC4833 38.7± 1.3 −0.61±0.05 7.7 29.0± 1.0 −0.59±0.05 7.5NGC5024 19.6± 0.4 −0.74±0.03 11.4 14.8± 0.3 −0.73±0.03 11.2NGC5053 127.9±40.6 −0.24±0.62 1.9 81.5±20.0 −0.39±0.37 1.9NGC5139 102.2± 4.9 −0.06±0.06 3.0 76.7± 3.6 −0.06±0.06 2.9NGC5272 24.4± 0.5 −0.88±0.04 20.3 18.4± 0.4 −0.87±0.04 20.2NGC5286 13.0± 0.4 −1.96±0.06 19.2 10.1± 0.3 −1.83±0.06 23.4NGC5466 76.6± 8.8 −0.31±0.15 3.8 56.4± 6.5 −0.29±0.15 3.2NGC5904 26.7± 0.7 −0.91±0.05 26.5 20.1± 0.5 −0.90±0.05 26.2NGC5927 27.4± 0.6 −1.00±0.05 7.7 20.6± 0.5 −0.98±0.05 7.1NGC5986 18.2± 0.5 −1.24±0.05 18.2 13.9± 0.4 −1.17±0.05 19.5NGC6093 12.6± 0.4 −1.17±0.05 61.8 9.5± 0.3 −1.13±0.05 61.0NGC6121 50.8± 3.2 −0.61±0.12 5.9 38.3± 2.4 −0.60±0.12 5.9NGC6144 40.4± 1.6 −0.49±0.06 10.6 30.5± 1.2 −0.49±0.06 10.5NGC6171 29.3± 1.3 −0.88±0.08 13.8 22.1± 1.0 −0.86±0.08 13.6NGC6205 35.3± 0.8 −0.65±0.05 13.6 26.5± 0.6 −0.64±0.05 13.3NGC6218 37.5± 1.8 −0.76±0.08 10.7 28.2± 1.4 −0.74±0.08 10.4NGC6254 30.7± 1.0 −0.64±0.04 12.4 23.1± 0.7 −0.62±0.04 12.2NGC6304 16.4± 0.7 −2.18±0.10 6.6 13.7± 0.3 −1.85±0.03 16.4NGC6341 16.9± 0.5 −1.12±0.04 26.4 12.8± 0.3 −1.10±0.04 26.7NGC6352 36.6± 3.8 −1.36±0.26 3.6 28.6± 1.7 −1.16±0.07 3.7NGC6362 70.3± 7.0 −0.68±0.19 8.9 52.1± 5.1 −0.60±0.19 6.5NGC6366 96.0±23.1 −0.57±0.41 2.1 109.3±11.4 −0.30±0.16 16.8NGC6388 13.7± 0.4 −2.39±0.04 183.3 12.3± 0.2 −1.84±0.03 79.0NGC6397 28.7± 2.2 −1.01±0.13 6.5 21.7± 1.7 −0.98±0.13 6.2NGC6441 12.9± 0.2 −2.18±0.04 544.3 10.6± 0.2 −1.82±0.04 261.1NGC6535 9.7± 1.2 −2.47±0.10 8.6 9.5± 1.1 −1.89±0.10 5.1NGC6541 14.3± 0.6 −1.44±0.07 47.7 10.9± 0.4 −1.41±0.07 47.9NGC6584 18.2± 0.7 −1.26±0.06 11.0 13.8± 0.5 −1.20±0.06 11.4NGC6624 10.6± 0.4 −1.82±0.07 13.0 8.0± 0.3 −1.80±0.07 12.3NGC6637 14.9± 0.5 −1.40±0.08 6.6 11.3± 0.4 −1.35±0.08 4.2NGC6652 7.2± 0.5 −2.63±0.16 10.3 5.6± 0.3 −2.51±0.15 6.4NGC6656 46.8± 1.4 −0.60±0.05 13.0 35.2± 1.1 −0.59±0.05 12.7NGC6681 13.0± 0.5 −1.08±0.05 16.8 9.8± 0.4 −1.06±0.05 16.9NGC6715 12.2± 0.5 −1.76±0.09 20.1 8.0± 0.4 −2.10±0.11 71.7NGC6717 12.1± 1.6 −1.72±0.19 13.4 9.2± 1.2 −1.66±0.19 12.9NGC6723 32.2± 1.1 −0.75±0.06 9.4 24.2± 0.8 −0.71±0.06 9.0NGC6752 28.8± 1.3 −0.99±0.07 58.1 21.8± 1.0 −0.99±0.07 57.7NGC6779 20.0± 0.6 −1.02±0.04 18.9 15.1± 0.4 −0.99±0.04 19.0NGC6809 77.5± 8.3 −0.57±0.16 6.3 57.4± 6.2 −0.50±0.16 5.2NGC6838 72.3±14.8 −1.47±0.47 1.1 54.5±11.1 −0.94±0.47 0.5NGC6934 11.2± 0.5 −1.35±0.05 33.9 8.5± 0.4 −1.31±0.05 33.9NGC6981 22.8± 0.8 −0.66±0.06 13.8 17.0± 0.6 −0.63±0.06 12.9NGC7078 17.1± 0.5 −0.98±0.04 163.0 12.9± 0.4 −0.97±0.04 161.9NGC7089 18.3± 0.6 −1.02±0.06 45.9 13.8± 0.4 −1.01±0.06 45.5NGC7099 16.1± 0.7 −1.04±0.05 80.4 12.2± 0.5 −1.03±0.05 79.724Chapter 3An Empirical Measure of the Rate ofWhite Dwarf Cooling in 47 Tucanae3.1 IntroductionThe majority of stars that form in our Galaxy will eventually end up as white dwarfs. Given enoughmass to fuse hydrogen and helium, but not enough mass to eventually form a neutron star or ablack hole, the range of initial masses for white dwarf stars encompasses much of the lower end ofthe initial mass function [Fontaine et al., 2001b]. The study of such objects can therefore reveal awealth of information about the population that formed them. A predictive model of white dwarfcooling can also be used to derive precise ages for stellar populations [Oswalt et al., 1996, Hansenet al., 2007a]. Additionally, white dwarfs themselves can serve as laboratories for constraining thephysics that determine their cooling rates.The aim of this chapter is to obtain an empirical measurement of the rate at which white dwarfscool. Theoretical work in this field began with Mestel’s 1952 paper, in which he laid out a simplemodel for the rate at which white dwarfs should lose energy. This model assumes a completelydegenerate and isothermal core surrounded by a thin non-degenerate envelope. Thermal motion ofthe non-degenerate ions releases energy from the core of the white dwarf, while the rate of energyloss from the star is determined by the opacity in the envelope [Mestel, 1952]. The end result is aprediction that the temperature of a white dwarf should decrease with time since formation, t, ast−0:4. This basic model has been gradually improved upon since then, with the introduction of moreand more physical considerations. Summaries of progress in this field can be found in D’Antona andMazzitelli [1990], Hansen [2004], and most recently Althaus et al. [2010]. Hansen also presents acomparison of some modern cooling models. However, comparisons of white dwarf cooling modelsto observable populations tend to focus on the cool end of the sequence, which is the most relevantfor dating. We present here a detailed comparison of a number of current models to measurablecooling rates up to temperatures of 30,000K.In this chapter, we will analyze a group of white dwarfs that have all evolved from the sameinitial population. These objects are part of the Galactic globular cluster 47 Tucanae (NGC 104).All of the objects in this cluster formed at the same time, or at least formed over a period of timethat is very small when compared with the current age of the cluster. The fact that our entiresample of white dwarfs comes from a single population allows us to draw a number of simplifying25assumptions.The first is that all of the white dwarfs are of a similar mass. Because all of the stars shouldhave formed around the same time, all of the objects on the main sequence have had roughly thesame amount of time to evolve. The population of the cluster also has a roughly uniform chemicalcomposition, and so the lifetime on the main sequence depends entirely on the mass of the star.Together, these imply that the variance in the masses of objects that are observed leaving the mainsequence should be small. Initial final mass relations [Kalirai et al., 2009] then suggest that variancein the mass of the white dwarfs, which these turn-off stars will eventually become, is also small.This agrees with the analysis of white dwarf spectra in globular clusters, which suggest that all ofthe objects are of similar mass [Kalirai et al., 2009, Moehler et al., 2004]. It is important to notethat this assumption will not be true for very old white dwarfs, but since our entire sample shouldbe younger than a few Gyrs, this assumption is valid.The second assumption is that objects arrive on the white dwarf sequence at a constant rate.This follows from the assumption that the rate at which objects are leaving the main sequence isconstant. And the time spent in post main sequence evolution for these stars is also a constant. Thisis supported by a comparison of the luminosity function of the lower giant branch with a modelluminosity function with a constant evolution rate (see Section 3.5). Given that this rate is constant,and all objects leaving the main sequence will eventually become white dwarfs, the rate at whichobjects are arriving onto the white dwarf sequence must also be constant.These two points allow us to put all of our white dwarfs on the same cooling sequence by fittingmodels to each object’s spectral energy distribution (SED) to determine a temperature, and thenrequiring that the white dwarfs are uniformly distributed in time along the sequence. This approachallows us to measure the rate at which the temperature of the white dwarfs is changing in a way thatis independent of any white dwarf cooling models.3.2 DataOur data set consists of 121 orbits of Hubble Space Telescope (HST) observations. During theseorbits, exposures were taken simultaneously using the Advanced Camera for Surveys (ACS) andWide Field Camera 3 (WFC3). The current chapter focuses on the WFC3 UVIS data in filtersF390W and F606W. These data are comprised of 26 orbits covering 13 adjacent fields (1 orbit ineach filter). The description of the methodology in this chapter will refer to that data set. Similarmethods were used to produce a cooling curve from the ACS data, which is taken in filters F606Wand F814W. The WFC3 data cover over 60 arcmin2 over a range of 6.5–17.9 pc from the core of47 Tuc (in projection). The observations and data reduction are described in detail in Kalirai et al.[2012]. The footprint of these observations can be seen in Figure 3.1.3.3 Selecting our White Dwarf SampleSelection began with a catalogue generated from the final drizzled images as described in Kaliraiet al. [2012]. Drizzling is a method of combining many overlapping images of a single field into asingle pixel grid. This is done by rotating and shifting to align all of the images and then “drizzling”the signal from the original pixels down onto a new higher resolution grid. The signal from eachpixel of the input images is fractionally distributed to the new pixels on which it would fall. Whenthis process is done with many images, a combined image that is higher resolution than any of the26Figure 3.1: The footprint of the ACS (central star shaped region) and WFC3 (surroundingswath) observations.individual input images can be achieved. This is the case because the cameras on Hubble are under-sampling compared to the limits imposed by the optics of the telescope. From the drizzled image,the point-spread-function (PSF) was determined using DAOPHOT II [Stetson, 1987] (a standard bitof publicly available software for determining the point-spread-function of an image). The point-spread-function is model of the charge distribution that results on the detector from a point sourceof light. This includes effects from the optics and the detector. Root-mean-square (RMS) residualsof the PSF fit, SHARP (defined as the ratio between the best fitting two-dimensional Gaussian andthe best fitting two dimensional delta-function) from ALLSTAR [Stetson, 1994], and photometricerrors from DAOPHOT were used to cut down the full catalogue to a higher quality sample. Starsthat made the cut were required to have errors of less than 0.1 magnitudes in F390W and F606W, aSHARP value within 2s of the mean SHARP value, and an RMS of the PSF fit less than 2.5%. TheCMD after selection cuts are made is shown in Figure 3.2.The region around the white dwarf sequence between roughly 22nd and 27th magnitudes inF606W was bounded manually. The objects inside this region were used as our final sample. Al-though some real objects are almost certainly lost due to these cuts, this can be properly accountedfor during the analysis of artificial star tests. The same cuts are used to determine which objects are“found” from our input artificial stars.Uncertainty also arises from manually bounding our sequence. We need to restrict the regionin which we consider objects to be white dwarfs to avoid including objects from the SMC in ourwhite dwarf sample. This restrictive region is also likely to exclude a small number of stars thatshould end up in our sample, but were too far from the centre of the white dwarf sequence (due27Figure 3.2: The CMD after making selection cuts in error, sharpness, and c2 of the PSF fit.The final WD selection region is shown in the enclosed area (red). 887 objects are foundwithin this enclosed region.to photometric error) to be included in the chosen region. To account for this loss of stars, weanalyze the distribution in colour along the white dwarf sequence and model the spread in colouras Gaussian. If we integrate this Gaussian over our enclosed area and compare this to the totalintegral, we find that we expect to lose 14 objects due to this bounding region. This means thatof our final sample of 887, each object should be weighted as (14+887)/887 objects due to thewhite dwarfs excluded from our selection region (prior to incompleteness corrections). This factoris approximately a 1.5% correction, and is included in the error of the time scaling of our object list.This scaling is discussed in more detail in Section 3.5.3.4 Fitting White Dwarf TemperaturesTo determine an effective surface temperature (Teff) for each object in our sample, we utilize spec-tral models from Tremblay et al. [2011a]. The models are for thick atmosphere (the H envelopemass fraction is approximately qH = 10−4) DA white dwarfs (white dwarfs with pure hydrogen at-mospheres), and are parametrized over a grid of surface gravity and Teff values ranging from 6.0 to10.0 and 6,000 K to 120,000 K, respectively. The DB (white dwarfs with pure helium atmospheres)28fraction is assumed to be negligible, which is supported by analysis from Woodley et al. [2012b].These models (Fmod) are then combined with the filter throughput curves (Sl ), which are availableonline from STScI [2009], to generate a magnitude in a particular filter for a given set of modelparameters:Mmod =−2:5log10(∫ ¥0 lFmodElSl dl∫ ¥0 lF0Sl dl)+5log10(dR): (3.1)Fmod is the spectral model from Tremblay et al. [2011a]. This is in units of power per wavelength(l ), and depends on Teff and MWD (the mass of the white dwarf). El is the interstellar extinctionmodel. It describes the fraction of power transmitted between the object and observer as a functionof l . The model comes from Fitzpatrick [1999b]. The Fitzpatrick reddening curve is parametrizedby E(B−V ), and RV . These are defined as follows. AV is the extinction inV and AB is the extinctionin B (the fraction of power lost in filter V and filter B respectively). E(B−V ) = AB−AV ; andRV = AVAB−AV . The symbols B and V refer to magnitudes in the B and V filters (centred at 445nm and551nm) in the photometric standard system of filters Johnson and Morgan [1953]. We fix RV = 3:1which is found to apply in most cases [Fitzpatrick, 1999a]. This value also gives good agreementbetween our white dwarf model sequence and our data, suggesting that it is appropriate for 47 Tucas well.The parameter d is the distance to the cluster, R is the white dwarf radius. The spectral energymodel Fmod gives the power at the surface of the star, so the last term in the equation appropriatelyscales this to an observer at distance d.The model magnitudes are fit to our sample of 887 white dwarfs using a maximum likelihoodfitting method. The likelihood of the data given the model is:L(data | parameters) =2Õi=1exp[−(datai−modi)22s2i](3.2)The product here is over the two filters in which our magnitudes are measured. Here the modelparameters are (m−M)0, E(B−V ), MWD, and Teff. (m−M)0 is the true distance modulus. Itis the distance in magnitude between the inherent luminosity and observed flux of an object dueonly to the geometric distance. This parameter has a one-to-one relation to d. The parameter Rfrom Equation 3.1 does not need to be specified as it is uniquely determined by MWD and Teff. Inthe white dwarf model grid these three parameters are related so that specifying any two of themuniquely determines the third.To condense this to the likelihood of Teff given our data, we integrate L(data | parameters) over(m−M)0, E(B−V ), and MWD with Gaussian priors of 13:30±0:15, 0:04±0:02, and 0:53±0:02M⊙ respectively [Woodley et al., 2012b]. The prior distribution over distance modulus used here isthe mean of all values quoted in Woodley et al. [2012b], which do not depend on the white dwarfsequence. The maximum likelihood in Teff is found for each object. The uncertainty in the distanceis twice as large as the value quoted by Woodley et al. [2012b]. This more accurately reflects thedistribution of distance measurements made by various groups, as presented in that paper. The sivalues in Equation 3.2 are the estimated photometric errors of each data point. These values arederived from artificial star tests presented in Kalirai et al. [2012]. The error distributions of the twoWFC3 filters are shown in Figure 3.4.29Figure 3.3: AmodelWD spectrumwith an effective surface temperature of 10,000K and a sur-face gravity of log(g)=7.5. Throughput curves for our two filters are shown overplottedon a separate axis scale on the right. The spectrum shown is the result of 3D hydrody-namic simulations of white dwarf atmospheres. Qualitatively, it is a blackbody spectrumwith hydrogen absorption features corresponding to the Balmer series transitions (n=2).The continuum suppression to the short wavelength side of 364nm is referred to as theBalmer Jump in stellar astrophysics. The underlying blackbody spectrum would followthe continuum seen to the right of this jump. The suppression to the left is due to contin-uum absorption. Photons of higher energy than 364nm can continue to ionize hydrogenatoms in the n=2 state, and so the continuum is suppressed beyond this limit.Figure 3.4: The standard deviation of input – output magnitude derived from artificial startests described in Kalirai et al. [2012] (for all fields combined).303.5 Determining the Timescale3.5.1 A Sorted List of TemperaturesAfter an effective surface temperature is determined from fitting spectral models to the photometryof each object, this list is sorted from hottest to coolest. This sorted list is the first step towarddetermining the rate at which the white dwarfs are cooling. Because the objects are arriving on thewhite dwarf sequence at a constant rate, the position in the sorted list should be proportional to thetime spent on the sequence. That is to say, the time between when the 100th and 200th objects arriveon the sequence should be the same as the time between the 200th and 300th.Completeness of the sample also needs to be taken into account. The cooler objects are alsofainter and thus less likely to be detected in a crowded field. Using artificial star tests as describedin Kalirai et al. [2012], we can assign a completeness correction to each object by first binning theartificial star data by magnitude and position. The completeness correction or fraction in some binis thenCcorr =NinputNfoundor Cfrac =NfoundNinput: (3.3)This can be interpolated to the magnitude of each object in our sample and the distance fromthe cluster center. To calculate the position of each object in a completeness corrected sorted list weusex(n) =nåi=0Ccorr(i): (3.4)Each temperature corresponds to one white dwarf in our sample, and each white dwarf hasan assigned completeness fraction dependent on it’s magnitude and location in the field. This isdetermined by finding the fraction of artificial stars that are recovered by our detection procedurewith similar magnitudes and positions in the field. This is Cfrac in the above equation. The Teffvalues are sorted from hottest to coolest. The completeness correction values are put in the sameorder. The completeness correction values are the weights that should be assigned to each starCcorr.For example, if we believe that our methods have only recovered 50% of stars similar to a group weare considering, then those stars should be weighted as 2 stars each to account for the unobservedobjects. If the first two stars in our sorted list each have completeness fractions of 50%, then the firsttwo positions in the completeness corrected list would be 2 and 4 from Equation 3.4. In practice, thecompleteness corrections are the smallest for the hottest stars (brightest) and largest for the cooleststars (faintest).If the entire sample were 100% complete, then this would be an array of positions in the sortedlist (1, 2, 3, etc.). In our sample the completeness corrections are small. Even the faintest objectsare more than 75% complete, and overall the sample is 89% complete. The completeness fractionas a function of magnitude is shown in Figure 3.5. The artificial star tests were actually analyzedfield by field. Fields closer to the centre of the cluster are obviously more dense, and have largerincompleteness. The results in Figure 3.5 are the total for all fields. Incompleteness was alsoconsidered in both filters. Finding a star is always limited by the filter in which it is fainter comparedto the surrounding objects. So, it was found to make no difference whether incompleteness valuesfrom a single filter, or a combination of both filters, was used. Since the corrections were so small,31this had no effect on our final fitted parameters. Completeness corrections for the ACS data werecalculated in the same way. The entire ACS sample is greater than 90% complete.Figure 3.5: Completeness fraction of our WD sample as a function of F606W magnitude (av-eraged over all fields).A sorted, completeness corrected list is shown in log-log space in Figure 3.6. This relationis well fit with a broken power law. The slopes and “break” of such a fit are also shown in thefigure. These parameters do not depend at all on the absolute scaling of the x-axis, since thisonly acts as an additive constant in log-space (log10 [∆t=object]). These power-law parameters andtheir uncertainties were derived by iterating our temperature fitting analysis. In each iteration weassume parameters of (m−M0, E(B−V ), and MWD with values drawn from the prior distributionsdiscussed above. After each iteration, the temperatures are sorted and fit with this broken powerlaw. The distributions of these fit parameters define the best fitting value, and the uncertainties. Thisprocess was also repeated for the ACS field with 292 WDs. The results are shown in the lowerpanel of Figure 3.6. All fitted parameters are consistent with the parameters from the WFC3 data.The uncertainties in these fit parameters take into account the uncertainty in the photometry of ourindividual objects, as well as the uncertainties in all of our marginalized parameters.3.5.2 The Rate at Which Stars are Leaving the Main SequenceIn order to put the completeness corrected list into units of time, so that our cooling curve can becompared with theoretical models in both dimensions, we need to find the rate at which objects arearriving onto the white dwarf sequence in our field. This rate should be the same as the rate at whichobjects are leaving the main sequence. One could imagine that dynamical relaxation of white dwarfscould cause problems with this assumption. The white dwarfs are considerably less massive thantheir giant branch progenitors and should be expected to relax on timescales similar to the relaxationtime for the field. However, this relaxation time is long. At the half-light radius the relaxation timeis 3:5 Gyrs [Harris, 1996]. Using the velocity dispersion profile of 47 Tuc calculated by Lane et al.[2010], together the footprint of our observations, and the density profile from our own data, we32Figure 3.6: The fitted Teff as a function of object number for our WDs. The object numberis proportional to the time on the white dwarf sequence, and is calculated as shown inequation 3.4. The y-axis is just the temperature fit to each object’s SED, as discussed insection 3.4. The results from both data sets agree well within uncertainties. Residualsare shown as Data−ModelData (the fractional residual).33calculate an expected mean relaxation time for our WFC3 field to be 10 Gyrs. We can compare thisto the expected time since the white dwarfs in our sample were the mass of current turn-off stars,which is < 2:5 Gyrs for stars brighter than 27.8 in F390W (utilizing models from Fontaine et al.[2001a]). Given that these objects have only existed for a small fraction of a relaxation time at theircurrent mass, there should be no significant mass segregation. Additionally, comparing the radialdistributions of the white dwarfs and the giant branch stars, no significant difference is apparent, asshown in Figure 3.7. This indicates that the white dwarfs have not yet relaxed, and their distributionin the field should be the same as their progenitors. To determine the age scaling for our coolingcurve we need to analyze the giant branch of the CMD. The upper giant branch is saturated in ourexposures, so we use the catalogue from Stetson [2000] of an overlapping field in B and V, selectingand completeness correcting only those objects that would be found within the region of our WFC3observations.Figure 3.7: A comparison of the cumulative radial distributions of the entire WD sample andthe lower RGB stars. A KS test indicates that there is a 73% chance that two samples ofthis size drawn from the same underlying distribution would differ by at least this much.We therefore find no evidence of significantly different radial distributions for the WDand RGB stars.On this portion of the CMD, stars evolve almost vertically as their outer layers expand accom-panied by a very moderate decrease in effective surface temperature. To determine an evolution ratefrom these objects, we use theoretical models from Dotter et al. [2008a]. Our input parameters forthe model are [Fe/H], (m−M)0, current cluster age (A), and an age range (∆A). This age rangecorresponds to the time it takes an object to evolve from one point on the lower giant branch toanother point (point P1 to point P2 in Figure 3.8). In our case, all of these parameters except theage range are well determined. We calculate the likelihood of our magnitude range resulting fromthis model over this whole parameter space using Equation 3.5. We marginalize out all parametersexcept the age range, using priors that describe the measured values of [Fe/H], (m−M)0, and A:L(∆Mr | parameters) = exp[−(∆Mm−∆Mr)22s2]: (3.5)34Here ∆Mr is the selected magnitude range (shown in Figure 3.8), while ∆Mm is the range pre-dicted by the model given a particular choice of input parameters. Chaboyer and Krauss [2002]show that there is an uncertainty of order 3% in the age of such models. To first order, object ageand magnitude are linearly related in this region, and so this translates to an uncertainty of 3% in theoutput magnitude range. This uncertainty is used as s in Equation 3.5. The priors used in Equation3.5 are all Gaussian with the following means and standard deviations: [Fe=H] = −0:75± 0:08,(m−M)F606W = 13:42±0:15, and A= 11:5±0:75 Gyr [McWilliam and Bernstein, 2008, Woodleyet al., 2012b, Dotter et al., 2009].Figure 3.8: The CMD from Stetson [2000], including only objects that fall within the regionof our WFC3 data. The region of the CMD between points P1 and P2 is fit with modelsfrom Dotter et al. [2008a] to determine the time it takes an object to evolve between thepoints.The result is a likelihood distribution of ∆A for our chosen magnitude range that takes intoconsideration uncertainties in the input model parameters, as well as uncertainty in the model itself.This distribution is shown in Figure 3.9.Our model prediction for the timespan corresponding to this chosen magnitude range is 496+95−83Myrs (2s uncertainties). The uncertainties on this value are entirely due to uncertainties in themodel as well as the input parameters. To determine the rate at which objects are evolving off the35Figure 3.9: Likelihood of ∆A fit to our selected magnitude range, after marginalizing out allother model parameters. ∆A corresponds to the time it takes an object to evolve frompoint P1 to point P2 in Figure 3.8.main sequence, the fit time interval must be divided by the completeness corrected number of objectsin this region of the CMD (shown in Figure 3.8) within our field. The relative uncertainty in thenumber of objects must be 1√N(where N is the number of objects before completeness correctionsare applied), resulting in a final rate of 1.7 ± 0.2 Myrs per object, an age scaling uncertainty of12%. This value has no physical meaning for the cluster as a whole. It is dependent entirely onour field. However, as this is the same field as our WD sample, it is appropriate to assume that thisrate is the same as the rate at which objects are arriving on the WD sequence. This time/object isthen multiplied by our ordered and completeness corrected list, resulting in an empirically measuredrelation between the temperature of the white dwarfs, and their age on the sequence. This methodalso assumes that the rate of evolution off the main sequence is constant. This assumption is wellsupported by analyzing the cumulative distribution of magnitudes between the two points in Figure3.8. This distribution is shown compared against a model distribution assuming a constant evolutionrate in Figure 3.10.The distribution of stars on the giant branch shown in Figure 3.10 strongly suggests that theformation rate appears constant over the last few hundred Myrs. However, we are interested in starsthat became white dwarfs up to 2 Gyrs ago. To derive the formation rate over time we can combinetwo terms. The first isdNdMµM−a−1; (3.6)which is the derivative of the Salpeter mass function [Salpeter, 1955]. The second isdMdtµM4: (3.7)This is the derivative of the relation between the stellar lifetime and mass from Duric [2004].36Figure 3.10: A cumulative distribution of Stetson’s RGB magnitudes between P1 and P2 inFigure 3.8 compared to the expected distribution assuming a constant evolution rate.A KS test suggests that there is a 40% chance that a sample drawn from the modeldistribution would differ by at least as much as this data.Combining these equations gives the formation rate:dNdt=dNdMdMdtµM3−a : (3.8)So the mass function needs to be quite steep (a = 3) in order for the assumption of a constantformation rate to hold over long times. Assuming the standard value of a = 2:3 and using therelation between mass and time from Duric [2004] leads todNdtµ t−0:233: (3.9)The time in this equation is actually the main sequence age of the stars. On the cooling curve,the stars that are at time zero are stars that just became white dwarfs and so t = t0 where t0 is the ageof the clsuter. Stars that are at an age of 109 yrs in Figure 3.13 are stars that became white dwarfs1 Gyr ago. So, to parametrize the formation rate in terms of the time on the white dwarf coolingsequence we havedNdtµ (t0− t)−0:233; (3.10)where now t is the time on the white dwarf cooling sequence. The current formation rate isderived to be 1.7 Myrs per object, so normalizing the formation history based on this leads todNdt=11:7Myrs(t0− tt0)−0:233: (3.11)Equation 3.11 is the historical formation rate of white dwarfs in this region of the cluster as-suming the Salpeter IMF. The difference between this calculated rate and the assumed constant rate37over the 2 Gyr range of our data is shown in Figure 3.11.Figure 3.11: The assumed constant formation rate is shown in red. The formation rate shownin black is the same as derived in Equation 3.11.Integrating this rate over time and rearranging for t(N) gives the times that should be assigned toeach numbered white dwarf in the sorted list shown in Figure 3.6. The difference in these integratedrates is shown in Figure 3.12.Figure 3.12: The cumulative number of white dwarfs formed in our field over the given timeis shown. The red line results from assuming a constant rate of formation. The blackline results from the assumptions of a Salpeter IMF outlined above.Overall, these more complicated assumptions lead to assigned ages to our white dwarfs thatare 1:8% younger by the end of the list. This approach does not lead to any greater agreement38between the data and models described in the following sections. In practice, this assumption leadsto negligibly better agreement for two of the four models tested, and negligibly worse agreementfor the remaining two models. In the interest of simplicity, and to avoid bringing additional modelparameters into our analysis we consider the white dwarfs as forming at a constant rate over the lasttwo Gyrs. The analysis above shows that this more thorough model of their formation history is notnecessary.3.6 The Cooling CurveThe cooling curve determined from our data is shown in Figure 3.13. The shape of the coolingcurve is well constrained, as can be seen in the residuals of Figure 3.6. The uncertainties quotedthere represent the 1s uncertainties in the upper and lower slopes taking into consideration both therandom uncertainties in the magnitudes of each object, as well as the uncertainty in the magnitudecalibration [Kalirai et al., 2012]. The contours shown in Figure 3.13 display the systematic uncer-tainty in the x-direction that arises from our uncertain time-scaling factor. This is the dominantsource of uncertainty in our time scaled cooling curve.3.6.1 The ModelsWe compare our data to models from four groups. The models of Fontaine et al. [2001a] havepure C cores with H and He envelope mass fractions of qH = 10−4 and qHe = 10−2. The Hansenet al. [2007a] models also have envelope mass fractions of qH = 10−4 and qHe = 10−2, and a mixedC/O core with profiles from Salaris et al. [1997]. Models from Lawlor and MacDonald [2006b] areevolved from the main sequence. By the time they reach the white dwarf stage, the envelope massfractions are qH = 9:5×10−5 and qHe = 3:1×10−2. The core is uniformly mixed with C and Omassfractions of 0.1 and 0.9 respectively. Lastly, we use the models of Renedo et al. [2010], specificallythose made publicly available with Z=0.001. The envelope thickness of these models varies withtime and is shown in Figure 4 of Renedo et al. [2010]. The core composition is discussed in section3.6 of that same paper with the C/O mass fractions being roughly 0.3 and 0.7.3.6.2 Model ComparisonAs can be seen in Figure 3.13, the general shape of the empirical cooling curve agrees reasonablywell with all four models. The models all show a linear region with a slope of approximately−0:20followed by a transition region, and then a section that agrees with Mestel’s value of −0:4. Fromour simple broken power law fit we find that the centre of our power-law break occurs at an age of(104 ± 26) Myrs, and a temperature of (19900 ± 1600) K. Fitting the same broken power law tothe models yields transition ages of 190 Myrs (Fontaine et al.), 129 Myrs (Lawlor & MacDonald),107 Myrs (Hansen), and 250 Myrs (Renedo et al.). This makes the time of our power law transitionconsistent with both Lawlor & MacDonald and Hansen, within 2s , but inconsistent with Fontaineet al. and Renedo et al. at greater than 3s . However, the temperatures at which this transitionoccurs in the models seems to differ from our data. These are 14,100K (Fontaine et al.), 16,600K(Lawlor & MacDonald), 16,200K (Hansen), and 13,500K (Renedo et al.). These differ at around2.5s for Lawlor & MacDonald and Hansen, and close to 4s for Fontaine et al. and Renedo et al.It is difficult to draw conclusions about how well these models fit the data by condensing our entire39Figure 3.13: The measured 1 and 2 s uncertainties of our cooling curve (the shaded blueregions). The shape in this region is well constrained as shown in Figure 3.6, but theage scaling leads to a large systematic uncertainty in the x-direction. This is the regionrepresented by these contours. Models of white dwarf cooling rates from Fontaine et al.[2001a], Lawlor and MacDonald [2006b], Hansen et al. [2007a], and Renedo et al.[2010] are shown for comparison.cooling curve down into a few statistical quantities. A more robust method involves comparing ourobserved cumulative distribution of temperatures to those that would be predicted by the models.In such an analysis, the scale of the x-axis is entirely irrelevant since only the relative time spent ata given temperature compared to any other temperature will determine the shape of the cumulativedistribution. This approach also compares the distribution over the entire temperature range, ratherthan condensing our curve down to a few statistics, such as the location of the power law break.From these cumulative distributions we can calculate the Kolmogorov-Smirnov (KS) separationstatistic “D” [Kolmogorov, 1933], indicating the largest separation between the cumulative distri-bution of our data and that of the model. Next, we draw samples of 813 objects from each modeldistribution 500,000 times. Here we use a sample size of 813 rather than our full sample of 887,because only 813 of our objects have temperatures that fall between 7,000 K and 30,000 K. Eachtime a random sample of temperatures is drawn from the model, each temperature is given a ran-dom shift corresponding to our measured temperature fit error distribution. All of the temperaturesin each sample are also given the same shift corresponding to our systematic temperature error re-sulting from an uncertain distance to the cluster. The uncertainty in distance leads to a systematicuncertainty in all of our absolute magnitudes and therefore a systematic uncertainty in all of ourfit temperatures. The errors added to the model samples correspond to a 1s distance modulus un-certainty of 0.15. This is a conservative estimate, and is actually 50% greater than the standarddeviation of the values quoted in Woodley et al. [2012b]. In each iteration we calculate the valueof D between the sample and the model it was drawn from. We then build a distribution of theseD values. This distribution is the range of D values expected for a sample the size of our data thatis drawn from each model. It is then possible to ask how likely it would be to draw a sample from40Figure 3.14: This plot shows measured cumulative distribution of our objects found to havetemperatures between 7,000 K and 30,000 K (813 objects WFC3, 188 objects ACS), aswell as those predicted by the four models shown in Figure 3.13.Figure 3.15: Cumulative distribution of our fit temperatures compared to the models whenassuming values at the extreme ends of literature distance measurements for 47 Tuc.41these models that produce a value of D greater than or equal to that of our data. This is done byintegrating the calculated distribution from the real data value of D to infinity. These distributionsare shown in Figures 3.16 and 3.17.Figure 3.16: Results from our iterative error analysis using data fromWFC3. Each distributionis a histogram of D values found when repeatedly comparing samples from each modelback to the theoretical model itself. The blue line indicates the value of D calculatedwhen comparing our real data to each model. Integrating from the blue line to infinitygives the probability that data drawn from each model would differ by at least as muchas our real data does. All of these cooling models are inconsistent with our data at the 3slevel. This is true even after considering a conservatively large systematic uncertaintyin distance and a random uncertainty in temperature.Even with this added uncertainty, all of the cooling models are ruled out as being consistent ourdata at around 3s . Hansen and Renedo et al. are consistent at roughly 3s , Lawlor & MacDonaldat around 3.5s , and Fontaine et al. at a bit less than 3s . This is due largely to the sharpness ofthe transition in our cooling curve, and the temperature at which the transition occurs. None ofthe models show such an abrupt transition to the Mestel slope, and so the predicted cumulative42Figure 3.17: Results from our iterative error analysis using data from the ACS field. Theseshow the same thing as figure 3.16, but using the ACS data rather than the WFC3 data.distributions differ significantly from the observations. It is important to note that the age scaledetermined from the giant branch in Section 3.5.2 can only add further constraint. The cumulativedistributions are normalized to run between 0 and 1 on the y-axis before comparisons are made. Theage scale determined from the giants will not affect the cumulative distributions presented above inany way, and cannot possibly make the data more consistent with the models.We also used other cooling models available from Bergeron et al. [2011]. We carry out the samecomparison described above with the mixed C/O thick atmosphere models, as well as the mixed C/Othin atmosphere models, in addition to the pure C thick atmosphere models used above. We foundthat, given our sample size and temperature range, we cannot distinguish between these models, andthey are all inconsistent with our data. Unfortunately, from this comparison we gain no insight intothe reasons for the discrepancy between the data and the models. It appears this difference cannotbe easily explained by varying the core composition or H atmosphere thickness.In addition, we consider the effects of varying the DB fraction on the expected cumulativedistribution of our temperatures and magnitudes (our original approach assumes our sample is 100%43DA white dwarfs). We find that we cannot make the data significantly more consistent with themodels of Fontaine et al. [2001a] by allowing the DB fraction to vary freely.A comparison of the temperature distribution from the ACS field to each model is shown inFigure 3.17. Even though the difference between the model distributions and ACS data is almostidentical to the difference between the WFC3 data and the models, due to the smaller sample size(188 vs. 813) all of the models are consistent with the data at around 2s .We also demonstrate in Figure 3.15 that the inconsistency between the models and our datacannot be explained in terms of a significantly different distance. We repeat our fitting procedurewhile assuming the two most extreme distance measurements from those summarized in Table 1 ofWoodley et al. [2012b]. We then build a cumulative distribution of fit temperatures and comparethis to the models. Even centring our distance prior 2s from our mean distance of 13.30, our resultsare inconsistent with any of the models at just below 3s .Finally, we perform a simple sanity check to show that the difference between the data and themodels is not a result of bad data in one or more of our filters. We calculate the expected cumulativedistribution of object magnitudes in each of our four filters for each model. We then compare thisto the observed cumulative distributions. These distributions are shown in Figure 3.18. All of thesedistributions show the data weighted toward brighter magnitudes than the models. This agrees withthe results we find when making the comparison in temperature space. We also show the comparisonto the luminosity function in F390W in Figure 3.19. We also observe a discrepancy between themodels and data in this space as the luminosity function is just the derivative of the cumulativedistribution. Or, more naturally perhaps, the cumulative distribution is the integral of the standardluminosity function as shown in Equation 3.12,C(M) =∫MM1F(M′)dM′∫M2M1 F(M′)dM′: (3.12)In this equation the cumulative distribution is calculated over the magnitude range [M1;M2] andnormalized to 1. The term F(M) is the luminosity function (dN/dM). Calculating the cumulativedistribution directly avoids the need to bin the data in either magnitude or temperature.3.7 DiscussionWe have presented an empirical measurement of the cooling rate of white dwarfs in the cluster 47Tucanae. Our results suggest two linear regions (in log-space) with a transition occuring around20,000K. In the upper region, between approximately 35,000K and 20,000K, the temperature of thewhite dwarfs changes like t−0:19. This slope agrees with all of the models within uncertainties. In thecooler region, between approximately 19,000K and 7,000K, the temperature decreases like t−0:42.This is close to Mestel’s original value of −0:40 and in fact agrees precisely with more thoroughmodern calculations, such as those summarized in Hansen [2004]. This region of the cooling curveagrees with all of the models as well. However, there is a clear discrepancy between the abrupttransition seen in the data and the more gradual transition seen in the models.Our measurements of both the white dwarf temperatures and their age on the sequence use arobust maximum likelihood fitting method and a thorough treatment of all the uncertainties involved.The disagreement of the available models with our data, as well as inconsistency between the modelsthemselves, suggests differences in the input physics involved in the transition region between the44Figure 3.18: Cumulative distributions of our object magnitudes in each filter compared withthose predicted by the models.two power-law regimes. The results of our analysis can therefore be used to constrain the inputphysics of these models.Such cooling models are used to predict the luminosity function of white dwarf populations ofvarious ages, and are used to date stellar populations, such as the Galactic disk [Oswalt et al., 1996],old globular clusters [Hansen et al., 2007a], and younger open clusters [von Hippel, 2005]. Theseresults therefore have potentially widescale implications for the age of the Galaxy and its variouscomponents.The unbinned method outlined in the next chapter could in principle provide different probabil-ities when comparing the data and models discussed in this chapter. However, it would not lead toany better agreement. The unbinned likelihood is a more strict test than the approach used in thischapter, and it is more difficult for data and models to match using that approach. Given that the dataand models disagree with the approach discussed in the previous chapter, it is highly likely that thisdiscrepancy would persists when using an even stricter analysis such as the unbinned likelihood.45Figure 3.19: The luminosity function of our WFC3 data compared to the models over theF390W range [ 21 , 26.5 ].46Chapter 4White Dwarf Cooling at HighTemperatures: The Effect of NeutrinoProduction on Cooling Rates4.1 IntroductionWhite dwarfs are often utilized as a way to indirectly test physical properties that are not directlyobservable themselves. White dwarfs have been used to measure the distances [Renzini et al., 1996,Zoccali et al., 2001, Woodley et al., 2012a], and ages [Hansen et al., 2007b, 2013] of globularclusters. They have also been used to constrain the initial mass function and age of the galaxy[Wood, 1992]. The physics of crystallization of dense white dwarf interiors has been indirectlyconstrained by Winget et al. [2009]. The white dwarf cooling sequence has even been used tomeasure dynamical relaxation in the globular cluster 47 Tucanae [Heyl et al., 2015]. White dwarfsin globular clusters in particular have proven useful as a testing ground for many types of models,because all of the objects in a cluster come from an initial population with roughly the same age,extinction, distance, and typically the same initial chemical composition.In this chapter we use the most well-populated white dwarf sequence ever measured in a globularcluster (47 Tuc) to indirectly constrain the thickness of the Hydrogen layer, and the rate of neutrinoproduction in the cores of white dwarfs. This will build on our previous work from Goldsbury et al.[2012], in which we constrained white dwarf cooling models from various groups against multi-band HST photometry in 47 Tuc. These models are again compared to the new data set in Section4.4.4.Our data set is obtained from nine HST orbits that cover the core of 47 Tuc in ultraviolet andvisible filters. In total more than 90% of the inner 5 arc-minutes of the cluster was observed, andover 3500 white dwarfs were detected. This data set is described in detail in Section 4.2. Theanalysis discussed in this chapter complements the work in Hansen et al. [2015]. In that paper,similar white dwarf cooling models were fit to the ACS data, shown as the surrounding annulus inFigure 4.1. These models were fit by binning the data in two dimensions (colour and magnitude)and using the c2 statistic to address the quality of the fit.Our analysis combines data in separate regions from both WFC3 and ACS images. To takefull advantage of the constraining power of this data set, we use an unbinned maximum likelihood47Table 4.1: Exposure times in seconds, for each of the nine visitsCamera Filter Exposure 1 Exposure 2WFC3 F225W 380 700F336W 720 485ACS F435W 290 690F555W 660 360method. The data occupy a three dimensional parameter space comprised of magnitudes in two fil-ters and a radial coordinate in the observed field. Theoretical models are transformed into this spaceusing the measured properties of each detection instrument and the reduction procedures. The dataare modified as little as possible and our assumptions about detection efficiency and photometricerror are built directly into the model itself. This process is described in detail in Section 4.3.4.2 DataOur dataset results from ten total HST orbits during Cycle 21 (GO-12971, PI: Richer), one ofwhich was rejected due to loss of guide stars. Observations were taken between November 2012and September 2013. In each orbit, the cluster was simultaneously observed using the Wide FieldCamera 3 (WFC3) and the Advanced Camera for Surveys (ACS). All of the observed fields areshown in Figure 4.1. Each orbit was split into two exposures for each of the four filters. Theseare listed in Table 4.1. Splitting each observation in two allows us to eliminate cosmic ray strikes.The exposure lengths are chosen to fit nicely around the required buffer memory dumps from thetelescope. The observations are designed to take advantage of the maximum possible observationtime during each orbit. The split observation also allows for a dither between each exposure. Adither is a small shift in the pointing of the telescope between observations.The raw WFC3 images were combined by drizzling onto a single reference frame [Fruchteret al., 2009]. The ACS images were reduced to photometry independently. In total 107,000 objectswere detected in the WFC3 fields, and 228,000 objects were detected in the ACS fields. The CMDsresulting from these catalogues are shown in Figure 4.2.4.3 Analysis4.3.1 Artificial Star TestsThe process by which we generate and measure artificial stars is one of many steps that are integralto the translation of the white dwarf cooling model from theory-space into data-space. By “theory-space” we mean the quantities that come out of the theoretical white dwarf calculations, such astheir physical size, their atmospheric composition, their surface temperatures, etc. By “data-space”we mean the dimensions of the space in which the data points exist. The WFC3 data that we areconsidering are the F225W magnitude, the F336W magnitude, and the radial distance in pixelsfrom the cluster centre of every white dwarf in the sample. So the points in the data set all lie inthis three-dimensional space. Similarly, the ACS data-space is made up of F435W, F555W, andradial distance. Any comparison of model to data necessarily requires the model and the data to be48Figure 4.1: The pattern of the observed fields in the core of 47 Tuc. Fields imaged by WFC3are shown in blue, while those observed by ACS are in pink. In the final data set onlynine of the ten orientations for both WFC3 and ACS were used, due to loss of guide starsduring one exposure.49Figure 4.2: The colour magnitude diagrams of the two fields. The central UV data are shownon the left. The visible CMD from the surrounding annulus is shown on the right. Ineach plot there are three distinct sequences, although in the left plot the central sequenceis quite faint and sparse. From left to right these are the white dwarf sequence of 47 Tuc,the main sequence of the Small Magellanic Cloud, and the main sequence of 47 Tuc.transformed to common quantities that are then compared. In our previous paper on white dwarfcooling [Goldsbury et al., 2012] we took an approach that involved modifying the data to produce atime-temperature relation that was then compared to the model. In this chapter we do the opposite,which is to leave the data as untouched as possible and performing all of the transformations on themodel itself.The artificial star tests for the WFC3 drizzled image (see Section 3.3) are described below. Theprocess is virtually identical for the ACS data, but is only performed on one of the nine visits.We make the assumption that the cluster’s projection is azimuthally symmetric, so the error andcompleteness measured across one visit are the same for the eight other ACS visits.First, 1000 values were chosen, evenly spanning the range of 18 to 28 magnitudes in F225W.1000 values were also chosen, evenly spanning 16 to 26 in F336W. A radial coordinate grid was thendefined in the drizzled image frame. This grid consists of concentric rings centred on the core ofthe cluster with a spacing of 50 pixels between each. Points on each annulus are evenly distributedwith a 50 pixel spacing between adjacent artificial stars on the same annulus. This means there50are approximately 25,000 points in the drizzled frame. This spacing ensures that the density ofartificial stars has no effect on the measurements of the artificial stars. Although the artificial starscan not interfere with measurements of other artificial stars, the real stars can interfere with themeasurements of artificial stars, which is what is measured by these tests. The locations of theartificial stars in magnitude and position space are shown in Figure 4.3.Figure 4.3: The real data are shown in black and the artificial star input grid is shown in red.The magnitude plot is centred on the white dwarf sequence.For each of the 1000 defined magnitudes in each of the two filters, one modified image wasgenerated. Each of these images is identical to the original science image, but with the additionof 25,000 artificial objects, all of the same magnitude. So, in total, 2000 new WFC3 images weregenerated (1000 for each filter). This process was repeated to generate an additional 2000 newimages for the ACS data.Each of these 2000 images was run through the same reduction pipeline used to generate thecatalogue of real objects. Whether or not the artificial object was found during the reduction wasrecorded. If the object was found, the output photometry was recorded as well.This output artificial star catalogue is the basis for the error and completeness model. We canindex the input artificial F225W magnitudes with i, the input F336W magnitudes with j and theinput positions with k. The real catalogue only contains stars that are identified in both filters.Because the positions indexed by k are consistent across both filters, we can consider any value iand j and access the output measurements of those two values at a consistent position k. Whilemeasuring the completeness of the artificial stars in this way we assess i∗ j ∗ k = 2:5×1010 uniquestars in the data-space.Considering only stars that are recovered, each unique combination (i; j;k) also has a measured∆F225W and ∆F336W (the difference between the magnitudes measured by the reduction pipelineand those input), as well as a ∆R. From here on, the difference in input and output position is ignored(∆R = 0). Our analysis assumes that stars have no uncertainty in their position in the field. All ofthese data are then used to build a five dimensional array indexed by: F225Winput, F336Winput, R,∆F225W, and ∆F336W.We call this array E, more explicitly: E(F225Winput;F336Winput;R;∆F225W;∆F336W). E is51Figure 4.4: Completeness fraction in magnitude space for a fixed radial distance is shown onthe top. Completeness fraction along the two-dimensional plane of R and F225Winputwith a fixed value of F336Winput is shown on the bottom. Radial distance is in units ofWFC3 pixels, which correspond to 0.04 arc seconds each. These surfaces are the valuesfrom C described in Equation 4.1. C is indexed by three parameters. In each plot above,one of these three is fixed, and C is shown varying across the remaining two.52generated by counting the artificial stars found after passing through our reduction procedure anddividing the count in each bin by the number of input artificial stars. The completeness of objects ata given pair of magnitudes and position in the field can be calculated asC =∫ ¥−¥∫ ¥−¥E d(∆F225W) d(∆F336W): (4.1)How the completeness varies across the three-dimensional data-space isC(F225Winput;F336Winput;R).The completeness values on two orthogonal slices through this three dimensional space are shownin Figure 4.4. The odd structure in the bottom plot is real and not an effect of poor sampling. Ourobservational design leads to much larger total integration time in the core of the cluster than inthe surrounding annuli, and this effect contributes positively to completeness as you move towardthe core of the cluster. However, stellar density also increases dramatically toward the core, whichdecreases the completeness. These two effects together lead to the behaviour seen in the bottomfigure.Another way to slice E is to fix the first three parameters. The variations in E over the remainingtwo parameters (∆F225W and ∆F336W) are the error distribution of those input stars. Two suchdistributions are shown in Figure 4.5.The error distribution describes the expected measured differences between the input and out-put artificial star photometry. From the figures one can see that the error distribution is highlynon-Gaussian, and changes significantly across the space of the data. The five dimensional array Eis the key to transforming from predicted model magnitudes to the likelihood of actually measur-ing objects at various magnitudes in the catalogue. E is an empirically determined model of ourphotometric scatter and measurement biases caused by the crowding of the field, the instrumentand filters used, the observation pattern, and the reduction process. How this information is used inconjunction with our theoretical cooling models is described in Section 4.3.3.4.3.2 MESA ModelsTo calculate the cooling models we will compare to the data, we used Modules for Experiments inStellar Astrophysics or MESA [Paxton et al., 2011]. We used these modules to perform simulationsof stellar evolution starting with pre-main-sequence models with a metallicity of Z= 3:3× 10−3appropriate for the cluster 47 Tuc. Wind parameters in the RGB and AGB phases of evolution areadjusted to produce white dwarfs with varying hydrogen layer thickness. The parameter adjusted ish . The mass loss rate on the RGB is calculated from Reimers [1975]:M˙RGB = 4×10−13h LL⊙RR⊙M⊙M: (4.2)The mass loss rate on the AGB is calculated from Blo¨cker [1995]:M˙AGB = 4:83×10−9(MM⊙)−2:1( LL⊙)2:7M˙RGB: (4.3)Hydrogen layers with qH between 10−6 and 10−3 are produced, where qH is the fraction ofthe WD mass that is in the hydrogen layer. The thickness of the residual helium layer is −1:45 <log10(qHe)<−1:30 for all of the models that reach the white dwarf stage. Since we adjust the windparameters and not the thickness of these layers directly, the masses of hydrogen and helium vary53Figure 4.5: 1, 2, and 3s contours (enclosing 68%, 95%, and 99:7%) of the photometric errordistribution are shown at two points in the data-space. The top panel shows the errordistribution for faint stars close to the outside of the inner WFC3 field. The lower panelshows the error distribution for brighter stars closer to the centre. R is again in units ofpixels.54coincidentally as h is changed. This process is stochastic, and very similar wind parameters canproduce white dwarfs with hydrogen layers that differ by orders of magnitude. However, using thisrevision of MESA, it was not possible to produce helium layers that vary over such a range. Themass of the helium layer stays relatively constant over all of our simulations. As a result, we do notconsider the helium layer mass as a separate parameter in our fitting. The fit distributions that wefind can be considered as the result after marginalizing over this parameter.During the evolution of these models, neutrino production is multiplied by a factor An rangingfrom 0.1 to 3.18. We use SVN revision 5456 of MESA and start with the model 1M pre ms t wdin the test suite. We changed the parameters initial mass and initial z of the star andadjusted the parameter log L lower limit to−6 so the simulation would run well into the WDcooling regime. We defined the starting time of our cooling models to be the last local maximum inluminosity. After the time we define as t0, the models only decrease in brightness.Initially it was thought that the mass of the white dwarf would be a dominant parameter. In fact,using the MESA framework, it is rather difficult to produce white dwarfs that differ significantlyfrom the canonical globular cluster value of 0:53M⊙. In order to change this mass by a largeamount, the main sequence evolution time scale and the turn-off mass need to be adjusted to valuesthat are entirely inconsistent with other observations in 47 Tuc and other clusters. Additionally,even when these models were produced, it was found that varying the mass does not produce a largeeffect on the predicted distribution of white dwarfs on our CMD. For this reason, we have chosen tofix this parameter and instead look to the hydrogen layer thickness and the neutrino production asour fit parameters. Both of these parameters have a much larger effect on the white dwarf coolingdistribution than the mass does.4.3.3 Moving the Model to Data-SpaceThe result of these calculations is a grid of white dwarf cooling sequences over varying qH and An .Each combination of hydrogen layer thickness and neutrino factor produces one cooling sequencethat describes how the physical properties of the white dwarf change with time. The ones weare interested in are effective surface temperature T , surface gravity g, and radius r. These threeparameters can be thought of as functions of time t but also qH and An . The state of the modelresulting directly from calculations in MESA isT (t;qH;An); (4.4a)g(t;qH;An); (4.4b)r(t;qH;An): (4.4c)These parameters exist in theory-space and need to be transformed to quantities in data-space.The first step is to determine how the model objects look through the HST filters. To do this weuse spectral models for DA white dwarfs described in Tremblay et al. [2011b] and the publiclyavailable HST filter throughputs available from STScI [2009]. We also model the reddening usingthe interstellar reddening model from Fitzpatrick [1999b].The spectral model grid provides the spectral flux density (energy per unit time per unit wave-length), F , at the surface of the object as a function of T and g. Explicitly we can write this asF(l ;T;g).55The Fitzpatrick reddening curve is parametrized by E(B−V ), and RV . These are defined asfollows. AV is the extinction in V and AB is the extinction in B (the fraction of power lost in filterV and filter B respectively). E(B−V ) = AB−AV ; and RV = AVAB−AV . The symbols B and V referto magnitudes in the B and V filters (centred at 445nm and 551nm) in the photometric standardsystem of filters Johnson and Morgan [1953]. We fix RV = 3:1 which is found to apply in mostcases [Fitzpatrick, 1999a]. This value also gives good agreement between our white dwarf modelsequence and our data, suggesting that it is appropriate for 47 Tuc as well. We do not vary RV , andit remains fixed at the standard value of 3.1 for the remainder of our analysis. The reddening curveis called k(l ;e), and is the fraction of the power transmitted as a function of wavelength betweenthe object and the observer.The filter throughput curves describe the fractional power transmittance as a function of wave-length for a given filter. This will be labeled as Si(l ) where i is an arbitrary index to indicate whichfilter is being used.All of these objects are then combined as shown in Equation 4.5 to generate a predicted magni-tude:m′i(T;g;r;e;d) =−2:5log10(∫ ¥0 lFSikdl∫ ¥0 lF0Sidl)+5log10(dr): (4.5)Here F0 is the zero-point flux density appropriate for each filter (the flux density that correspondsto a magnitude of zero). The new parameter d is introduced as the distance to the cluster. If wesubstitute Equations 4.4 into 4.5, we getm′i(t;qH;An ;e;d): (4.6)At this point m′ is used to describe the magnitude rather than m, as the photometric biasesdescribed in Section 4.3.1 have not been taken into account yet. The two filters used in the WFC3observations are F225W and F336W. Magnitudes in these filters will be referred to as m2 and m3,respectively.The thing we wish to calculate is the probability density, f , of finding a white dwarf at a givenset of magnitudes (M2 and M3) and a given position in the field (R). This is f = dNdm2dm3dR . In orderto manipulate Equation 4.6 into the form we need to compare to the data, we first solve for t as afunction of everything else, and then take two derivatives:dtdm′2dm′3=dtd(m′2)2dm′2dm′3: (4.7)Next we multiply by the formation rate of objects in the field, N˙ = dNdt :dNdm′2dm′3=dtd(m′2)2dm′2dm′3dNdt: (4.8)Now we have to bring back the error distribution, E, from Section 4.3.1. Before we refered tothis as E(F225Winput;F336Winput;R;∆F225W;∆F336W), but now we can see that F225Winput =m′2,F336Winput = m′3, ∆F225W= m2−m′2, and ∆F336W= m3−m′3. So what we have is:E(m2−m′2;m3−m′3;m′2;m′3;R): (4.9)56We can perform an integral to transform Equation 4.8 into the form that we want, which is f :f =dNdm2dm3dR= r(R)∫ ¥−¥∫ ¥−¥dNdm′2dm′3E(m2−m′2;m3−m′3;m′2;m′3;R)dm′2dm′3: (4.10)This step is a convolution integral of the error distribution with the theoretical two dimensionalluminosity function. The error distribution changes with position in data-space. The assumed radialdensity distribution, r(R), must also be included at this point. Neglecting this term is equivalentto an assumption of a density distribution that is constant with radius. We use a King-Michiemodel [King, 1966, Michie, 1963] to calculate the projected density distribution. The assumed Kingradius and tidal radius values are taken from Goldsbury et al. [2013] and Harris [1996], respectively(r0 = 32′′ and rt = 3800′′). Note that r(R) in Equation 4.10 is not the three dimensional densitydistribution, it is the projected and azimuthally integrated distribution with units of objects per pixel.As discussed in Chapter 2, it is not appropriate to model a population of stars with varying massesusing a single King-Michie model. The white dwarfs we are considering are all within 1% of0:53M⊙, and so their density distribution can be described by a single model. The model used isappropriate for giant stars or turn-off stars in the cluster. Even though the white dwarfs are muchless massive (0:53 compared to 0:8M⊙), their distribution has not relaxed dynamically since theylost this mass.At this point we have a distribution f that gives the probability density of finding an objectat any point in the data-space (F225W;F336W;R) for any combination of our model parameters(qH;An ;d;e; N˙). We will use a test statistic referred to as the “unbinned likelihood” to determinehow well the data are predicted by this model.4.3.4 The Unbinned LikelihoodA quick derivation of the unbinned likelihood starts by considering a bin in the data-space withdimensions ∆R by ∆m2 by ∆m3. This bin is centred at at the location (R j;m2k;m3l). From the modeldiscussed above, the predicted number of objects in this bin is given as f (R j;m2k;m3l)∆R∆m2∆m3.The probability of finding n objects in this bin is given by the Poisson distribution:P(n; f (R j;m2k;m3l)) =[ f (R j;m2k;m3l)∆R∆m2∆m3]ne− f (R j;m2k;m3l)∆R∆m2∆m3n!: (4.11)The limit is then taken as the bin size approaches zero. This leaves a number of bins equal tothe number of data points with one object in them, and an infinite number of bins with zero objectsin them. Then there are two possible forms for the probability to take:P(1; f (R j;m2k;m3l)) = f (R j;m2k;m3l)∆R∆m2∆m3e− f (R j;m2k;m3l)∆R∆m2∆m3 (4.12a)P(0; f (R j;m2k;m3l)) = e− f (R j;m2k;m3l)∆R∆m2∆m3 (4.12b)57The log-likelihood is the sum of the logarithm of the probabilities of all of the bins:logL=åilog( f (Ri;m2i;m3i))+åilog(∆R∆m2∆m3)−åjklf (R j;m2k;m3l)∆R∆m2∆m3: (4.13)The first two sums are only performed over i, which is the index of the data points, since theseterms only exist in the case of Equation 4.12a. The third sum is carried out for all bins, even thosewith no data in them, since this term is common to Equations 4.12a and 4.12b. The second sumdoes not depend in any way on the data or our chosen model. It will be an additive constant to thelog-likelihood and so can be ignored. Given the limit as the bins become infinitesimal, we can alsosee that the third sum becomes an integral of the probability density over the entire data-space. Withthese two modifications we havelogL=åilog( f (Ri;m2i;m3i))−∫∫∫data−spacef dRdm2dm3; (4.14)logL=åilog( f (Ri;m2i;m3i))−Npred: (4.15)The integral in Equation 4.14 is the total number of objects that our model predicts in the regionof data-space in which we are doing the comparison. To simplify the notation, we have written thisas Npred in Equation 4.15. The bound region in data-space does not have to be rectangular or anyother particular shape; as long as the integral in Equation 4.14 is performed over exactly the sameregion from which we select the data, this works.The description up until this point has been specific to the central WFC3 data, but a virtuallyidentical approach is used for the ACS data. All nine ACS fields were reduced separately. Sincewe assume that the cluster distribution and also the incompleteness are azimuthally symmetric, weonly perform artificial star tests on one of these nine images. These tests are still parametrized byradial position R and two magnitudes (in this case F435W and F555W ). We can use these tests tobuild a second probability density function f specific to the ACS data.Although it is left out of the likelihood equations above, the probability density functions f , andtherefore the likelihood as well, both depend on the choice of model parameters (qH;An ;d;e; N˙).The first four parameters are consistent across both the ACS data and the WFC3 data. This isbecause all of the stars should have the same hydrogen layer thickness, have the same neutrinophysics, lie at the same distance, and experience the same reddening. The fifth parameter N˙ is therate at which white dwarfs are being produced in each field, so will not be the same for both datasets.To understand the N˙ parameter more intuitively, imagine first selecting the region shown inFigure 4.6. If we could watch this CMD change over hundreds of millions of years, we would seeobjects moving down the white dwarf sequence as they cool and leaving the region through thebottom. We would also see new white dwarfs arriving through the top of the region as they contractand cool coming from the horizontal branch. N˙ would be the flux into the fit region as observed inthis way. This flux changes on long timescales, as the cluster turn-off moves down the mass functionto lower mass stars, but it will change slowly enough over the age range of white dwarfs in the fitregion that we can assume that it is constant.In total we now have six parameters; (qH;An ;d;e; N˙WFC3; N˙ACS). We can combine the two58likelihood functions to generate a likelihood over parameter space for all of the data together.4.3.5 Degenerate Parameters and Prior KnowledgeA number of the parameter constraints are partially degenerate. Distance and reddening are partiallydegenerate, due to the fact that a larger reddening requires a larger extinction in each filter, producinga similar effect to putting the cluster farther away. The formation rate is degenerate with distance,which is because the cooling rate roughly follows a power law with time. For a perfect power law itis impossible to distinguish between a vertical and horizontal shift in log-space. Because the modelcooling rate is similar, changing the timescale looks roughly the same as making all of the objectsbrighter or fainter.In addition to these degeneracies, the data set is not particularly useful for obtaining constraintson these parameters. Without using any prior distributions on other parameters, the fit distributionsfor distance or reddening are 2-3 times wider than current constraints. This isn’t worrying however,since the parameters we are most interested in studying are qH and An , and the data are able toprovide useful constraints on these parameters. So for these reasons, we choose to use priors for allof our model parameters except for qH and An . This leads to a constraint on the neutrino productionrate, dependent on what we already know about the other cluster parameters.Table 4.2: Values for N˙ are determined as described in the paragraph following this table.They are provided in units of objects per Myr.Parameter Mean Std. Dev. Reference(m−M)0 13.30 0.13 Woodley et al. [2012a]E(B−V ) 0.04 0.02 Harris [1996]N˙WFC3 (Myr−1) 8.2 0.3N˙ACS (Myr−1) 2.61 0.07We use Gaussians for all of the prior distributions. The parameters of these are shown in Table4.2. The values for the formation rate of white dwarfs in each field (N˙) are estimated in the same wayas in Goldsbury et al. [2012]. We use data on the giant branch in each field, along with completenesscorrections from the artificial star tests, and main sequence evolution models from Dotter et al.[2008b] to model the rate at which stars are leaving the main sequence and evolving up the giantbranch. We make the assumption that this is the same rate at which stars are arriving onto the whitedwarf sequence.4.3.6 Restricting the Fitting RegionThe region in magnitude space over which we fit our model is restricted in both fields to the regionimmediately around the white dwarf sequence. These regions are shown in blue in Figure 4.6. Thereare two main reasons for the requirement that the fitting be restricted to these regions. In both theWFC3 and ACS data, the lower white dwarf sequence begins to overlap with the main sequence ofthe Small Magellanic Cloud, which lies in the background of 47 Tuc.We do not include the SMC in our model and so fitting over this region would be unreliable.In the case of the WFC3 data, the white dwarf models also deviate significantly from the measuredsequence at fainter magnitudes. The cause of this discrepancy is unclear. After a thorough analysis,59Figure 4.6: The colour magnitude diagrams of the two fields. Overplotted are the region inwhich we fit our model in blue, and the white dwarf cooling sequence in red for themean prior values given in Table 4.2, before convolution with the error model.we find no apparent systematic offset introduced by the photometric reduction process. The pho-tometry of the sequence is slightly biased by diffuse background light and crowding of the innerfield, but this could account for only 10% of the discrepancy observed at most. Additionally, theerror model described in Section 4.3.1 fully accounts for the discrepancy attributable to diffuse back-ground light. The model offset even persists through entirely different reductions performed withand without drizzling. Avoiding the lower WFC3 sequence does have a considerable impact on thetotal sample size, but if the fitting region were not restricted as shown in Figure 4.6, all of the whitedwarf models would be strongly ruled out. Using the region shown, there are regions of parameterspace that are consistent with the data. Although the plot shows this region in colour-magnitudespace, fitting is done in magnitude-magnitude space.4.3.7 Markov Chain Monte Carlo SamplingTo explore the model parameter space, we use a Markov chain Monte Carlo (MCMC) method.Specifically, we use the Metropolis-Hastings algorithm [Metropolis et al., 1953, Hastings, 1970]with Gibbs sampling [Geman and Geman, 1984]. The total chain used to generate the parameterdistributions in this analysis is 100,000 points long. The process of thinning, in which only every60Figure 4.7: The auto-correlation of each chain parameter as a function of the lag betweenpoints in the chain. No memory of the past position in parameter space is retained beyond40 steps.Nth point in a Markov chain is kept, is often used to reduce the correlation between points withinthe chain. We do not follow this approach, since these small-scale correlations are averaged out bythe comparably long length of the total chain. This is explained in detail in Link and Eaton [2012].We can look at the auto-correlation of each of our model parameters to see that beyond 20 stepsthe auto-correlation drops to effectively zero (Figure 4.7). It is possible to break the total chain intomany smaller chains, each of which gives the same result within statistical uncertainties (Figure4.8). These simple tests indicate that the total chain will be not be biased by the small length scalecorrelations within it.4.4 Results4.4.1 Hydrogen Layer Thickness (qH) and Neutrino Scaling (An )The only two parameters fit without priors are An and qH. The joint fit distribution for these pa-rameters is shown in Figure 4.9. The distribution for qH alone is shown in Figure 4.10. This is thehistogram of the qH values in the total chain. It is also the likelihood distribution of this parameter61Figure 4.8: A comparison of the likelihood distribution for each of our model parametersshown for the entire chain in red. The thinner black lines are independent segments 1000steps long. Each of these separate chains is much longer than the length scale on whichcorrelations persist in the chains. All of them average to the same broad distribution thatwe find with the entire chain, but with more statistical noise.after marginalizing over all of the other parameters. The range of 2:2×10−5 to 8:9×10−5 contains95% of the distribution. The distribution for An alone is shown in Figure 4.11. This is the histogramof An values in the total chain. It is also the likelihood distribution of this parameter after marginal-izing over all other parameters. The range of 0.83 to 1.22 contains 95% of the distribution. Thedefault amount of neutrino production in the MESA models corresponds to An = 1. Our findingsare consistent with this. Keep in mind that these parameters are fit over the range 0.1 to 3.18 for Anand 10−6 to 10−3 for qH. If another peak exists outside this range, our analysis will not find it.Figures 4.12 and 4.13 show how the cooling models change with varying An and qH respectively.From these figures it is easy to see that changing An has a much larger effect than changing qH overthe range of values we are considering. This results in a tighter constraint on An .4.4.2 Neutrino SpeciesThe MESA software uses the results of Itoh et al. [1996] to calculate the neutrino emissivity. Theneutrino production in white dwarfs in this temperature range is dominated by the plasmon neutrino62Figure 4.9: Contours of the joint fit distribution of An and qH. 1, 2, and 3s contours (enclosing68%, 95%, and 99:7%) are shown.process. The term “plasmon” refers to a photon moving through a plasma. Unlike photons in avacuum, plasmons travel more slowly than c due to their interaction with the plasma. Some of theenergy of the photon is transferred to the free electrons in the plasma. This allows the interactingplasmon to decay into a neutrino-antineutrino pair while still conserving momentum: g → n + n¯[Winget et al., 2004, Kantor and Gusakov, 2007]. The factor we refer to as An is a multiplicationof Qplasma from equation 4.1 of Itoh et al. [1996]. The plasmon neutrino rate is calculated in theMESA code as described in Itoh et al. [1996], and this rate depends on the central density of thewhite dwarf, core temperature of the white dwarf, the weak mixing angle (more often parametrizedas sin2qW), and the effective number of non-electron neutrino species (n). We multiply the rate thatcomes out of this calculation by the factor An . The rate Qplasma can be affected by any of the fourparameters listed above, and so the factor An could be thought of as a change to any combinationof these four. However, interpreting An as a change to either the central density or temperaturewould be inconsistent, since these parameters were not consistently altered elsewhere in the code.It therefore only makes sense to interpret a constraint on this parameter as either a constraint onsin2qW and n.The relation between An , sin2qW, and n isAn µ (C2V +nC′2V ); (4.16)63Figure 4.10: A histogram of the posterior distribution of qH values from the Markov chain.This is the likelihood distribution of qH after marginalizing over all other parameters inthe model with prior distributions on each, except for An . The range of 2:0× 10−5 to8:8×10−5 contains 95% of the distribution (shown in red).whereCV and C′V are defined asCV =12+2sin2qW (4.17)andC′V = 1−CV : (4.18)Using these relations, the likelihood in Figure 4.11 can be mapped onto the two-dimensionalspace of sin2qW and n. This is shown in Figure 4.14. From the particle data group [Olive et al.,2014] sin2qW = 0:23155. Looking at this value in Figure 4.14, one can see that changing n willnot produce a large change in An . This means that our constraint on An does not allow us to put areasonable constraint on n.4.4.3 Fit QualityWhile these distributions in parameter space inform us which values lead to better or worse fits,they cannot indicate which fits are actually consistent with the dataset. This is possible in the same64Figure 4.11: A histogram of the posterior distribution of An factor values from the Markovchain. This is the likelihood distribution of An after marginalizing over all other param-eters in the model with prior distributions on each, except qH. The range of 0:80 to 1:21contains 95% of the An fit distribution (shown in red).way that one could determine a best-fitting line for data that are clearly not linear. A maximum inlikelihood will be found regardless of whether or not the fit is good. It is possible that even thoughwe have found constrained likelihood distributions for qH and An , the best fitting model would stillbe a poor fit to the data.To evaluate whether this is the case we use our best-fitting model to create simulated datathrough Monte Carlo sampling. The model being referred to here is Equation 4.10, so these simu-lated samples take full account of the empirically determined incompleteness as well as photometricbias and scatter. We generate millions of simulated data sets from the best-fitting model and in eachcase calculate the unbinned likelihood of the simulated data against the model. These values makeup the distribution of likelihood we would expect if data such as ours were drawn from our best-fitting model. We can then ask whether or not the likelihood value we find from comparing the realdata to the model is an outlier in this distribution. If the real likelihood value falls well outside thedistribution of likelihood values from the simulated data, then we would conclude that it is unlikelythe real data could be drawn from our best-fitting model. The distribution of likelihood values fromsimulated samples, as well as the likelihood value of the real data compared to the model, are shownin Figure 4.15.Integrating from the red dashed line to the left tail of the distribution in Figure 4.15 will enclose65Figure 4.12: The relation between luminosity and time shown for varying An values while qHis held fixed.4.9% of the total area under the black line. This indicates that 4.9% of the time we would expect todraw data from our model that differs at least as much as the real data do.4.4.4 Cooling Models From Other GroupsWe compared the four cooling models used in Goldsbury et al. [2012] to this data set as well. Thesemodels come from Fontaine et al. [2001b], Hansen et al. [2015], Lawlor and MacDonald [2006a],and Renedo et al. [2010]. The model from Hansen et al. [2015] is more recent than the versionfit in Goldsbury et al. [2012], while the other models are the same. These models do not haveany intrinsic parameters to be varied, so the fits are only performed over the parameter space (withpriors) defined by the four parameters in Table 4.2. The quality of fit assessment for each modelis done in the same way as for our MESA models (summarized in Figure 4.15. Simulated data isrepeatedly drawn from each best-fitting model and compared back to that model. The distributionof resulting likelihood values is then integrated from the value obtained when comparing the realdata to the model. The results are shown in Table 4.3. With this data set we cannot strongly rule outany of these models.66Figure 4.13: The relation between luminosity and time shown for varying qH values while Anis held fixed.Table 4.3: The “fit probability” is defined as the probability that data drawn from the givenmodel would produce a likelihood value less than the real data do when compared back tothat model.cooling model fit probabilityFontaine et al. [2001b] 0.032Hansen et al. [2015] 0.025Lawlor and MacDonald [2006a] 0.078Renedo et al. [2010] 0.118this work (MESA) 0.0494.5 ConclusionsWe have described the implementation of a fitting statistic, known as the unbinned maximum like-lihood, which is common in some fields but relatively unused in astronomy. This approach does notrequire that the data be binned in any way, or an assumption of Gaussian error bars, both of whichare necessary for a c2 fit. Given the same data and model, an unbinned approach will allow fortighter constraints on the same parameters than a binned fit will. In the case we have described, theerror model is combined directly with the theoretical model to produce a continuous distribution for67Figure 4.14: The likelihood of sin2qW and n. For the accepted value of sin2qW (0.23155 fromOlive et al. [2014]),C′V is very small. This means that the dependence of An on n is veryweak. As a result, our constraint on the neutrino production (An ) does not reasonablyconstrain n.where points are expected to be found in data-space. The data can be kept unmodified from theirraw state, and all assumptions about parameters and errors are built into the model itself.We used this approach to test white dwarf cooling models against the largest sample of whitedwarfs ever collected in a single globular cluster (47 Tuc). Fitting in six-dimensional parameterspace was performed with an MCMC sampler. We found constraints on the thickness of the Hydro-gen layer and the amount of neutrino production in the white dwarf cores.The constraint on qH is shown in Figure 4.10. The data prefer thicker hydrogen layers (qH =3:2×10−5) and we can strongly rule out thin layers (qH = 10−6). The 95% range of this parameteris 2:0× 10−5 to 8:8× 10−5. Our constraint on neutrino production is much better than for qH,with 0:80 < An < 1:21 at 95% confidence. In Hansen et al. [2015] the best-fitting hydrogen layerthickness was found to be qH = 4× 10−4. Our analysis strongly rules out that value, but we alsofind better fitting hydrogen layers at the thick end of the range usually considered by modellers[Bergeron et al., 2001]. The 95% confidence range on fs in Hansen et al. [2015] (which is equivalentto our parameter An ) was found to be 0:6 < fs < 1:7, which agrees with what we have found. Thedifferences in these values can be attributed both to the inclusion of additional data in our analysis,and to different statistical assumptions made in fitting the models.68Figure 4.15: The distribution of likelihood values expected when sample data are drawn fromthe best fitting cooling model is shown as the black distribution. The red line indicatesthe log-likelihood value found when comparing the real data to the model. We wouldexpect to draw data from our best fitting model that differs at least as much as the realdata do from that model 4.9% of the time.69Chapter 5Conclusion5.1 Cluster Mass SegregationIn this thesis I have presented research that probes two separate phenomena in globular clusters.In Chapter 2, I have shown how the measurable effects of mass segregation in clusters lead to adiscrepancy between the surface brightness profile of a cluster and the mass distribution of thecluster. A bias in measured cluster parameters follows directly from this discrepancy and can beseen in the values currently used in the literature. A quantitative analysis of these effects lead usto a simple power-law model to account for them. The index of this power law quantifies howstrongly the stars are segregated by mass within a particular cluster. This approach was used toquantify mass segregation in 54 clusters using public catalogues [Sarajedini et al., 2007]. Thedimensionless power-law index correlates strongly with other dynamical cluster parameters, suchas central velocity dispersion, cluster relaxation time, and the physical size of the cluster.A table of the determined fit values modeling core radius as a function of mass was published.Using this table of values, it is possible to roughly calculate the expected distribution of stars at agiven mass. The Harris catalogue [Harris, 1996] is the most commonly used resource for globularcluster parameters. Using the core radius value from this catalogue to describe the distribution ofany stars other than the most massive in the cluster would be inappropriate. The fits summarized inChapter 2 allow one to determine a core radius specific to the type of stars being considered, for anymass.Future work in this area would focus on deriving a physically motivated model for quantifyingmass segregation in clusters. The main benefit of the work summarized in Chapter 2 was to bring tothe attention of other researchers the source of this systematic bias, and to provide a straightforwardway to approximately correct it.5.2 White Dwarf CoolingIn Chapter 3, we used data to measure a rate of cooling for white dwarfs in 47 Tuc. We foundthat the cooling observed was consistent with a broken power law. The transition between the twoslopes was found to occur at around 20,000K. When the white dwarfs are young and hot, between35,000K and 20,000K, their effective surface temperature changes with time like t−0:19. The slopein this region agrees with all of the models tested within the uncertainties of our data. Once the white70dwarfs cool past this transition, between 20,000K and 7,000K, their effective surface temperaturechanges with time like t−0:42. This value is quite close to Mestel’s predicted rate of t−0:40, and itagain agrees with all tested models in this region.Despite the agreement with models in both the hot region and the cool region separately, whentested against the entire data set, all of the models are ruled out. This can be attributed to the slightlyhotter and much sharper transition between the two slopes that is seen in the data. Additionally, noneof the models tested are consistent with each other at the level of our data uncertainties. That is tosay, if we were to draw a sample of white dwarfs from one of the theoretical models and add ouruncertainties, this simulated data would be inconsistent with the other theoretical models.White dwarf cooling models are used to determine ages of various stellar populations, such asthe disk of the galaxy [Oswalt et al., 1996], old globular clusters [Hansen et al., 2007a], and youngeropen clusters [von Hippel, 2005]. The findings of Chapter 3 indicate much larger uncertainties inthe physics of the cooling models themselves than is commonly considered by independent groupsfitting their own models to data.In Chapter 4 we probed the physics of white dwarf cooling further, but with entirely different(and much more robust) techniques from those used in Chapter 3. Using a statistical approachknown as the “unbinned maximum likelihood”, we fit models to our data in a way that does notrequire data to be binned, and we do need to assume that our errors are Gaussian. We build ourknowledge of biases and errors in the photometry directly into the model we are fitting, while thedata are left alone, and we predict the probability density of white dwarfs continuously across ourentire data space.This approach was used to constrain various theoretical white dwarf cooling models against thelargest sample of white dwarfs ever observed in a single globular cluster. This analysis stronglysuggests that most of the white dwarf models that are publicly available from various groups areinadequate at predicting the distribution of white dwarfs in the globular cluster 47 Tuc. We measurethe likelihood distribution across our six dimensional parameter space using an MCMC sampler.This analysis provided new constraints on the thickness of the hydrogen layer and the neutrinoproduction rate in white dwarf cores. The direct comparison between our results and those ofHansen et al. [2015] suggest that the unbinned likelihood leads to tighter constraints on the samedata set.The most important aspect of this thesis is not any result in particular, but the method describedin Chapter 4. The unbinned likelihood allows one to take advantage of the full dimensionality ofthe available data in a statistically appropriate way. Chi-square is still the statistic of choice for thisfield of research, but this statistic is no longer appropriate when the number of objects in each binis small. Small counts in each bin is a natural result of higher dimensional data. The unbinnedlikelihood method described here could easily be adapted to a larger number of filters. White dwarfcooling models have not been test in more than two filters simultaneously in a single field. Thiswould be a logical next step, and would add further constraints to the models.71BibliographyL. G. Althaus, A. H. Co´rsico, J. Isern, and E. Garcı´a-Berro. Evolutionary and pulsationalproperties of white dwarf stars. 18:471–566, Oct. 2010. → pages 25J. Anderson, A. Sarajedini, L. R. Bedin, I. R. King, G. Piotto, I. N. Reid, M. Siegel, S. R.Majewski, N. E. Q. Paust, A. Aparicio, A. P. Milone, B. Chaboyer, and A. Rosenberg. The acssurvey of globular clusters. v. generating a comprehensive star catalog for each cluster. TheAstronomical Journal, 135(6):2055–2073, 2008. → pages 12P. Bergeron, S. Leggett, and M. T. Ruiz. Photometric and spectroscopic analysis of cool whitedwarfs with trigonometric parallax measurements. The Astrophysical Journal SupplementSeries, 133(2):413, 2001. → pages 68P. Bergeron, G. Fontaine, P.-E. Tremblay, and P. M. Kowalski. Synthetic colors and evolutionarysequences of hydrogen and helium-atmosphere white dwarfs, Oct. 2011. URLhttp://www.astro.umontreal.ca/∼bergeron/CoolingModels. → pages 43J. Binney and S. Tremaine. Galactic Dynamics. Princeton University Press, 1987. → pages 13T. Blo¨cker. Stellar evolution of low and intermediate-mass stars. i. mass loss on the agb and itsconsequences for stellar evolution. Astronomy and Astrophysics, 297:727, 1995. → pages 53M. Catelan, P. B. Stetson, B. J. Pritzl, H. A. Smith, K. Kinemuchi, A. C. Layden, A. V. Sweigart,and R. M. Rich. Deep HST Photometry of NGC 6388: Age and Horizontal-Branch Luminosity.651:L133–L136, Nov. 2006. doi:10.1086/509720. → pages 19B. Chaboyer and L. M. Krauss. Theoretical uncertainties in the subgiant mass-age relation and theabsolute age of w centauri. The Astrophysical Journal, 567:45–48, 2002. → pages 35F. D’Antona and I. Mazzitelli. Cooling of white dwarfs. Annual Review of Astronomy andAstrophysics, 28:139–181, 1990. → pages 25A. Dotter, B. Chaboyer, D. Jevremovic´, V. Kostov, E. Baron, and J. W. Ferguson. The dartmouthstellar evolution database. The Astrophysical Journal Supplement Series, 178(1):89–101, 2008a.→ pages viii, 12, 16, 19, 34, 35A. Dotter, B. Chaboyer, D. Jevremovic´, V. Kostov, E. Baron, and J. W. Ferguson. The dartmouthstellar evolution database. The Astrophysical Journal Supplement Series, 178(1):89, 2008b. →pages 8, 5972A. Dotter, J. Kaluzny, and I. B. Thompson. Globular cluster ages from main sequence fitting anddetached, eclipsing binaries: The case of 47 tuc. The Ages of Stars, Proceedings of theInternational Astronomical Union, IAU Symposium, 258:171–176, 2009. → pages 35A. Dotter, A. Sarajedini, J. Anderson, A. Aparicio, L. R. Bedin, B. Chaboyer, S. Majewski,A. Marı´n-Franch, A. Milone, N. Paust, G. Piotto, I. N. Reid, A. Rosenberg, and M. Siegel. TheACS Survey of Galactic Globular Clusters. IX. Horizontal Branch Morphology and the SecondParameter Phenomenon. The Astrophysical Journal, 708:698–716, Jan. 2010.doi:10.1088/0004-637X/708/1/698. → pages 19N. Duric. Advanced astrophysics. Cambridge University Press, 2004. → pages 36, 37E. L. Fitzpatrick. Correcting for the effects of interstellar extinction. The Publications of theAstronomical Society of the Pacific, 26:63–75, 1999a. → pages 29, 56E. L. Fitzpatrick. Correcting for the effects of interstellar extinction. The Publications of theAstronomical Society of the Pacific, 26:63–75, 1999b. → pages 29, 55G. Fontaine, P. Brassard, and P. Bergeron. The potential of white dwarf cosmochronology. ThePublications of the Astronomical Society of the Pacific, 113(782), 2001a. → pages ix, 8, 9, 34,39, 40, 44G. Fontaine, P. Brassard, and P. Bergeron. The potential of white dwarf cosmochronology1.Publications of the Astronomical Society of the Pacific, 113(782):409–435, 2001b. → pages 25,66, 67D. A. Forbes and T. Bridges. Accreted versus in situ Milky Way globular clusters. Monthly Noticesof the Royal Astronomical Society, 404:1203–1214, May 2010.doi:10.1111/j.1365-2966.2010.16373.x. → pages 19A. Fruchter, M. Sosey, W. Hack, L. Dressel, A. Kockemoev, J. Mack, M. Mutchler, and N. Pirzkal.The multidrizzle handbook. STScI, Baltimore, 2009. → pages 48D. Geisler, G. Wallerstein, V. V. Smith, and D. I. Casetti-Dinescu. Chemical Abundances andKinematics in Globular Clusters and Local Group Dwarf Galaxies and Their Implications forFormation Theories of the Galactic Halo. The Publications of the Astronomical Society of thePacific, 119:939–961, Sept. 2007. doi:10.1086/521990. → pages 19S. Geman and D. Geman. Stochastic relaxation, gibbs distributions, and the bayesian restoration ofimages. Pattern Analysis and Machine Intelligence, IEEE Transactions on, (6):721–741, 1984.→ pages 60R. Goldsbury, H. B. Richer, J. Anderson, A. Dotter, A. Sarajedini, and K. Woodley. The acs surveyof globular clusters. x. new determinations of centers for 65 clusters. The Astronomical Journal,140(6):1830–1837, 2010. → pages 18R. Goldsbury, J. Heyl, H. Richer, P. Bergeron, A. Dotter, J. Kalirai, J. MacDonald, R. Rich,P. Stetson, P.-E. Tremblay, et al. An empirical measure of the rate of white dwarf cooling in 47tucanae. The Astrophysical Journal, 760(1):78, 2012. → pages iii, 47, 50, 59, 6673R. Goldsbury, J. Heyl, and H. Richer. Quantifying mass segregation and new core radii for 54milky way globular clusters. The Astrophysical Journal, 778(1):57, 2013. → pages iii, vii, 20,57R. G. Gratton, E. Carretta, and A. Bragaglia. Multiple populations in globular clusters. TheAstronomy and Astrophysics Review, 20(1):1–49, 2012. → pages 4B. Hansen. The astrophyics of cool white dwarfs. Physics Reports, 399:1–70, 2004. → pages 25,44B. Hansen, J. Anderson, J. Brewer, A. Dotter, G. Fahlman, j. Hurley, j. Kalirai, I. King, D. Reitzel,H. Richer, R. M. Rich, M. Shara, and P. Stetson. The white dwarf cooling sequence of ngc 6397.The Astrophysical Journal, 671:380–401, 2007a. → pages ix, 25, 39, 40, 45, 71B. Hansen, J. Kalirai, J. Anderson, A. Dotter, H. Richer, R. Rich, M. Shara, G. Fahlman, J. Hurley,I. King, et al. An age difference of two billion years between a metal-rich and a metal-poorglobular cluster. Nature, 500(7460):51–53, 2013. → pages 47B. M. Hansen and E. S. Phinney. Stellar forensicsi. cooling curves. Monthly Notices of the RoyalAstronomical Society, 294(4):557–568, 1998. → pages 8B. M. Hansen, J. Anderson, J. Brewer, A. Dotter, G. G. Fahlman, J. Hurley, J. Kalirai, I. King,D. Reitzel, H. B. Richer, et al. The white dwarf cooling sequence of ngc 6397. TheAstrophysical Journal, 671(1):380, 2007b. → pages 8, 47B. M. Hansen, H. Richer, J. Kalirai, R. Goldsbury, S. Frewen, and J. Heyl. Constraining neutrinocooling using the hot white dwarf luminosity function in the globular cluster 47 tucanae. TheAstrophysical Journal, 809(2):141, 2015. → pages 47, 66, 67, 68, 71W. E. Harris. Spatial structure of the globular cluster system and the distance to the galactic center.The Astronomical Journal, 81:1095–1116, 1976. → pages 1W. E. Harris. Globular cluster systems in galaxies beyond the local group. Annual Review ofAstronomy and Astrophysics, 29:543–579, 1991. → pages 1W. E. Harris. A catalog of parameters for globular clusters in the milky way. The AstronomicalJournal, 112:1487, 1996. → pages vii, 1, 11, 13, 15, 16, 18, 19, 21, 22, 32, 57, 59, 70W. K. Hastings. Monte carlo sampling methods using markov chains and their applications.Biometrika, 57(1):97–109, 1970. → pages 60J. Heyl, H. B. Richer, E. Antolini, R. Goldsbury, J. Kalirai, J. Parada, and P.-E. Tremblay. Ameasurement of diffusion in 47 tucanae. The Astrophysical Journal, 804(1):53, 2015. → pages4, 47N. Itoh, H. Hayashi, A. Nishikawa, and Y. Kohyama. Neutrino energy loss in stellar interiors. vii.pair, photo-, plasma, bremsstrahlung, and recombination neutrino processes. The AstrophysicalJournal Supplement Series, 102:411, 1996. → pages 7, 62, 6374H. Johnson and W. Morgan. Fundamental stellar photometry for standards of spectral type on therevised system of the yerkes spectral atlas. The Astrophysical Journal, 117:313, 1953. → pages29, 56J. S. Kalirai, D. Saul Davis, H. B. Richer, P. Bergeron, M. Catelan, B. M. S. Hansen, and R. M.Rich. The masses of population ii white dwarfs. The Astrophysical Journal, 705:408–425, 2009.→ pages 26J. S. Kalirai, H. B. Richer, J. Anderson, A. Dotter, F. Greg, B. Hansen, J. Hurley, I. King,D. Reitzel, R. M. Rich, M. Shara, P. Stetson, and K. Woodley. A deep, wide-field, andpanchromatic view of 47 tuc and the smc with hst: observations and data analysis methods. TheAstronomical Journal, 143(1), 2012. → pages viii, 26, 29, 30, 31, 39E. Kantor and M. Gusakov. The neutrino emission due to plasmon decay and neutrino luminosityof white dwarfs. Monthly Notices of the Royal Astronomical Society, 381(4):1702–1710, 2007.→ pages 63I. R. King. The structure of star clusters. III. Some simple dynamical models. The AstronomicalJournal, 71:64, Feb. 1966. doi:10.1086/109857. → pages 4, 13, 14, 57M. Koleva, P. Prugniel, P. Ocvirk, D. Le Borgne, and C. Soubiran. Spectroscopic ages andmetallicities of stellar populations: validation of full spectrum fitting. Monthly Notices of theRoyal Astronomical Society, 385:1998–2010, Apr. 2008.doi:10.1111/j.1365-2966.2008.12908.x. → pages 19A. N. Kolmogorov. Sulla determinazione empirica di una legge di distribuzione. na, 1933. →pages 40L. M. Krauss and B. Chaboyer. Age estimates of globular clusters in the milky way: Constraintson cosmology. Science, 299(5603):65–69, 2003. → pages 1W. Kutta. Beitrag zur na¨herungweisen integration totaler differentialgleichungen. 1901. → pages14D. Q. Lamb and H. M. Van Horn. Evolution of crystallizing pure 12C white dwarfs. TheAstrophysical Journal, 200:306–323, 1975. → pages 7R. R. Lane, L. L. Kiss, G. F. Lewis, R. A. Ibata, A. Siebert, T. R. Bedding, and P. Sze´kely. TestingNewtonian gravity with AAOmega: mass-to-light profiles and metallicity calibrations from 47Tuc and M55. Monthly Notices of the Royal Astronomical Society, 401:2521–2530, Feb. 2010.→ pages 32T. Lawlor and J. MacDonald. The mass of helium in white dwarf stars and the formation andevolution of hydrogen-deficient post-agb stars. Monthly Notices of the Royal AstronomicalSociety, 371(1):263–282, 2006a. → pages 8, 66, 67T. M. Lawlor and J. MacDonald. The mass of helium in white dwarf stars and the formation andevolution of hydrogen-deficient post-agb stars. Monthly Notices of the Royal AstronomicalSociety, 371:263–282, 2006b. → pages ix, 39, 4075W. A. Link and M. J. Eaton. On thinning of chains in mcmc. Methods in Ecology and Evolution, 3(1):112–115, 2012. → pages 61D. E. McLaughlin and R. P. van der Marel. Resolved Massive Star Clusters in the Milky Way andIts Satellites: Brightness Profiles and a Catalog of Fundamental Parameters. The AstrophysicalJournal Supplement Series, 161:304–360, Dec. 2005. doi:10.1086/497429. → pages 11A. McWilliam and R. A. Bernstein. Globular cluster abundances from high-resolutionintegrated-light spectra. i. 47 tuc. The Astrophysical Journal, 684:326–347, 2008. → pages 35L. Mestel. On the theory of white dwarf stars. i. the energy sources of white dwarfs. MonthlyNotices of the Royal Astronomical Society, 112:583, 1952. → pages 5, 25N. Metropolis, A. W. Rosenbluth, M. N. Rosenbluth, A. H. Teller, and E. Teller. Equation of statecalculations by fast computing machines. The Journal of Chemical Physics, 21(6):1087–1092,1953. → pages 60R. W. Michie. The dynamics of spherical stellar systems, IV. Monthly Notices of the RoyalAstronomical Society, 126:499, 1963. → pages 4, 13, 57A. Milone, G. Piotto, L. Bedin, I. King, J. Anderson, A. Marino, A. Bellini, R. Gratton, A. Renzini,P. Stetson, et al. Multiple stellar populations in 47 tucanae. The Astrophysical Journal, 744(1):58, 2012. → pages 4P. Miocchi, B. Lanzoni, F. R. Ferraro, E. Dalessandro, E. Vesperini, M. Pasquato, G. Beccari,C. Pallanca, and N. Sanna. Star count density profiles and structural parameters of 26 Galacticglobular clusters. ArXiv e-prints, July 2013. → pages 22S. Moehler, D. Koester, M. Zoccali, F. R. Ferraro, U. Heber, R. Napiwotzki, and A. Renzini.Spectral types and masses of white dwarfs in globular clusters. Astronomy and Astrophysics,420:515–525, 2004. → pages 26K. A. Olive et al. Review of Particle Physics. Chin. Phys., C38:090001, 2014.doi:10.1088/1674-1137/38/9/090001. → pages x, 64, 68T. D. Oswalt, J. A. Smith, M. A. Wood, and P. Hintzen. A lower limit of 9.5 gyr on the gae of thegalactic disk from the oldest white dwarf stars. Nature, 382:692–694, 1996. → pages 25, 45, 71B. Paxton, L. Bildsten, A. Dotter, F. Herwig, P. Lesaffre, and F. Timmes. Modules for experimentsin stellar astrophysics (mesa). The Astrophysical Journal Supplement Series, 192(1):3, 2011. →pages 8, 53W. H. Press, S. A. Teukolsky, W. T. Vetterling, and B. P. Flannery. Numerical Recipes: The Art ofScientific Computing. Cambridge University Press, 2007. → pages 16D. Reimers. Circumstellar absorption lines and mass loss from red giants. Memoires of the SocieteRoyale des Sciences de Liege, 8:369–382, 1975. → pages 53I. Renedo, L. G. Althaus, M. M. Bertolami, A. D. Romero, A. H. Co´rsico, R. D. Rohrmann, andE. Garcia-Berro. New cooling sequences for old white dwarfs. The Astrophysical Journal, 717(1):183, 2010. → pages 8, 66, 6776I. Renedo, L. G. Althaus, M. M. Miller Bertolami, A. D. Romero, A. H. Co´rsico, R. D. Rohrmann,and E. Garcı´a-Berro. New Cooling Sequences for Old White Dwarfs. The AstrophysicalJournal, 717:183–195, July 2010. → pages ix, 39, 40A. Renzini, A. Bragaglia, F. R. Ferraro, R. Gilmozzi, S. Ortolani, J. Holberg, J. Liebert,F. Wesemael, and R. C. Bohlin. The white dwarf distance to the globular cluster ngc 6752 (andits age) with the hubble space telescope. The Astrophysical Journal Letters, 465(1):L23, 1996.→ pages 47M. Salaris, I. Dominguez, E. Garcia-Berro, M. Hernanz, J. Isern, and R. Mochkovitch. TheCooling of CO White Dwarfs: Influence of the Internal Chemical Distribution. TheAstrophysical Journal, 486:413, Sept. 1997. → pages 39E. E. Salpeter. The luminosity function and stellar evolution. The Astrophysical Journal, 121:161,1955. → pages 36A. Sarajedini. An acs survey of galactic globular clusters, Jan. 2011. URLhttp://www.astro.ufl.edu/∼ata/public hstgc/. → pages 12A. Sarajedini, L. R. Bedin, B. Chaboyer, A. Dotter, M. Siegel, J. Anderson, A. Aparicio, I. King,S. Majewski, A. Marn-Franch, G. Piotto, I. N. Reid, and A. Rosenberg. The acs survey ofgalactic globular clusters. i. overviews and clusters without previous hubble space telescopephotometry. The Astronomical Journal, 133(4):1658–1672, 2007. → pages 12, 23, 70M. Schwarzschild. Structure and evolution of the stars. 1958. → pages 7P. B. Stetson. Daophot - a computer program for crowded-field stellar photometry. ThePublications of the Astronomical Society of the Pacific, 99:191–222, 1987. → pages 27P. B. Stetson. The center of the core-cusp gobular cluster m15: Cfht and hst observations, allframereductions. The Publications of the Astronomical Society of the Pacific, 106:250–280, 1994. →pages 27P. B. Stetson. Homogeneous photometry for star clusters and resolved galaxies. ii. photometricstandard stars. The Publications of the Astronomical Society of the Pacific, 112:925–931, 2000.→ pages viii, 34, 35STScI. Multi-cycle treasury programs, Nov. 2009. URLhttp://www.stsci.edu/institute/org/spd/mctp.html/. → pages 29, 55S. C. Trager, I. R. King, and S. Djorgovski. Catalogue of Galactic globular-clustersurface-brightness profiles. The Astronomical Journal, 109:218–241, Jan. 1995.doi:10.1086/117268. → pages 11P. E. Tremblay, P. Bergeron, and A. Gianninas. An improved spectroscopic analysis of da whitedwarfs from the sloan digital sky survey data release 4. The Astrophysical Journal, 730:128,2011a. → pages 28, 29P. E. Tremblay, P. Bergeron, and A. Gianninas. An improved spectroscopic analysis of da whitedwarfs from the sloan digital sky survey data release 4. The Astrophysical Journal, 730:128,2011b. → pages 5577H. Van Horn. Cooling of white dwarfs. In White Dwarfs, pages 97–115. Springer, 1971. → pages6, 7T. von Hippel. From Young and Hot to Old and Cold: Comparing White Dwarf Cooling Theory toMain-Sequence Stellar Evolution in Open Clusters. The Astrophysical Journal, 622:565–571,Mar. 2005. → pages 45, 71D. Winget, D. Sullivan, T. Metcalfe, S. Kawaler, and M. Montgomery. A strong test of electroweaktheory using pulsating db white dwarf stars as plasmon neutrino detectors. The AstrophysicalJournal Letters, 602(2):L109, 2004. → pages 63D. E. Winget, C. Hansen, J. Liebert, H. Van Horn, G. Fontaine, R. Nather, S. Kepler, and D. Lamb.An independent method for determining the age of the universe. The Astrophysical Journal,315:L77–L81, 1987. → pages 7D. E. Winget, S. O. Kepler, F. Campos, M. H. Montgomery, L. Girardi, P. Bergeron, andK. Williams. The physics of crystallization from globular cluster white dwarf stars in ngc 6397.The Astrophysical Journal Letters, 693(1):L6, 2009. → pages 47M. Wood. Constraints on the age and evolution of the galaxy from the white dwarf luminosityfunction. The Astrophysical Journal, 386:539–561, 1992. → pages 47K. A. Woodley, R. Goldsbury, J. Kalirai, H. Richer, P.-E. Tremblay, J. Anderson, P. Bergeron,A. Dotter, L. Esteves, G. Fahlman, et al. The spectral energy distributions of white dwarfs in 47tucanae: The distance to the cluster. The Astronomical Journal, 143(2):50, 2012a. → pages 47,59K. A. Woodley, R. Goldsbury, J. S. Kalirai, H. B. Richer, P. E. Tremblay, J. Anderson, P. Bergeron,A. Dotter, L. Esteves, G. G. Fahlman, B. M. S. Hansen, J. Heyl, J. Hurley, R. M. Rich, M. M.Shara, and P. B. Stetson. The spectral energy distributions of white dwarfs in 47 tucanae: Thedistance to the cluster. The Astronomical Journal, 143, 2012b. → pages 29, 35, 40, 44M. Zoccali, A. Renzini, S. Ortolani, A. Bragaglia, R. Bohlin, E. Carretta, F. Ferraro, R. Gilmozzi,J. Holberg, G. Marconi, et al. The white dwarf distance to the globular cluster 47 tucanae and itsage. The Astrophysical Journal, 553(2):733, 2001. → pages 4778

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.24.1-0228563/manifest

Comment

Related Items