diff options
Diffstat (limited to 'noao/onedspec/doc/fitprofs.hlp')
| -rw-r--r-- | noao/onedspec/doc/fitprofs.hlp | 403 |
1 files changed, 403 insertions, 0 deletions
diff --git a/noao/onedspec/doc/fitprofs.hlp b/noao/onedspec/doc/fitprofs.hlp new file mode 100644 index 00000000..ed21e7b1 --- /dev/null +++ b/noao/onedspec/doc/fitprofs.hlp @@ -0,0 +1,403 @@ +.help fitprofs Mar92 noao.onedspec +.ih +NAME +fitprofs -- Fit 1D profiles to features in image vectors +.ih +USAGE +fitprofs input +.ih +PARAMETERS +.ls input +List of input images to be fit. The images may be one dimensional +spectra (one or more spectra per image) or long slit spectra. Other +types of nonspectral images may also be used and for two dimensional +images the fitting direction will be determined from either the keyword +DISPAXIS in the image header or the \fIdispaxis\fR parameter. +.le +.ls lines = "" +List of lines, columns, or apertures to be selected from the input image +format. The default empty list, "", selects all vectors in the images. +The syntax is a list of comma separated numbers or ranges, where a range +is a pair of hyphen separated numbers. +.le +.ls bands = "" +List of bands for 3D images. The empty list, "", selects all bands. +.le +.ls dispaxis = ")_.dispaxis", nsum = ")_.nsum" +Parameters for defining vectors in 2D and 3D images. The +dispersion axis is 1 for line vectors, 2 for column vectors, and 3 for band +vectors. A DISPAXIS parameter in the image header has precedence over the +\fIdispaxis\fR parameter. The default values defer to the package +parameters of the same name. +.le + +The following are the fitting parameters. +.ls region = "" +Region of the input vectors to be fit specified as a pair of space +separated numbers. The coordinates are defined in terms of the linear +image header coordinate parameters. For dispersion corrected spectra this +is usually wavelength in Angstroms and for other data it is usually pixels. +A fitting region must be specified. +.le +.ls positions = "" +File of initial or fixed profile positions and (optional) peaks, profile +types, and widths. The +format consists of lines with one or more whitespace separated fields. +The fields are the position, peak relative to the continuum with +negative values being absorption, profile type of gaussian, lorentzian, +or voigt, and the gaussian and/or lorentzian full width at half maximum. +Trailing fields may be missing and fields to be set from default parameters +or the image data (the peak value) may be given as INDEF. +Comments and any additional columns are ignored. The positions and +widths are specified in the coordinate units of the image, usually +wavelength for dispersion corrected spectra and pixels otherwise. +.le +.ls background = "" +Background values defining the linear background. If not specified the +single pixel values nearest the fitting region endpoints are used. +Otherwise two whitespace separated values are expected. If a value is +a number then that is the background at the lower or upper end of the +fitting region (ordered in pixel space not wavelength). The special +values "avg(w1,w2,z)" or "med(w1,w2,z)" (note that there can be no +whitespace) may be specified, where w1 and w2 are dispersion values, and z +is a multiplier. This will take the average or median of pixels within the +specified range and multiply the result by the third argument. The +dispersion point used for that value in computing the linear background is +the average of the dispersion coordinates of the pixels used. +.le +.ls profile = "gaussian" (gaussian|lorentzian|voigt) +Default profile type to be fit when a profile type is not specified in +the positions file. The type are "gaussian", "lorentzian", or "voigt". +.le +.ls gfwhm = 20., lfwhm = 20. +Default gaussian and lorentzian full width at half maximum (FWHM). +These values are used for the initial and/or fixed width when they are +not specified in the position file. +.le +.ls fitbackground = yes +Fit the background? If "yes" a linear background across the fitting region +will be fit simultaneously with the profiles. If "no" the background will +be fixed. +.le +.ls fitpositions = "all" +Position fitting option. This may be "fixed" to fix all positions at their +initial values, "single" to fit a single shift to the positions while +keeping their separations fixed, or "all" to independently fit all the +positions. +.le +.ls fitgfwhm = "all", fitlfwhm = "all" +Profile width fitting options. These may be "fixed" to fix all widths +at their initial values, "single" to fit a single scale factor to the initial +widths, or "all" to independently fit all the widths. +.le + +The following parameters are used for error estimates as described +below in the ERROR ESTIMATES section. +.ls nerrsample = 0 +Number of samples for the error computation. A value less than 10 turns +off the error computation. A value of ~10 does a rough error analysis, a +value of ~50 does a reasonable error analysis, and a value >100 does a +detailed error analysis. The larger this value the longer the analysis +takes. +.le +.ls sigma0 = INDEF, invgain = INDEF +The pixel sigmas are modeled by the formula: + +.nf + sigma**2 = sigma0**2 + invgain * I +.fi + +where I is the pixel value and "**2" means the square of the quantity. If +either parameter is specified as INDEF or with a value less than zero then +no sigma estimates are made and so no error estimates for the measured +parameters is made. +.le + +The following parameters determine the output of the task. +.ls components = "" +All profiles defined by the position file are simultaneously fit but only +a subset of the fitted profiles may be selected for output. A profile +or component is identified by the order number in the position file; +i.e. the first entry in the position file is 1, the second is 2, etc. +The components to be output are specified by a range list. The empty +list, "", selects all profiles. +.le +.ls verbose = yes +Print fitting results and record of output images created on the +standard output (normally the terminal). +The fitting information is printed to the logfile so there is normally +no need to redirect this output. The output may be turned off when +the task is run as a background task. +.le +.ls logfile = "logfile" +Logfile for fitting results. If not specified the results will not be +logged. +.le +.ls plotfile = "plotfile" +File to contain plot output. The plots show the image vector with +overplots of the total fit, the individual components, and the residuals. +The plotfile may be examined and manipulated later with tools such as +\fBgkimosaic\fR. +.le +.ls output = "" +List of output images. If not specified then no output images are created. +If images are specified the list is matched with the input list. +.le +.ls option = "fit" (fit|difference) +Image output option. The choices are "fit" to output the fitted image +vector which is the sum of the fitted profiles (without a background), +or "difference" to output the data with the profiles subtracted. +.le +.ls clobber = no, merge = no +Clobber or modify any existing output images? If clobbering is not +enabled a warning is printed and any existing output images are not +modified. If clobbering is enabled then either new images are created +if merge is "no" or the new fits are merged with the existing images. +Merging is meaningful when only a subset of the input is fit such +as selected lines or apertures. +.le +.ih +DESCRIPTION +\fBFitprofs\fR fits one dimensional profile functions to image vectors +and outputs the fitting parameters, plots, and model or residual +image vectors. This is done noninteractively using a file of initial +profile positions and widths. Interactive profile fitting may be +done with the deblending option of \fBsplot\fR or +\fBstsdas.fitting.ngaussfit\fR. + +The input consists of images in a variety of formats. These include +all the spectral formats as well as standard images. For two dimensional +images (or the first 2D plane of higher dimensional images) either the +lines or columns may be fit with possible summing of adjacent lines or +columns to increase the signal-to-noise. A subset of the image apertures, +lines, or columns may be specified or all image vectors may be fit. + +The fitting parameters consist of a fitting region, a list of initial +positions, peaks, and widths, initial background endpoints, the fitting +function, and the parameters to be fit or constrained. The coordinates and +units used for the positions and widths are those defined by the standard +linear coordinate header parameters. For dispersion corrected spectra +these are generally wavelengths in Angstroms and otherwise they are +generally pixels. A fitting region must be specified by a pair of +numbers. + +The background parameter may be left empty to select the pixel values at +the endpoints of the fitting region for defining the initial linear +background. Or values at the endpoints of the fitting region may be given +explicitly in pixel space order (i.e. the first value is for the edge of +the fitting region which has smaller pixel coordinate0 Values can also be +computed from the data using the functions "avg(w1,w2)" or "med(w1,w2)" +where w1 and w2 are dispersion coordinates. The pixels in the specified +range are average or medianed and the dispersion point for the linear +background is the average of the dispersion coordinates of the pixels. + +The position list file consists of one or more columns. +The format of this file has +one or more columns. The columns are the wavelength, the peak value +(relative to the continuum with negative values being absorption), +the profile type (gaussian, lorentzian, or voigt), and the +gaussian and/or lorentzian FWHM. End columns may be missing +or INDEF values may be specified to use the default parameter +values (the profile and widths) or determine the peak from the data. +Below are examples of the file line formats + +.nf + wavelength + wavelength peak + wavelength peak (gaussian|lorenzian|voigt) + wavelength peak gaussian gfwhm + wavelength peak lorentzian lfwhm + wavelength peak voigt gfwhm + wavelength peak voigt gfwhm lfwhm + + 1234.5 <- Wavelength only + 1234.5 -100 <- Wavelength and peak + 1234.5 INDEF v <- Wavelength and profile type + 1234.5 INDEF g 12 <- Wavelength and gaussian FWHM +.fi + +where peak is the peak value, gfwhm is the gaussian FWHM, and lfwhm is +the lorentzian FWHM. This format is the same as used by \fBsplot\fR +and also by \fBartdata.mk1dspec\fR (except in the latter case the +peak is normalized to a continuum of 1). + +The profile parameters fit are the central position, the peak amplitude, +and the profile widths. The fitting may be constrained in number of ways. +The linear background may be fixed or simultaneously fit with the +profiles. The profile positions may be fixed, the relative separations +fixed but a single zero point shift fit, or all positions may be fit +simultaneously. The profile widths may also be fixed, the relative ratios +of the widths fixed while fitting a single scale factor, or all widths fit +simultaneously. The profile amplitudes are always fit. + +The fitting technique uses a nonlinear iterative Levenberg-Marquardt +algorithm to reduce the Chi-square of the fit. The execution time +increases rapidly with the number of profiles fit so there is an +effective limit to the number of profiles that can be fit at once. + +The output includes a number of formats. The fitted parameters are +recorded in a logfile (if specified) and printed on the standard +output (if the verbose flag is set). This output includes the date, +image vector, fitting parameters used, and a table of fitted or +derived quantities. The parameters included some quantities relevant to +spectral lines but others apply to any image data. The quantities are +the profile center, the background or continuum at the center of the +profile, the integral or flux of the profile (which is negative for +profiles below the background), the equivalent width, the profile peak +amplitude or core value, and the profile full width at half +maximum. Pure gaussian and lorentzian profiles will have one of +the widths set to zero while voigt profiles will have both values. + +Summary plots are recored in a plotfile (if specified). The plots +show the data with the total fit, individual profiles, and residuals +overplotted. The plotfile may be examined and printed using the +task \fBgkimosaic\fR as well as other tasks which interpret GKI metacode. + +The final output consists of images in the same format as the input. +The images may be of the total fit (sum of profiles without background) +or of the difference (residuals) of the data minus the model. +.ih +ERROR ESTIMATES +Error estimates may be computed for the fitted parameters. +This requires a model for the pixel sigmas. Currently this +model is based on a Poisson statistics model of the data. The model +parameters are a constant Gaussian sigma and an "inverse gain" as specified +by the parameters \fIsigma0\fR and \fIinvgain\fR. These parameters are +used to compute the pixel value sigma from the following formula: + +.nf + sigma**2 = sigma0**2 + invgain * I +.fi + +where I is the pixel value and "**2" means the square of the quantity. + +If either the constant sigma or the inverse gain are specified as INDEF or +with values less than zero then no noise model is applied and no error +estimates are computed. Also if the number of error samples is less than +10 then no error estimates are computed. Note that for processed spectra +this noise model will not generally be the same as the detector readout +noise and gain. These parameters would need to be estimated in some way +using the statistics of the spectrum. The use of an inverse gain rather +than a direct gain was choosed to allow a value of zero for this +parameters. This provides a model with constant uncertainties. + +The error estimates are computed by Monte-Carlo simulation. The model is +fit to the data (using the noise sigmas) and this model is used to describe +the noise-free spectrum. A number of simulations, given by the +\fInerrsample\fR, are created in which random Gaussian noise is added to +the noise-free spectrum based on the pixel sigmas from the noise model. +The model fitting is done for each simulation and the absolute deviation of +each fitted parameter to model parameter is recorded. The error estimate +for the each parameter is then the absolute deviation containing 68.3% of +the parameter estimates. This corresponds to one sigma if the distribution +of parameter estimates is Gaussian though this method does not assume +this. + +The Monte-Carlo technique automatically includes all effects of +parameter correlations and does not depend on any approximations. +However the computation of the errors does take a significant +amount of time. The amount of time and the accuracy of the +error estimates depend on how many simulations are done. A +small number of samples (of order 10) is fast but gives crude +estimates. A large number (greater than 100) is slow but gives +very good estimates. A compromise value of 50 is recommended +for many applications. + +.ih +EXAMPLES +1. The following example creates an artificial spectrum and fits it. +It requires the \fBartdata\fR and \fBproto\fR packages be loaded. + +.nf + cl> mk1dspec test slope=1 temp=0 lines=testlines nl=20 + cl> mknoise test rdnoise=10 poisson=yes + cl> fields testlines fields=1,3 > fitlines + cl> fitprofs test reg="4000 8000" pos=fitlines + # Jul 27 17:49 test - Ap 1: + # Nfit=20, background=YES, positions=all, gfwhm=all, lfwhm=all + # center cont flux eqw core gfwhm lfwhm + 6832.611 1363.188 -13461.8 9.875 -408.339 30.97 0. + 7963.674 1507.641 -8193.58 5.435 -395.207 19.48 0. + 5688.055 1217.01 -7075.11 5.814 -392.006 16.96 0. + 6831.3 1363.02 -7102.01 5.21 -456.463 14.62 0. + 7217.335 1412.323 -10110. 7.158 -427.797 22.2 0. + 6709.286 1347.437 -4985.06 3.7 -225.346 20.78 0. + 6434.317 1312.319 -7121.03 5.426 -342.849 19.51 0. + 6130.415 1273.506 -6164. 4.84 -224.146 25.83 0. + 4569.375 1074.138 -3904.6 3.635 -183.963 19.94 0. + 5656.645 1212.999 -8202.81 6.762 -303.617 25.38 0. + 4219.53 1029.458 -5161.64 5.014 -241.135 20.11 0. + 4551.424 1071.845 -3802.61 3.548 -139.39 25.63 0. + 4604.649 1078.643 -5539.15 5.135 -264.654 19.66 0. + 6966.557 1380.294 -11717.5 8.489 -600.581 18.33 0. + 4259.019 1034.501 -4280.38 4.138 -213.446 18.84 0. + 5952.958 1250.843 -8006.98 6.401 -318.313 23.63 0. + 4531.89 1069.351 -712.598 0.6664 -155.197 4.313 0. + 7814.418 1488.579 -2926.49 1.966 -164.891 16.67 0. + 5310.929 1168.846 -10132.2 8.669 -487.502 19.53 0. + 5022.948 1132.066 -7532.8 6.654 -325.594 21.73 0. + +.fi + +2. Suppose there is no obvious continuum level near the fitting +region but you want to specify a flat continuum level as the average +of pixels in a specified wavelength region. The background region +would be specified as + +.nf + background = "avg(4250,4425.3) avg(4250,4425.3)" +.fi + +Note that the value must be given twice to get a flat continuum. +.ih +REVISIONS +.ls FITPROFS V2.11.3 +Modified to allow a more general specification of the background. +.le +.ls FITPROFS V2.11 +Modified to include lorentzian and voigt profiles. The parameters and +positions file format have changed in this version. A new parameter +controls the number of Monte-Carlo samples used in the error estimates. +.le +.ls FITPROFS V2.10.3 +Error estimates based on a simple noise model are now computed. +.le +.ls FITPROFS V2.10 +This task is new. +.le +.ih +TIME REQUIREMENTS +The following CPU times were obtained with a Sun Sparcstation I. The +number of pixels in the fitting region and the number of lines fit +were varied. The worst case of fitting all parameters and a background +was considered as well as the constrained case of fitting line positions +and a single width with fixed background. + +.nf + Npixels Nprofs Fitbkg Fitpos Fitsig CPU(sec) + 100 5 yes all all 1.9 + 100 10 yes all all 3.3 + 100 15 yes all all 5.6 + 100 20 yes all all 9.0 + 512 5 yes all all 4.7 + 512 10 yes all all 10.0 + 512 15 yes all all 17.6 + 512 20 yes all all 27.8 + 1000 5 yes all all 8.0 + 1000 10 yes all all 18.0 + 1000 15 yes all all 31.8 + 1000 20 yes all all 50.2 + 1000 25 yes all all 72.8 + 1000 30 yes all all 100.2 + 512 5 no all single 2.8 + 512 10 no all single 5.3 + 512 15 no all single 8.6 + 512 20 no all single 12.8 +.fi + +Crudely this implies CPU time goes as the 1.4 power of the number of profiles +and the 0.75 power of the number of pixels. +.ih +SEE ALSO +splot, stsdas.fitting.ngaussfit +.endhelp |