Model Fitting

The properties of an ensemble of faint galaxies can be defined by the distribution of the parameters of the models that best fit each individual object. This is only practical for faint images if the number of adjustable parameters to be fitted is kept to a minimum; a possible initial set would include total magnitude, half-light radius, ellipticity, orientation and centroid position. For brighter images, where more signal is available, we can fit additional shape parameters, such as bulge-to-disk ratio. The model also includes instrumental parameters, such as the sky background (estimated locally) and amplitude of the noise. The PSF and the distribution of faint cosmic rays are measured separately. It should be pointed out that the modelling approach for extended objects is relatively insensitive to small-scale errors in the description of the PSF.

For any set of model parameters, a smooth galaxy image is generated and convolved with the expected PSF. The expected sky brightness is added and the image is multiplied by the detector sensitivity (flat field). The probability of the observed value is then computed for each pixel using the adopted error distribution. Since the parameters are estimated as the maximum likelihood model that fit the observations, we can ignore defective pixels without any need to interpolate them. The likelihood function is defined as the sum of the natural logarithm of the probability for individual pixels, including both galaxy pixels and neighboring sky pixels (in fact, pixels do not have to be separated into `galaxy' and `sky'). For a Gaussian error distribution, the likelihood function is proportional to the of the distribution; however, the method allows for a different error distribution can be used, to take into account both the uneven probability distribution of different data numbers due to the ADC error, and the wide positive tail due to faint cosmic rays and warm pixels.

The choice of models is based on standard galaxy models and include an exponential disk and an law. Because the method involves a non-linear optimization, the choice of initial guesses for the parameters is quite important; we generally use initial values estimated from moments of the image, as well as the results of the deconvolution if available. Optimization uses the Numerical Algorithm for Maximum Likelihood Estimation (NAMaLiE) developed by Ratnatunga and Casertano (1991). The information on errors of the estimated parameters and their correlation can be obtained from the covariance matrix, which is the inverse of the Hessian (matrix of second order derivatives) at the peak of the likelihood function. Extensive numerical simulations have shown that both parameter estimates and their errors are unbiased.

This approach also requires an error map for the calibrated image, which is not generated by the current HST WFPC pipeline calibration. We use a simple noise model to estimate the rms noise in each pixel. Special consideration must be given to noise in the calibration images (bias, dark and flat field), which, unlike image noise, does not decrease as multiple exposures are stacked together, unless the exposures are shifted with respect to one another. In some cases involving multiple exposures, the noise could be estimated from the internal scatter between individual exposures, plus the contribution of the calibration files.

The contribution of faint cosmic rays to the noise is found to be described adequately by a Weibull distribution with power law index of 0.25; convolved with a Gaussian distribution for shot and read noise, this is a good representation for the total `noise' distribution over the full range in data numbers. The two additional parameters required to describe the distribution of cosmic rays are derived from a single fit for the whole CCD, and therefore do not increase the number of parameters to be fitted for each image. Empirically, the scatter between different exposures is found to be small.