The visual lightcurve of comet C/1995 O1 (Hale-Bopp) from 1995-1999

The long-term brightness evolution of the great comet C/1995 O1 (Hale-Bopp) presented a remarkable opportunity to study the behavior of its coma over four years. We used approximately 2200 total visual magnitudes published in the International Comet Quarterly taken from 17 observers during the period of 1995 July - 1999 September to create a secular lightcurve. In order to account for observer differences, we present a novel algorithm to reduce scatter and increase precision in a lightcurve compiled from many sources. It is implemented in a publicly available code, ICQSPLITTER. This code addresses the differences among observers by using a self-consistent statistical approach, leading to a sharper lightcurve, and improving the precision of the measured slopes. To first order, the comet's lightcurve approximately follows a r$^{-4}$ response for both pre- and post-perihelion distances. Interestingly, the pre-perihelion data are better fit with a fifth-order polynomial with inflection points at 4.0, 2.6, 2.1 and 1.1 au. We analyze these specific regions and find that they are associated with physical phenomena in the comet's evolution. Contrary to other reports, the lightcurve shows no evidence for the comet having been in outburst at discovery. Afrho values derived from the visual lightcurve data are consistent with a r$^{-1.5}$ dependence on heliocentric distance, which is similar in shape to those derived from spectroscopy and narrow-band photometry. We present correlation equations for visual magnitudes and CO and H2O production rates, which are consistent with the pre-perihelion visual magnitudes increasing almost entirely due to CO outgassing until a heliocentric distance of about 2.6 - 3.0 au. We also present two correlation equations that should prove highly useful for observation planning and data analysis, and can be generalized to be applicable to other comets.


INTRODUCTION
Comet nuclei are among the best-preserved icy-rocky remnants from the solar system's formation. Far from the Sun, when comets are inactive, the heliocentric dependence of their bare nuclei's brightness will follow an r −2 response due to the reflection of sunlight off the surface, assuming a constant phase angle, θ, surface area and albedo ( Fig. 1). Figure 1. A comet's geocentric distance, ∆, and phase angle, θ, can both affect its apparent brightness.
When comets get close enough to the Sun, the ices incorporated in the nucleus sublimate and generate comae of gas and dust. Once they become active, the heliocentric dependence of their brightness changes due to a variety of parameters about which we know very little, including dust particle size distribution, relative contributions of gas and dust, and possibly other sources of energy (Fulle 2004). Coma brightness may change when the mass loss rate changes, or large particles or grains in the coma fragment into smaller particles without changing the actual mass loss rate. Moreover, the nucleus's spin pole orientation will affect a comet's heating from the Sun (Whipple 1980).
Over human history, most reports of comets are from visual observations of their comae, including some reports of comet 1P/Halley ∼ 800 hundred years ago, e.g., Choi et al. (2018). Estimates of a comet's visual magnitude are often used to assess its level of activity. For comets with high amounts of dust, lightcurves can also provide useful information about the bulk behavior of the dust coma, which when coupled with more rigorously derived values of the dust, and images from infrared, will provide observational constraints for nucleus composition models. Changes in slope may pinpoint a significant transition in the comet nucleus's processing, such as the initiation of distant activity, ramp up of water-ice sublimation, breakup event or outbursts .
One of the most continuous datasets on any comet was the visual secular lightcurve of Halley's coma during its 1986 perihelion compiled from thousands of apparent magnitudes and corrected for geocentric distance of the comet (Green & Morris 1987). The lightcurve showed many significant changes in slope, peaked at ∼ 60 days past perihelion, and continued brighter post-perihelion. A similar opportunity arose with the "great" comet C/1995 O1 (Hale-Bopp), because it was exceptionally bright upon its discovery at ∼ 7 au and easily observed for ∼ 18 months with the naked eye, and much longer with binoculars and small telescopes. From 1995From -1999, dozens of observers dedicated nearly every clear night to recording the visual magnitude of Hale-Bopp and submitted brightness estimates to the International Comet Quarterly (ICQ) 1 to make them available for further analysis by others.
From the ∼ 15,000 magnitudes of Hale-Bopp archived with the ICQ, we studied a subset of 2,240 magnitudes from 17 highly experienced observers. We then created a final lightcurve using a reduced number of 785 magnitudes from 13 observers. In this paper we describe our procedures for removing points, correcting for geocentric distance and phase angle, and our new self-consistent statistical approach for combining different observers, and its efficacy. This reduction process increases the accuracy and precision of a comet lightcurve.
We present this algorithm through the analysis of the comet Hale-Bopp, which was the subject of the largest observational campaign in history. We also highlight some of the major features and implications of the refined lightcurve. We present our open source software package ICQSplitter (Curtis et al. 2020), which was used to develop these results and is discussed in Appendix 1. This serves as the companion paper to our archive submission with the NASA Planetary Data Systems Small Bodies Node archive (Womack et al. 2018).
2. OBSERVATIONS AND REDUCTIONS 2.1. Pruning data points from the original dataset The original apparent visual magnitudes, m app , used in this study were from the ICQ issues 96-100 and provided by Daniel W. E. Green (personal communication), who selected observers with good coverage for both hemispheres (see Table 1). Initially, all 2240 of the visual magnitudes were used. We removed 751 points following the best practice guidelines from Green & Morris (1987) and Green 2010 (personal communication), which included removing values when any of the following occurred: • No magnitude method was listed • Methods of obtaining magnitudes that were not the Vsekhsvyatskij-Stevenson-Sidgwick (S), Van Biesbroeck-Bobrovnikoff-Meisel (B), Modified Out-Out, Infocus, or Extrafocal-Extinction methods

• Poor observing conditions
• Under 20 degrees elevation with no extinction correction applied • Telescopes used for m app = 5.5 and brighter • Binoculars used for m app = 1.5 and brighter • CCD or photoelectric detectors used • Poor quality star comparison catalogs used • More than one magnitude per day submitted by an observer (if this occurred, we kept the measurement made with the smallest aperture instrument, which is believed to be most accurate) Some observers, such as Herman Mikuz, submitted a mixture of visual and CCD magnitude estimates. When used properly, CCD and photoelectric detectors can provide very accurate values. However, they do not always agree well with visual magnitudes determined by eye, and it is not straightforward to combine data obtained with the two techniques (see Kidger (2002); Sosa & Fernández (2009);de Almeida et al. (2009)). Thus, we removed all CCD measurements from this dataset. Since CCDs were not widely used during the passage of Hale-Bopp, this affected relatively few points.
The remaining apparent magnitudes are the starting point for our analysis and are shown in Fig. 2a, which plots magnitude vs. heliocentric distance. The other panels b, c and d of these figures show the lightcurve after correcting for geocentric distance, ∆, phase angle, θ, and observer bias, respectively, as we describe in the next sections. To facilitate comparison with other published studies the data are also shown in slightly different format in Fig. 3a, which plots magnitudes as a function of "days from perihelion."

Orbital Parameters that Affect the Apparent Brightness of a Comet
The changing orbital configuration of Hale-Bopp, Earth and the Sun sometimes caused the comet's brightness to change over time, which if uncorrected, would be mistakenly attributed to an intrinsic change in output. A comet's apparent visual luminosity, L app , is often written in terms of its absolute luminosity at 1 au from the Sun and Earth, L 0 , scaled with a power-law based on its distance to the Earth (∆ in au) and Sun (r in au), where k describes the geocentric distance dependence, which is usually assumed to be "2" due to reflected sunlight, and n, which is called the "activity index," and is meant to quantify the comet's behavior with respect to heliocentric distance. For a relatively new comet, n = 4 is often used, but the value may change during a comet's orbit and can vary between 1 -8 (Meisel & Morris 1976;. In this paper we assume k = 2, and fit models to the lightcurve data to obtain n.
Using the luminosity-magnitude equivalency equation: and assuming k = 2, we can rewrite Eq. 1 in terms of magnitudes m app = m 0 + 5 log 10 ∆ + 2.5n log 10 r, where m 0 is referred to as the absolute magnitude, or how bright a comet would appear if at 1 au from the Earth and Sun for the case where the comet's phase angle, θ, is ignored (Fig. 1). The visible light from the dust coma includes both reflection and scattering of sunlight off of dust particles in the coma, as well as a much smaller contribution from the comet nucleus. Light scattering off dust particles is a complicated process that depends on the physical properties of the particles in the coma, and is wavelengthdependent (Kolokolova et al. 2004). When observers do correct visual magnitudes for scattering, they often assume a linear response, which is appropriate for comets with small phase angles (less than ∼ 30 degrees) (Meech & Jewitt 1987;Kolokolova et al. 2004). A further complication of correcting for phase angle scattering of visual observations is that many comae observed at visual wavelengths have significant contributions from C 2 , CO + and CN gaseous emission, and this light cannot be separated from contributions from the dust without using spectroscopy or narrow-band filters. Hale-Bopp was a very dust-dominated comet, and thus we assumed that all of the light reported by the observers was from dust reflection and scattering. We used the phase correction fuction, φ[θ], normalized to 0 • , from (Schleicher & Bair 2011) 2 , which is a composite of a previous fit for comet Halley photometry (Schleicher et al. 1998), for smaller angles, and the phase function derived in Marcus (2007a,b) for mid-and largephase angles, which is appropriate for the large heliocentric range of the Hale-Bopp lightcurve. Including the phase angle correction term leads to m app = m 0 + 5 log 10 ∆ + 2.5n log 10 r + 2.5 log 10 (φ[θ]). (4) Magnitudes corrected to a geocentric distance, ∆, of 1 au are referred to as m helio and are calculated with m helio = m app − 5 log 10 ∆ using the geocentric distance of the comet at the time the observation is made. Magnitudes are further corrected to a phase angle, θ, of 0 degrees, which we call m phase , using We plot the lightcurve for each step in individual panels of Figs. 2 and 3.

Statistical analysis of differences among the observers
In analyzing visual data from multiple observers, the questions inevitably arise of which data to reject, and under what justification, and whether combining data from observers, each with his or her own systematic errors, leads to a biased result. Without instrumental calibration, there is no certain answer to these questions, but as discussed by (Kidger 2002;Mousis et al. 2014), such calibration is itself problematic, and in any case is not available for the observations discussed here.
We offer a systematic approach to combining data from multiple observers yielding a self-consistent consensus fit. In application to comet Hale-Bopp, the procedure does not significantly affect the grossest measure, the activity index n, but, applied to data already corrected for geocentric distance and for phase, does reduce the statistical error bars.
We assume three categories of errors: 1. Every observer reports the brightness of an object on a scale that is shifted up or down from other observers, but by the same number of magnitudes, δ obs ,   . Fifth-order polynomial fits to the Hale-Bopp lightcurve data for each step of data reduction. Correcting for geocentric distance (red) and phase angle (green) makes noticeable changes to the lightcurve shape and generally increases the comet's brightness.
independent of distance or brightness. Without instrumental calibration, we can best estimate δ obs as that observer's mean deviation from a consensus fit to the data (that is, an average).
2. Some observers may have a slope bias, underestimating the brightness of dimmer objects and overestimating those of brighter ones, or vice-versa. While it is difficult to correct for such error without calibration, the bias can be detected (relative to the consensus fit), and that observer's data discarded. Table 1. List of observers whose ICQ data were analyzed. We also include the final shift, −δ obs , applied to their magnitudes, based on the statistical analysis discussed in section 2.3. Cells marked "n/a" signify that the observer too few observations for statistical analysis or that the observations were determined to bias the slope of the lightcurve (against log distance) relative to the consensus and so were discarded in the final iteration. 3. Finally, some observers may have a great deal of scatter in their data but no bias. We can weight these observations less in the fits.
We seek a consensus fit to the data. The comet exhibits small but noticeable deviations from a power law (luminosity ∼ r −n ) on time (distance) scales larger than any outbursts. In particular, the pre-perihelion data between 2 and 1 au suggest a positive curvature in the graph of (minus) magnitude against log distance (See Fig. 5.) A straight-line fit would penalize observations that report this feature accurately. Based on the total number of apparent features between the closest (0.91 au) and furthest (9.4 au post-perihelion) observation, we fit to fifth-order polynomials. The effect is to reflect the smoothed consensus.
The deviations between an observer's measurements and the consensus fit at the same distance are considered noise. The set of all such deviations by one observer defines the noise distribution for that observer, characterized by a mean δ obs , variance σ 2 obs , skewness, excess kurtosis, etc.
To compute the self-consistent fit, we iterate the following until convergence to a fractional tolerance of 0.0001 (absolute tolerance if any fitting parameter is less than 0.0001) of all six polynomial fitting parameters: i. We fit a fifth-order polynomial through all the data of the (possibly shifted) magnitudes against log distance by the method of least squares (Legendre 1805), weighting each observation inversely as the observer's variance of deviations, σ 2 obs . Initially, the magnitudes have been corrected for geocentric distance and for phase angle, but have not been shifted. Also initially, the weights are all equal, because we do not know the distributions of deviations.
ii. Using the same data as in i, we can also fit a straight line, recovering the activity index n and statistical error bars on that fit based on the shifted, weighted data.
iii. For each observer, we consider the distribution of deviations between the observations and the fit at the given distance (i.e., in the graph, the vertical vectors between data points and the polynomial fit). For each observer separately, we treat the distribution as noise and estimate the mean, δ obs , and variance, σ 2 obs . iv. For each observer separately, we shift all magnitudes by −δ obs as computed in iii.
For the present data set, these iterations converge to the specified tolerance after between ten and twelve iterations. We now have a self-consistent fit and set of shifted data. Row 4a of Table 2 shows small changes in best-fit slopes after the procedure, with the statistical error of the straight-line fit cut roughly in half. Note that nonstraight-line features in the consensus polynomial fit limit how far the statistical error bars in a straight-line fit can shrink. Table 2. Best-fit activity indices n for pre-perihelion and post-perihelion observations before and after corrections and the self-consistent method. (The slope of magnitude against log 10 (r) is equal to 2.5n; see Eq. (3).) Steps 4a and 4b use self consistently shifted magnitudes and weight each observer's data inversely as σ 2 obs , the variance of the noise. Error bars are standard errors of fit (Press et al. 2007) and would be considerably smaller, preperihelion, if slow variations ( Fig. 5) were subtracted. The fitting method reduces statistical error but does not significantly alter the activity index estimates.

pre-perihelion post-perihelion Correction
Step -7.0 to -0.9 au 0.9 to 8.0 au 1. raw data ( The self-consistent procedure has eliminated the need to discard data arbitrarily (e.g., points more than some number of standard deviations above or below a consensus fit, which otherwise would throw off least-squares fits) by weighting points inversely as the observer's variance. However, as noted above, an observer whose systematic error changes with brightness would still affect the activity index adversely.
We can detect such a systematic non-stationarity by applying Student's t-test to each observer's data set, comparing the mean deviation from the consensus fit in the first half of the observer's data (sorted by distance) to the mean deviation in the second half. Since the variances of the two halves may not be equal, we normalize by the "pooled variance" to get an approximate t-statistic (Press et al. 2007). Assuming approximately Gaussian noise, we calculate the p (probability) value that the t-statistic would be as large as observed or larger under the null hypothesis that first and second halves of the data were drawn from the same distribution, i.e., that the observer did not contribute bias to the slope relative to the consensus (stationarity).  For a p value less than 0.05, we reject the null hypothesis and conclude, on the basis of the t-test, that the observer's data bias the slope relative to the consensus. We then discard such observers from the data set and repeat the self-consistent iteration (i-iv), starting from the original (phase-corrected, unshifted) magnitudes but without the discarded observations. If repeating the t-test at the end of second set of iterations results in more observers biasing the slope (relative to the consensus) based on the t-test, we then repeat the procedure until no new observers are discarded. For the pre-perihelion data, one additional observer failed the t-test after the second set of iterations, so we iterated a third time without that observer's data. At the end of the third set of iterations, no new observers failed the t-test.
As shown by Line 4b of Table 2, throwing out observers results in small changes in the slopes. For the pre-perihelion data, seven observers out of 17 were discarded, accounting for 430 of the original 1,003 data points. For post-perihelion data, four out of 12 observers were discarded, accounting for 232 out of the original 486 data points.
For comet Hale-Bopp, the self-consistent procedure resulted in adjustments to the slope roughly comparable to the original statistical error bars, while cutting the error bars approximately in half, as reflected in visibly smaller scatter in the lightcurve. (See Table 2 between lines 3 and 4b and Fig. 2 between panels (c) and (d).) These changes are far smaller than those associated with correcting for geocentric distance and phase. The observations in the present work were taken by people using generally reliable star comparison catalogs. Future work may rely on less homogeneous amateur networks, in which case a self-consistent method for combining and weighting magnitudes, and discarding subsets with systematic slope bias, could prove useful.
The shifts, −δ obs , to the data resulting from this self-consistent procedure are listed in Table 1 and are equivalent to the last step in the equation m app = m 0 + 5 log 10 ∆ + 2.5n log 10 r + 2.5 log 10 (φ[θ]) + δ obs , and we define or as it is written in complete form: Using Eq. 9, we can rewrite Eq. 7 as: which we use to measure n. The m shift values represent the final corrected magnitudes of the lightcurve and are shown in Figs. 2d and 3d, as well as a close-up version in Fig. 5.
Correcting for geocentric distance and phase angle makes noticeable changes to the lightcurve shape and increases the comet's brightness. By design, making selfconsistent shifts to the data does not change the lightcurve's shape, but reduces the uncertainty. For ease of comparison, fifth-order polynomial fits to the data at each correction step are shown in Fig. 4.
The original apparent magnitudes and those at all steps outlined in this paper underwent full panel review with the NASA Planetary Data Systems Small Bodies Node and are publicly archived (Womack et al. 2018). All data reduction was performed using the ICQSplitter program, which is available for download on GitHub. 3 This paper describes the functionality of ICQSplitter Version 3.0 as of 28 January 2020. Since ICQSplitter is still evolving, refer to the on-line documentation for up-to-date information. Appendix A discusses ICQSplitter's functionality and its subroutines in more detail.

Comparison with other Hale-Bopp lightcurves
Several other groups produced visual lightcurves of Hale-Bopp that spanned many years; however, they were not always derived the same way. For example, some did not account for the effect of geocentric distance, even though this had up to a 42% impact on the activity index, n, for Hale-Bopp. Fig. 2 shows an illustration of this effect, where a "shoulder" feature, seen in panel 2a at ∼ 3.5 -4.0 au pre-perihelion, is effectively removed once we account for the geocentric distance in panel 2b. This type of artifact might invite unwarranted physical interpretation. Table 3 lists some of these longterm lightcurve data sets for Hale-Bopp. Our lightcurve is the only one in Table 3 that corrected for light scattering in the coma (see Sec. 2.2). Our calculations show that scattering can be significant for Hale-Bopp's lightcurve and at times contributed up to a 92% change to the calculated brightness (e.g., compare Figs. 2b and 2c and Fig. 4). This result underscores our motivation in sharing the ICQSplitter code, so that this potentially important correction step may be more accessible to other researchers.
The comet magnitude websites from the Comet OBServational Database 4 (COBS) and Seiichi Yoshida 5 provide data and interactive tools to analyze visual and CCD magnitudes submitted to the International Comet Quarterly (ICQ) and the Minor Planet Center (MPC). Both Yoshida and COBS present only apparent magnitudes (m app ), but include correction for geocentric distance. However, when fitting slopes to their lightcurves, they do not correct for phase angle variations, or include any adjustments for observer differences (Yoshida, S. 2020, personal communication, Zakrajsek & Mikuz (2018); Filonenko & Churyumov (2001)).
The COBS online tool provides greater flexibility by allowing the user to choose which starting and ending dates to use in a fit, whereas Yoshida provides fits to preselected temporal ranges. Based on our analysis here, we recommend that caution should be used when interpreting any changes in slopes derived from Yoshida and COBS lightcurves. They present uncorrected apparent magnitudes, which might  (2001) 1995-2000 Corrected for ∆, fit determined by eye Kidger et al. (1997Kidger et al. ( ) 1995Kidger et al. ( -1997 Corrected for ∆, least squares fit to all data Ferrín (2010Ferrín ( ) 1993Ferrín ( -2009 Corrected for ∆, θ (nucleus), least squares fit to brightest magnitudes from many observers include observational artifacts, as discussed in the above mentioned example for Hale-Bopp.
Liller (2001) and Kidger et al. (1997) both applied geocentric distance corrections, but did not correct for phase angle. Liller (2001) constructed Hale-Bopp lightcurves using visual magnitudes from the IAU circulars and derived an activity index through "fitting by eye". There was no assessment of potential observer differences or fitting bias. Kidger et al. (1997) used visual observations from The Astronomer Group and The British Astronomical Association Comet Section. They applied a least-squares fit to their data, but did not consider any adjustments for observer differences.
Another well known lightcurve was published by Ferrín (2010), which includes geocentric distance correction and adjustment of visual magnitudes due to scattering. However, scattering was considered for the nucleus only. This constitutes a very small change for Hale-Bopp's nucleus, since it was surrounded by a significant coma.
Since all of the lightcurves in Table 3 were corrected for geocentric distance (which we refer to as Step 2), we briefly compare their activity index, n, for this Step of the lightcurve calculations, using m helio = m 0 + 2.5n log 10 r, which we derive from Eqs. 3 and 5 6 . Overall, the Step 2 lightcurves for the groups in Table 3 agree within the uncertainties. The online tools provided by Yoshida and COBS, and Liller (2001) lightcurves show the best agreement with our Step 2 values. This is not surprising considering all four teams started with ICQ measurements. Yoshida and COBS give activity indices that are within ∼ 10% of each other's values, Liller's "fit by eye" value, and our measurements. Our lightcurve is also within 2-6% agreement with Kidger et al. (1997), except for 2.6-2.1 au (1996 Oct 25 -Dec 4), when their calues are ∼ 30% higher than ours, and during 1.1-0.9 au (1997 Mar 1 -Apr 1) when their n is ∼ 60% of our index.
6 Two conventions are used to describe lightcurves and, unfortunately, both use the letter n. COBS, Ferrín (2010), Liller (2001), and this paper calculate activity indices, n, according to Eq. 3, whereas Yoshida and Kidger et al. (1997) use n for the combined term 2.5n in Eq. 11. When we mention results from Yoshida and Kidger et al. (1997), we have already divided their n values by 2.5.
However, the data have high disperion in both of these regions in Kidger et al. (1997) and they overlap with our values within their uncertainties. Ferrín (2010) present two possible activity indices within 6.3 au of the Sun that are ∼ 20-30% flatter than everyone else's value for this long range. Their lower index may be due to the fact that they only fit to the 5% brightest points and combine visual and CCD-derived magnitudes. Ferrín (2010) incorporate naked-eye visual magnitudes from the ICQ as well as CCD measurements, so many individual observers contributed results. As we have shown in this paper, some observers submitting to the ICQ are not selfconsistent, and some submit systematically brighter, or lower, magnitudes than the average. Thus, it may be that preferentially fitting the brightest points when multiple observers are involved contributed to biases in Ferrín (2010)'s results for Hale-Bopp.

Pre-perihelion Changes in the Lightcurve Slope
Given the significant effect that the phase angle correction can have on coma brightness, and the benefits of using a self-consistent method for adding data from different observers, the rest of our analysis is completed using our final lightcurve data (Step 4).
Comet Hale-Bopp's lightcurve provides information about the long-term response of the nucleus to solar heating and continuous thermal processing at different levels. Additional information about the nucleus, such as shape, rotation, or active regions is needed to constrain models in detail (Marshall et al. 2019). However, long-term lightcurves with higher resolution and accuracy are useful in constraining the parameter space for nucleus thermal models, and provide an important way of connecting models of coma activity and nucleus processing, especially over multiple orbits ot larger arcs of a given orbit. A change in the lightcurve's shape may indicate the onset of sublimation or outgassing of a new volatile, or pinpoint the timing of prominent outbursts.
To first order, the activity index of the entire pre-perihelion lightcurve is n = 3.63 ± 0.01 (Fig. 5). This is within ∼10% of the r −4 value that is often assumed for relatively new, inbound comets . It is also consistent with nucleus models that are under the influence of a CO-dominated outgassing regime with a possible water-ice phase change from amorphous to crystalline state (Prialnik 1997). As Fig. 5 shows, however, the lightcurve noticeably deviates from a straight line in several places, and a fifth order polynomial is a better fit with inflection points at r = 4.0, 2.6, 2.1 and 1.1 au. The locations of these changes in slope were also noted by Kidger et al. (1997). In Table 4 we list the power-law activity index for the five regions at each of the data reduction steps. Next, we discuss the lightcurve changes in the context of other coma phenomena during the same time.
3.2.1. Region A: 7.2 ≥ r ≥ 4.0 au (1995 July 17 -1996 June 23) Table 4. Hale-Bopp's pre-perihelion lightcurve showed changes in slope at r = 4.0, 2.6, 2.1, and 1.1 au. This table lists the activity index (n) and absolute magnitude (m 0 ) fits for each region. Corrections for geocentric distance and phase angle can produce significant changes in slopes. The preferred activity index for each heliocentric range are listed in bold font.

Heliocentric Dust coma description
Step n m 0 region (au) Several observers reported 1-2 magnitude outbursts soon after discovery at ∼ 7 au (see discussion in Prialnik 1997) and near 4.0 au (Liller 2001), but we see no evidence for either in our lightcurve. The activity index for region A is n = 4.46 ± 0.11. The comet was still far enough from the Sun that water-ice sublimation, the dominant activity mechanism for most comets, was inefficient. The coma was probably controlled by outgassing of CO, CO 2 , and other volatiles (Jewitt et al. 1996;Womack et al. 1997;Crovisier 1997;Biver et al. 2002a), and possibly moderated by crystallization of amorphous water ice (Prialnik 1997).
During this time, optical images show that Hale-Bopp's dust coma exhibited a single radial jet feature which persisted for several months and then transitioned to four or more jets at ∼ 4 au (Braunstein et al. 1997;Mueller et al. 1997). The appearance of the jets coincides with the next change in lightcurve slope and could indicate the onset of several active areas on the nucleus (Sekanina 1997), perhaps related to triggering of subsurface volatile pockets, or mass wasting events due to mechanically or thermally unstable regions on the nucleus.
3.2.2. Region B: 4.0 > r ≥ 2.6 au (1996 Jun 23 -October 25) Here, the slope flattens with n = 3.06, which it maintained for about four months. Previous analyses by Kidger et al. (1997) reported a much flatter r −1.5 brightening law for this region, which was attributed to the onset of water-ice sublimation near 3 au (Bockelée-Morvan & Rickman 1997). However, our calculations show that the flattening did not have a physical cause, but instead was due to a prominent and uncorrected phase angle effect. This is clearly seen if we compare the indices for the geocentric-corrected (n = 1.62) and phase-corrected (n = 3.11) lightcurves of region B in Table 4.
An outburst was reported in mid-September (3.1 au), when at least a few of the jets appeared disrupted (Braunstein et al. 1997;Tao et al. 2000;McCarthy et al. 2007). This outburst did not appear in our lightcurve, nor did it cause a change in the lightcurve, as the coma maintained a slope of n ∼ 3 for another month until ∼ 2.6 au (mid-late October 1996). This slope change is within the heliocentric distance where it is is commonly thought that water-ice sublimation becomes prominent. If rapid water-ice sublimation is initiated it can easily surpass outgassing of all other volatile species. Due to increased thermal flux and spatial distribution throughout the nucleus, water-ice quickly becomes the primary carrier of dust grains Prialnik (1997), which can significantly increase the reflectance behavior of the coma.
3.2.3. Region C. 2.6 >r ≥ 2.1 au (1996 October 25 -December 4) In this region, the lightcurve flattened even more with a r −2.36 trend. During this time, the comet maintained at least 4 nearly-radial jets in the inner coma, possibly from the same active areas that started at the end of region A. If indeed supported by subsurface pockets or newly excavated patches, this sustained jet activity can inform us about the volume of volatile species (mostly water ice) and how accessible it is to sublimation and distinct outward flow.
3.2.4. Region D. 2.1 > r ≥ 1.1 au (1996December 5 -1997 At ∼ 2 au, the lightcurve returns to a r −4 brightening law. No outbursts are detected. During this time the inner coma morphology transitioned once again. Between 1996 Nov 16 (2.33 au) and 1997 Feb 13 (1.23 au) the coma's dominant features changed from the nearly radial jets to curved jets consistent with radiation pressure and beginning of a "gull wing" feature (Braunstein et al. 1997). Modeling and simulations of many optical and infrared images indicated that the changes from jet-like to arc-like features were due to a changing observing perspective and increased angular resolution, and not due to physical changes in the nucleus or coma . This again can be considered consistent with a view of water-ice rich activated pockets that are present just below the surface.
3.2.5. Region E. 1.1 > r ≥ 0.9 au (1997 Mar 1 -Apr 1) Here we see the fastest rate of brightening of any region, with a power-law index of r −7.5 . The extreme brightening here compared to other regions has been interpreted to mean that an outburst occurred during the dates leading up to perihelion . The brightest points in the lightcurve occurred 3 days before perihelion, in agreement with previous analysis (Ferrín 2010;Kidger et al. 1997). During this time the optical images continued to show the inner coma features dominated by curved jets with shells, or arcs.

Post-perihelion lightcurve
Instead of a series of changing slopes, the post-perihelion data are well-fit by a r −3.94±0.02 power law superimposed with some outbursts. We did not notice the breaks claimed by Kidger et al. (1997), nor do we see any evidence that the comet faded more rapidly as water-ice sublimation turned off around 3 au, which they reported. The lack of measurable slope change over 8 au is consistent with the comet nucleus achieving a level of stability after its perihelion passage. This is slightly (∼8%) steeper than the overall change measured over the pre-perihelion lightcurve from 7 to 0.9 au.
Outbursts may have occurred at 2.15 ± 0.07 au and 7.36 ± 0.20 au (Fig. 5), both lasting a few weeks. The 2.15 au outburst in 1997 Jul led to a brightness change of 0.47 ± 0.08 magnitudes. The 7.36 au outburst in 1999 Jan brightened the coma by 0.68 ± 0.23 magnitudes, and was also observed by Liller (2001)'s CCD dataset, who measured an increase of ∆m ∼ 1.1 magnitude.

Absolute magnitude
Absolute magnitudes are a theoretical construct, representing what brightness the comet would have had if viewed at 1 au from both the Sun and Earth. Astronomers typically correct the apparent magnitude to 1 au from the Earth with Eq. 5, but the correction to 1 au from the Sun is often model-dependent. However, since Hale-Bopp passed and was observed twice through the heliocentric distance of 1.0 au, we have the opportunity to measure the absolute magnitude directly. This is valuable for constraining models of comet activity. Using all the m shift values that were within 0.05 au of 1.0 au, we calculate these average absolute magnitudes: Pre-perihelion: m 0 = -1.74 ± 0.13 Post-perihelion: m 0 = -1.92 ± 0.23 where the uncertainties are one standard deviation. Within these uncertainties we can see that the comet was equally bright pre-and post-perihelion at 1 au. These absolute magnitudes confirm that Hale-Bopp was one of the brightest comets ever recorded.

Af ρ Calculations
Af ρ is a quantity often used as a proxy for dust production rate, and it can be thought of as equivalent to the size of a hypothetical disk needed to approximate the light reflected from the dust at the comet's distance. It is derived from A, the albedo of the dust, f , the filling factor, and ρ, the coma radius (A'Hearn et al. 1984) where is used to convert magnitudes, m, to flux, F , and we assume m Sun = −26.7. We used a value of half the observers' reported coma diameter for the aperture radius, ρ.
To our knowledge, these are the first reports of Af ρ calculated from visual magnitudes for Hale-Bopp, or any comet. The quantity is generally not considered to be well-suited for visual magnitudes, since an aperture radius is usually not provided by the observers. Also, as discussed earlier, visual magnitudes recorded without a filter (which apply to all the data in our study) may be contaminated with gaseous emission, especially within ∼ 2 au of the Sun, when gas output increases significantly. Nonetheless, since we could not find evidence of this being done before, we calculated Af ρ values from the final shifted visual magnitudes (with the geocentric correction reversed, since Eq. 12 accounts for geocentric distance).
In Fig. 6 we plot our magnitude-determined Af ρ values along with the values derived from spectroscopic and filter-imaging techniques, using data from Fig. 4 from Weiler et al. (2003). As Fig. 6 shows, the values show a very similar trend in the preand post-perihelion data. When combined, our pre-perihelion and post-perihelion values follow a Af ρ ∼ 2x10 6 r −1.5 cm heliocentric dependence. No heliocentric dependence was provided by Weiler et al. (2003), but both Af ρ curves show broadly similar shapes with a steep increase within ∼ 3 au, and are much flatter beyond that. This behavior change at ∼ 3 au could be indicative of the water-ice sublimation turn-on at this distance. A rapid increase in sublimation-driven activity is due to both ejection of dust particles from the surface and dragging of dust particles from the subsurface (through porous flow or disruptive mass wasting). Typically, our Af ρ values are 2-3 times higher than those derived with spectroscopic or photometric instrumentation. This offset is most likely due to the relatively wide bandpass that visual observations use, which admit gaseous emission into the observer's field of view, such as from C 2 , NH 2 and CO + , especially closer to perihelion.
As a reminder, the magnitudes used for this work were all obtained from visual observations by eye, sometimes with binoculars or telescopes, but always without filters or detectors. These Af ρ values are likely to be contaminated by gas emission, which will lead to an overestimate of the comet's dust coma brightness, and thus, higher Af ρ values. Sky brightness, or less than ideal observing conditions would also lead to an underestimation of the coma diameter in a visual observation, which might contribute to the higher Af ρ resulting from visual observations. Also, CCD magnitudes tend to be fainter than visual magnitudes, because they are obtained with Professional pre-perihelion Professional post-perihelion Amateurs pre-perihelion Amateurs post-perihelion Figure 6. Af ρ values (in cm) derived from our final Hale-Bopp lightcurve data (blue circles) and Af ρ values derived from optical spectroscopy and R-filter images (orange squares) from Fig. 4 of Weiler et al. (2003). Open symbols are used for pre-perihelion data and filled symbols denote post-perihelion. While there is some overlap, in general, the Af ρ values derived from visual magnitudes are frequently 2-4 times higher, on average, than those derived using optical spectroscopy and R-filter imaging. Both sets show that the Af ρ values have a strong inverse relationship with heliocentric distance within 3 au, but appears to change very little beyond that.
short exposures to provide a good centroid for astrometry and avoid saturating the detector. A short exposure often means that the part of the coma will not be present in the photometric aperture thus leading to a lower magnitude. CCD observers also use telescopes with long focal length (vs. a visual observer's need for a short focal length) to achieve a low power and thus make the coma as small as possible to account for the low surface brightness that characterizes the outer part of comae. In contrast, the human eye is more likely to take in the whole coma + nucleus and hence the visual estimates tend to be brighter, leading to higher Af ρ values. With future observations, this can be addressed by using narrow-band filter-imaging or spectroscopy where the dust contribution can be more accurately measured and calibrated (e.g., Milani et al. (2013)).

Correlation of Water and CO production rates with visual magnitudes
We analyzed published CO and H 2 O production rates (Biver et al. 2002b), Q(CO) and Q(H 2 O), respectively, in the context of the fifth-order polynomial fit to our m shift magnitudes (Fig. 4). We wanted to test whether there was evidence for a possible effect of gas production on the visual magnitudes. In addition, this serves as an independent proof of concept that our approach can be readily combined with other data set to promote advanced analysis.
As Fig. 7 shows, the pre-perihelion water correlation slope has a noticeable break at m shift ∼ 2. Hale-Bopp was this bright predominantly when it was at r ∼ 2.6 au (see Figs. 2d,4 and 5).
We derive these correlation equations for the two pre-perihelion H 2 O regimes: m shift ≥ 2 (beyond r ∼ 2.6 au pre-perihelion): m shift ≤ 2 (within r ∼ 2.6 au pre-perihelion): A similar change in the water-magnitude correlation was reported by Bockelée-Morvan & Rickman (1997); however, they saw the break at ∼ 3 magnitudes, corresponding to r ∼ 3.1 au. We attribute the difference to the fact that they used m helio values for their analysis. This point is illustrated in Fig. 7 where we also show the correlation lines computed with our m helio values with dotted lines, which reproduces Bockelée-Morvan & Rickman (1997)'s results. This change in correlation relationship for water production rate and magnitude may be useful in constraining models of Hale-Bopp activity; however, to get the most accurate identification of the heliocentric distance when this change occurred, we recommend applying a phase angle correction to the magnitudes.
With the post-perihelion data, we find a single equation for H 2 O: which is similar to correlation relationships derived for other long-period comets (Sosa & Fernández 2011;Jorda et al. 2008).
In contrast, we find no break in the CO pre-perihelion correlation plot and derive a single equation: and a nearly identical relationship for post-perihelion: The pre-and post-perihelion correlation equations for CO (Eq. 17,18) are remarkably similar, so that we recommend one equation for CO at all heliocentric ranges in comet Hale-Bopp: The correlation data are consistent with the visual magnitudes increasing almost entirely due to CO outgassing until about 2.6 au, when the water also significantly contributes to the dust release. This coincides with the approximate heliocentric distance where water-ice sublimation began to dominate Hale-Bopp's activity.

Estimating CO and H 2 O Production Rates for Other Comets
We can now take the correlation equations we derived for Hale-Bopp and apply it to other comets. If we note our assumption that the variance of the overall response of different comets is smaller than our uncertainties, we can proceed with caution. We may now use visual magnitudes of other new comets to predict what could be their associated gas production rates. Following this approach may be useful for planning observations of H 2 O, CO and even other volatile species (OH, CN, HCN, CH 4 , H 2 CO, and CH 3 OH), if we assume a range of typical cometary ratios for these species relative to water or CO.
Hence, we present another correlation equation for CO and m app , which we derived from our Q(CO) and m helio correlation analysis for both post-perihelion and preperihelion Hale-Bopp data, and is plotted in red dotted line in Fig. 7: log 10 (QCO) = 29.9 − 0.24m helio = 29.9 − 0.24(m app − 5 log 10 ∆) Although not corrected for phase angle 7 , it should still be useful for a quick approximation of the CO production rate using using visual magnitudes, m app from the ICQ, COBS, or Yoshida, and the comet's geocentric distance, ∆. Of course, most comets will not be as productive as Hale-Bopp, but it is still often used as a comparison standard, especially for newly discovered or brighter than average comets. As a demonstration, we consider a comet that has generated a great deal of interest recently: C/2017 K2 (PANSTARRS). Upon discovery it had a coma at very large heliocentric distances (> 13 au), and its activity was hypothesized to be driven by CO outgassing (Meech et al. 2017;Jewitt et al. 2017). This Oort Cloud comet continues to brighten and is heading for a 2022 December perihelion. If we want to estimate its CO production rate for 2020 Jul 1, we use Eq. 20 and the reported apparent magnitude from that date, m app = 16.4 (Minor Planet Electronic Circular 2020-M114) and obtaining ∆ = 8.6 au (JPL Horizons), and predict that Q(CO) = 1.2x10 27 mol s −1 . A 3-σ detection limit for this value would probably require more time than is typically granted at infrared and mm/submm-wavelength telescopes.
If C/2017 K2 (PANSTARRS)'s CO production rate matches that of Hale-Bopp, then we predict that it will be detectable at infrared and mm/submm wavelengths with Q(CO) ∼ 10 28 mol s −1 (Senay & Jewitt 1994;Gunnarsson et al. 2003;Paganini et al. 2013;Wierzchos & Womack 2020), when it reaches m helio ∼ 7.9. Based on recent lightcurve data from COBS or Yoshida, this is likely to occur during the first half of 2021 as the comet reaches m app < 12.5 and ∆ = 8.5 au. It may be detected even sooner with more observing time. Of course, CO emission may still not be detectable in this comet, but such non-detections would also be valuable for constraining models of nucleus composition and activity. Thus, Eq. 20 may be useful for developing time estimates applicable in observing proposals for distantly active comets.
Similarly, we present this equation from our correlation analysis of the postperihelion water production rate and visual magnitudes for planning purposes: log 10 (QH 2 O) = 30.5 − 0.46(m app − 5 log 10 ∆).
This equation is plotted as a dotted blue line (post-perihelion) in Fig. 7. 7 We do not correct for phase angle here, because accurate information about phase correction is unlikely to be known soon after a comet is discovered

SUMMARY AND CONCLUSIONS
We present a detailed analysis of visual magnitudes for comet Hale-Bopp submitted to the International Comet Quarterly from more than a dozen observers, resulting in a lightcurve for each of these steps: 1) apparent magnitudes, 2) magnitudes corrected for geocentric distance only, 3) magnitudes corrected for both geocentric distance and phase angle, and 4) magnitudes corrected for geocentric distance, phase angle, and systematic observer variation.
Correcting Hale-Bopp's apparent magnitudes for geocentric distance brightens the lightcurve and changes its shape, especially around 3 -4 au pre-perihelion. Correcting for light scattering, as measured by the comet's phase angle, changes the lightcurve slope by an average of ∼ 30%. Not correcting for phase angle led to a reported flattening in the lightcurve from 4 to 2.6 au pre-perihelion, which was reported by other groups and attributed to the onset of water-ice sublimation near 3 au. However, we show that most of this reported flattening in this region was an artifact of not correcting for phase angle.
The self-consistent procedure resulted in adjustments to the slope roughly comparable to the original statistical error bars, while cutting the error bars approximately in half. The computer code used to derive the lightcurves, ICQSPLITTER, is available to the community for download on GitHub, for use on other comets.
To first order, the pre-perihelion and post-perihelion lightcurves follow a ∼ r −4 trend, which is consistent with typical trends often assumed for a relatively new comet. A closer look at the pre-perihelion data shows that the lightcurve deviates from a linear fit at these locations: 4.0, 2.6, 2.1 and 1.1 au. The 4.0 au inflection point coincided with the reported turn-on of four or more jets, or active areas, on the comet nucleus. The second change in slope at 2.6 au may coincide with the onset of vigorous sublimation of water ice from the nucleus. We are not aware of any other significant developments in Hale-Bopp's coma for the other changes in slope.
The post-perihelion lightcurve was more linear than the pre-perihelion data, and was remarkably constant out to ∼ 9 au.
The comet appeared to undergo an outburst a few days before perihelion, and two outbursts at r = 2.15 and 7.36 au, post-perihelion.
The final absolute magnitudes derived from the lightcurve are m o = -1.74 and -1.92, which confirms that Hale-Bopp was one of the brightest comets in recent history.
The Af ρ values derived from visual magnitudes in the final lightcurve show a similar power-law dependence on heliocentric distance, r −1.5 , as those derived using spectroscopic and imaging techniques, but the visually-derived Af ρ values are 2-4 times larger.
We calculated and analyzed correlation laws between the published gas production rates (for H 2 O and CO) and final shifted magnitudes. The data are consistent with the visual magnitudes increasing almost entirely due to CO outgassing until ∼ 2.6 au, when the water begins to significantly contribute to the dust release. This coincides with the approximate heliocentric distance where water-ice sublimation began to dominate Hale-Bopp's activity.
We present two equations to use for predicting water and CO production rates in comets using their apparent magnitude reports. Namely, log 10 (QH 2 O) = 30.5 − 0.46(m app − 5 log 10 ∆) log 10 (QCO) = 29.9 − 0.24(m app − 5 log 10 ∆) These relations have been validated for Hale-Bopp, which is a comet with prolonged and sustained activity that has been observed over a large heliocentric range. Further validation and improvement on the accuracy of such magnitude-production rate correlations will have to wait for future opportunities with similarly active objects. We suggest these relations be used as a quick reference for guiding modeling work and proposal planning. extract the heliocentric distance, geocentric distance, and phase angle of the target. This function will perform a single query of JPL HORIZONS over the range of dates provided in increments inputted by the user. The default time interval is 30 minutes increments. For instance, if your first date is 1996:01:19 00:00, final date is 1996:01:19 01:00, and your increment size is every 30 minutes then it will query JPL Horizons for the ephemerides of your object at 1996:01:19 00:00, 1996:01:19 00:30, and 1996:01:19 01:00. As JPL only allows users to pull 100,000 in a single query, if one wishes to analyze all of their data in a single compilation of ICQSplitter then small time increments should be avoided. ICQSplitter uses the ephermides data to perform heliocentric corrections following equations discussed in section 2.2.

-PHASE
Applies phase corrections to the given magnitudes. If -heliocentric and -phase are called at the same time then the phase angle corrections will be applied onto the heliocentric corrected magnitudes as discussed in section 2.2, else JPL will be queried for the first time and phase corrections will be applied to the raw magnitudes. The phase angles are cross referenced to Dave Schleicher's Composite Dust Function for Comets (Schleicher & Bair 2011).

-STATS
Performs the statistical analysis discussed in section 2.3. The program will automatically split any dataset into pre-and post-perihelion and perform the statistics on each set separately. ICQSplitter follows procedures for regression analysis through the methods of singular value decomposition using NumPy's Linear Algebra package. After a polynomial fit has been taken to convergence, Python's Statistics package is used to perform the Students t and probability tests on each observer's data. If one observer is found to fail the stationarity test in either epoch, then that observer is removed from the dataset and the procedure is repeated. The -stats command is always issued after -heliocentric and -phase (if those commands have also been given).

-PLOT
This command will produce plots with Pythons matplotlib package. For instance, if a user runs her code on data with the -phase and -stats commands then -plots will produce individual graphs of m app , m helio , m phase , and m shift .