Warning: This program is not maintained anymore. In the near future, it will be removed from the set of GIPSY tasks. An improved version with graphical user interface is called XGAUPROF. Program: GAUPROF Purpose: GAUPROF examines influence of parameters in GAUFIT by plotting Gaussian estimate and fit of a selected profile. Category: ANALYSIS, PROFILES, PLOTTING File: gauprof.c Author: M.G.R. Vogelaar Keywords: STARTNEW= Start improved version of this program? [Y]/N A new version with graphical user interface is called XGAUPROF. You can start it here. INSET= Give set name and operation axis ( & range ): Maximum number of subsets is 2048. ** SUBDIV= View surface subdivisions x,y: [1,1] It is possible to have more plots on the view surface. The number of plots can vary in both directions x and y. GRDEVICE= Graphics device [list of all graphics devices] ** PAPER= Give width (cm), aspect ratio: [0.0,1.0] Change the size of the (output) view surface. The aspect ratio is defined as height/width. The default is a calculated size and aspect ratio 1. ** LINEWIDTH= Give line width (1-21): [2] Only required if a hardcopy device is selected. ** XMINMAX= Give Xmin, Xmax in plot: [calculated] ** YMINMAX= Give Ymin, Ymax in plot: [calculated] BOX= Enter profile position(s) in ..., ...: [quit] Examples: BOX=3 4 Examine profile at (3,4) BOX=0 0 3 3 Examine profiles from (0,0) to (3,3) or: BOX=D 3 3 CPOS=0 0 Examine profiles from (-1,-1) to (1,1) PROFILE= (Q)uit, (S)ame, (P)revious, (B)ox, [N]ext profile 1) (Q)uit program 2) (S)ame profile (you can change parameters now) 3) (P)revious profile if it exists. 4) (N)ext profile Next keywords become hidden after their first specification. NGAUSS= Give max. number of gaussians in profile: [1] The maximum number is 5. The calculation of the initial estimates. ========================================= The fit routine in GAUFIT needs initial estimates for the parameters. These estimates are calculated by the program using a method described by Schwarz, 1968, Bull.Astr.Inst. Netherlands, Volume 19, 405. The associated routine reads a number of parameters. It uses NGAUSS=, and the hidden keywords: ESTRMS=, ESTCUTAMP= and ESTCUTDISP=. ** ESTRMS= Give rms noise level in profile: [0.0] Rms of noise in profile in units of the data. The value is used in an automatic window method to define the signal region in the routine that calcu- lates initial estimates for each profile. ** ESTCUTAMP= Give cutoff for profile: [0.0] Estimated gaussians below this amplitude will be discarded. The amplitude is in units of the data. ** ESTCUTDISP= Give critical dispersion of gaussian: [0.0] Estimated gaussians with dispersion smaller than the critical dispersion will be discarded. The critical dispersion is given in units of the operation axis. Q= Give smoothing parameter: [2] See description. Calculation of the fitted parameters. ==================================== The actual fitting is a least-squares fit of a function to a set of data points. The method used is described in: Marquardt, J.Soc.Ind.Appl.Math. 11, 431 (1963). This method is a mixture of the steepest descent method and the Taylor method. The weights for the data points are all equal to 1.0. TOLERANCE= Fractional tolerance for the chi square: [0.01] Fitting is stopped when successive iterations fail to produce a decrement in reduced chi-squared less than TOLERANCE= The value cannot be less than a certain minimum as set by the system. This means that maximum accuracy can be obtained with TOLERANCE=0.0 ** MAXITS= Maximum number of iterations in fit: [50] If this number is exceeded, no parameters could be fitted. ** LAB= Value for mixing parameter: [0.001] Mixing parameter in the least squares fit function. LAB= determines the initial weight of steepest descent method relative to the Taylor method. LAB= should be a small value (i.e. 0.001). Edge fitting ============ ** CENTVAL= Central value for edge fit in physical coordinates: [Central value on profile axis] The position in physical units, needed to determine in which direction data along profile is stored to be used in a least squares fit, now called 'edge' fit. ** EDGEFIT= Do you want an 'edge' fit? Y/[N] Default: skip the 'edge 'fitting'. Description: The decomposition of profiles into Gaussian components can be a powerful method in the study of the kinematics of neutral hydrogen. The program GAUFIT can make initial estimates of Gaussian parameters that can be stored or used as initial values for a least-squares fit. This program (GAUPROF) helps you to examine the influence of different variables like the ones specified with the keywords NGAUSS=, ESTCUTAMP=, ESTCUTDISP=, Q= and TOLERANCE=. There is a loop facility that enables you to change keywords for the same profile. In the observed profile the second derivative of a fitted second order polynomial is used to estimate the parameters of the Gaussian components. The polynomial is fitted at each pixel using q points distributed symmetrically around any pixel. The number q equals 2*Q+1 where Q is given by the smoothing parameter (Q=). The threshold values given by ESTRMS= (rms noise of the profile), and Q=, are used to discriminate against spurious components. Other criteria to select valid Gaussians are a critical amplitude (ESTCUTAMP=) and a critical dispersion (ESTCUTDISP=) in units selected by the program. Gaussians with amplitude smaller than ESTCUTAMP= or with dispersion smaller than ESTCUTDISP= will be discarded. It is assumed that the observed profile can be approxi- mated by the sum of a few Gaussian functions. The number of Gaussians you allow in the calculations is set by NGAUSS=. The maximum number is 5. (For more details: Schwarz, 1968, Bull.Astr.Inst. Netherlands, Volume 19, 405-413.) The fit of the estimated set of components could be improved by means of a least-squares analysis (Full line in the plot). Fitting of the profile stops when successive iterations fail to produce a decrement in reduced chi-squared less than TOLERANCE=. If its value is less than the minimum tolerance possible, it will be set to this value. This means that maximum accuracy can be obtained by setting TOLERANCE=0.0. The weights for all data points is set to 1.0. In order to fit Gaussians in profiles you first need to specify the input set and the subsets (INSET=). The number of subsets must be greater than 1 and less than 2048 and the subsets are specified by giving the name of the axis in the direction in which you want your profiles. Suppose your set AURORA has the axes: RA-NCP from -7 to 8 DEC-NCP from -7 to 8 FREQ-OHEL from 1 to 59 and you want an estimate and fit of Gaussian components of a profile in frequency direction at RA=2, DEC=4, specify: INSET=AURORA FREQ 1:59 BOX=2 4 The 'edge' fitting routine selects data from the position (in physical coordinates) of the estimated Gaussian with maximum amplitude ('center') to a profile edge and performs a 'lsq' fit. The edge is chosen in the following way: If the value in CENTVAL= is smaller than 'center', all data with axis values (in physical coordinates) greater than center are assembled for the 'lsq' fit etc. Sometimes the amplitude of the edge fit is larger than the maximum value of the used data. To adjust your plot use XMINMAX= and YMINMAX=. These keywords are hidden, so you have to specify them for instance at the same time as PROFILE= in the plot loop. Example: gauprof INSET=m8320 f 1:58 Set M8320 has 3 axes RA-NCP from -50 to 40 DEC-NCP from -50 to 20 FREQ-OHEL from 1 to 59 GRDEVICE=tek BOX=-3 -3 4 4 BOX range for set M8320 : RA-NCP from -3 to 4 DEC-NCP from -3 to 4 NGAUSS=2 ESTRMS=0.3 ESTCUTAMP=0.6 ESTCUTDISP=0.4 Q=2 CENTVAL= GAUPROF: ESTIMATE FIT =======================================================... gauss ampl. centre disp. ampl. centre ... W.U. KM/S KM/S W.U. KM/S ... =======================================================... 1 2.89 212.01 8.61 3.43 214.34... 2 0.45 192.04 3.56 0.61 192.63... NGAUSS=2 RMS=0.30 CUTAMP=0.60 CUTDISP=0.40 Q=... PROFILE=n 1 0.88 225.13 3.69 2.54 217.01... 2 0.58 192.50 3.49 0.42 195.15... PROFILE=q gauprof - +++ FINISHED +++ Updates: Jun 11, 1991: VOG, Document created. Oct 10, 1991: WZ, PGPLOT standard names implemented Mar 25, 1994: VOG, Rearranged keywords. May 30, 1995: VOG, Copied 'func' and 'derv' from GAUFIT. Added MAXITS= and LAB= keywords. Updated documentation