The CModel approach to model-fit galaxy photometry - also known as the "Sloan Swindle" - is an approximation to bulge+disk or Sersic model fitting that follows the following sequence:

Fit a PSF-convolved elliptical exponential (Sersic n=1) model to the data.
Fit a PSF-convolved elliptical de Vaucouleurs (Sersic n=4) model to the data.
Holding the positions and ellipses of both models fixed (only allowing the amplitudes to vary), fit a linear combination of the two models. In the limit of pure bulge or pure disk galaxies, this approach yields the same results as a more principled bugle+disk or Sersic fit. For galaxies that are a combination of the two components (or have more complicated morphologies, as of course all real galaxies do), it provides a smooth transition between the two models, and the fraction of flux in each of the two parameters is correlated with Sersic index and the true bulge-disk ratio. Most importantly, this approach yieled good galaxy colors in the SDSS data processing.

In this implementation of the CModel algorithm, we actually have 4 stages:

In the "initial" stage, we fit a very approximate PSF-convolved elliptical model, just to provide a good starting point for the subsequence exponential and de Vaucouleur fits. Because we use shapelet/Gaussian approximations to convolved models with the PSF, model evaluation is much faster when only a few Gaussians are used in the approximation, as is done here. In the future, we may also use a simpler PSF approximation in the initial fit, but this is not yet implemented. We also have not yet researched how best to make use of the initial fit (i.e. how does the initial best-fit radius typically relate to the best-fit exponential radius?), or what convergence criteria should be used in the initial fit. Following the initial fit, we also revisit the question of which pixels should be included in the fit (see CModelRegionControl).
In the "exp" stage, we start with the "initial" fit results, and fit an elliptical exponential profile.
In the "dev" stage, we start with the "initial" fit results, and fit an elliptical de Vaucouleur profile.
Holding the "exp" and "dev" ellipses fixed, we fit a linear combination of those two profiles. In all of these steps, the centroid is held fixed at a given input value (take from the slot centroid when run by the measurement framework).

Units

Unlike most measurement algorithms, CModel requires the Exposure it is given to have both a Wcs and a PhotoCalib. This is because it makes use of Bayesian priors, and hence it has to know the relationship between the raw units of the image (pixel and count) and the global units in which the priors are defined.

In fact, all of the nonlinear fits in CModel are done in a special, local coordinate system, defined by a Wcs in which the "pixels" have units of arcseconds (because we never create an image in this system, we don't have to worry about the size of the pixels) and the fluxes should be of order unity. In addition to allowing us to use priors, it also ensures that the parameters all have the same order of magnitude, which improves the behavior of the optimizer.

See Units and Coordinate Systems for more information.

Forced Photometry

In forced photometry, we replace the three nonlinear fits with amplitude-only fits, and then repeat the final linear fit, using the ellipses from the reference catalog in all casees. We do allow the relative amplitudes of the two components to vary in forced mode, though in the future we will add an option to hold this fixed as well as the ellipses.

Shapelet Approximations to the PSF

The CModel algorithm relies on a multi-shapelet approximation to the PSF to convolve galaxy models. It does not compute this approximation directly; for CModelAlgorithm methods that take inputs directly as arguments, the PSF must be supplied as a shapelet::MultiShapeletFunction instance. When using SourceRecords for input/output, CModel assumes that the ShapeletPsfApprox plugin has already been run (see psf.py), and uses the fields created by that plugin to retrieve the PSF approximation.

Code Organization

The CModel implementation consists of many classes, defined in this file and CModel.cc. These mostly fall into four categories:

Control structs: C++ analogs of Python Config classes, these structs contain the configuration parameters that control the behavior of the algorithm. These are nested; the CModelControl struct contains a PixelFitRegionControl instance and three CModelStageControl (one for each of "initial", "exp", and "dev"). The configuration for the final amplitude-only fit goes in CModelControl itself; because it is a simpler linear fit, it doesn't have much in common with the first three stages.
Result structs: while the algorithm has methods to use SourceRecord objects for input/output, it can also take inputs directly as arguments and return the outputs using these structs. Like the Control structs, the master CModelResult struct holds three CModelStageResult classes, for each of the three nonlinear fits.
Keys structs: these private classes (defined in an anonymous namespace in CModel.cc) hold the afw::table::Key and FunctorKey objects that provide a mapping from the Result structs to Schema fields. They also provide methods to transfer values from Results to Records, or the reverse. These are also split into a master CModelKeys struct, which holds three CModelStageKeys structs.
Impl classes: these private classes contain the actual algorithmic code. Once again, we have the master implementation class (CModelAlgorithm::Impl) and a class for the nonlinear fitting stages (CModelStageImpl). In addition to these categories, we also have the CModelAlgorithm class, which is the C++ public interface to all of this, and the CModelStageData class, a private class that aggregates per-source state and makes it easier to pass it around.