Comets formed billions of years ago in the outer parts of the solar system and are among the oldest and purest samples of the makeup of the solar nebula. By studying the gasses emitted by comets as they approach the sun allows us to study the chemical makeup of the early solar system. It also give insight into the chemistry of the early earth, since many comets impacted the surface of the earth during the great bombardment. Comets which I have observed are Hale-Bopp C/1995 O1, LINEAR C/2002 T7, and NEAT C/2001 Q4. See my publications section for details.

Outflow Modeling

Models used by optical, infrared, and single-dish radio telescopes to translate their measured quantities into real world values (column density, etc.) do not always apply to observations done by interferometers. Single-element telescopes are most sensitive to large-scale structures, which for comets translates to the outer coma where both temperature and outflow velocity are well behaved and are best described by the Haser model (Haser, L. 1957, Acad. R. Belgique Bull. Cl. Sci. Ser. 5, 43, 740). However, arrays are sensitive to the inner coma where both the temperature and outflow velocity can vary significantly. I developed the Variable Temperature and Outflow Velocity model (VTOV) to do for arrays what the Haser model does for single dish observations (see Friedel et al. 2005, ApJ, 630, 623. The VTOV model is written as a Miriad fortran routine and is described in the following text (the example refers to calculations for Comet LINEAR (C/2002 T7), but the numerical values can be adjusted to suit any cometary size and distance):

The variable temperature and outflow velocity (VTOV) model consists of concentric spherical shells centered on a cometary nucleus of 5 km. The shells have a thickness of 100 km (with the exception of the first shell, which was 95 km thick) out to a cometocentric radius of 10,000 km and a thickness of 2000 km for radii from 10,000 to 30,000 km. Thirty thousand kilometers is equivalent to a cometocentric radius of 94" on the sky for comet LINEAR. This distance was chosen as the outer limit of our model because structures of this size contribute less than 0.1% to the peak flux, and the population of HCN at this radius is less than 40% of that of the central shell due to photodissociation. Thus, by the combination of these factors the contribution from any shells beyond 30,000 km would be negligible.

While quantities such as the outflow velocity and temperature vary in a very nonlinear way across the coma of the comet, the changes across any given shell are small enough to be considered linear. Thus, the quantities calculated for each shell, unless otherwise specified, are calculated at the average radius of each shell, since they are very close to the average value for the entire shell.

For the model, calculate the following quantities for each shell i.

  1. The time, ti, each molecule spends in the shell is given by
    where Ri and Ri-1 are the outer and inner radii of shell i, respectively, and Vi is the average velocity of the molecule across the shell (from Combi et al. 1999, ApJ, 512, 961).
  2. The total time, ti,tot, it takes a molecule to travel from the nucleus to midway between Ri-1 and Ri is given by
  3. The number of molecules, ni, in each shell is given by
    where Q is the production rate, α is the photodissociation rate [assumed to be α(1 AU) = 1.3×10-5 s-1], and rH is the heliocentric radius of the comet. Note that at this point it does not matter what value we give Q, since it does not vary from shell to shell and thus is a constant for this part of the model.
  4. The BIMA array's flux recovery factor for the projected size of the shell on the sky is defined as the peak flux of a shell compared to that of a point source. This is a combination of a geometric effect (i.e., the flux is spread out over larger areas for larger shells), which is dominant for source sizes comparable to our synthesized beam size, and the fact that the array resolves out flux from sources that are much larger than the synthesized beam (array‐filtering effect). The recovery factor is found in the following way.
    1. Each source size (or shell) is modeled by first creating a Gaussian1 source image for each concentric shell, of the appropriate FWHM and at the position on the sky of the comet, with the MIRIAD task imgen. Every Gaussian shell is assumed to contribute equally to the flux recovery for a uniformly extended source.
    2. In order to remove geometric and filtering effects, each image is run through uvmodel to produce a u-v data set that is the array's response to the input source based on the u-v tracks of our actual observations.
    3. Each data set is then inverted, cleaned, and restored as usual to produce the final output map.
    4. The peak intensity (pi) of each map is compared to that of a point-source (pps) model of equal integrated intensity to give the array's recovery factor, ri, where
      Note that for an infinitely extended source, ri→0 for a point source, and if the source just fills the beam.

    The array's recovery factor is illustrated in this figure. The top panel is for comet LINEAR, and the bottom panel is for comet NEAT. The abscissa is the cometocentric radius, and the ordinate is ri. The solid line shows ri for all shells.
  5. The production rate Q of the comet is found by modifying the column density equation,
    It is noted that 4ln(2)/πab is the projected area on the sky of our Gaussian synthesized beam (in square centimeters), Ni is the beam-averaged column density of the shell, Wi is the integrated line intensity of the shell, and Ti is the rotation temperature for the shell (from Combi et al. 1999). Rearranging equation this equation and incorporating the fact that
    Incorporating equation (3) and rearranging equation (7) gives
    Since all factors on the right-hand side of equation (8) are known, Q can be found. Note that in equation (5) we divided ni by the projected beam area rather than the projected area of the shell, even for those shells with a projected area larger than our beam. This is done because although there are some molecules outside the synthesized beam, they are contributing flux to the area inside the synthesized beam.
  6. Now that Q has been found, we use it in equation (3) to find ni. Then we use each ni in equation (5) to find Ni for each shell. Then the total beamaveraged column density, NT, is given by the sum over all shells,
1In this step we approximate the shells as Gaussian sources rather than as actual shells because any structure on the sky that has sharp edges, such as a shell, will produce ringing in the u‐v plane and artifacts in the final map. Since the actual observations of the comets observed continuous distributions on the sky with no sharp edges, approximating the shell by a Gaussian is reasonable.