Program: XGAUPROF #begin section gettingstarted GETTING STARTED =============== XGAUPROF fits parameters of a number of Gauss related functions to data from a GIPSY set or to data from file. The functions include also fitting of skewness and kurtosis parameters and Lorentzian line widths. Enter the name of a GIPSY set in the 'Setname' input field. Do not specify any axes or ranges yet. If the set exists, you get a 'Profile axis' menu to the right of the 'Setname' field. Select the axis that represents the X-axis of your profile(s). A default grid range will appear. Enter a box. If you press in an empty box field, a default will appear. The first profile in the box will be plotted and if you see a green curve (initial guess curve) you can push the button to fit the profile. #end section gettingstarted #begin section buttons DESCRIPTION OF THE BUTTONS AND INPUT FIELDS =========================================== MAIN GUI Usually you will start the program without specifying any keywords. To get results you will have to enter the name of a GIPSY image first (if you want to import data from other sources, see description at ). Then the 'Profile axis' menu will become active. This menu contains the names of all the axes in your set. After this, the default limits in grid coordinates of your profiles is displayed in the input field 'Grid range'. This value can be changed at any time to exclude unsuitable data. The input field 'Phys. range' corresponds to the grid range. Here you can change the limits in physical coordinates. Pressing in the 'Box' input field will generate a default box. Note that the default generates the maximum possible box, so if you start with a big cube, then it can take a while before all data is read. Finally you get a default in the 'Position' input field. Usually this is the first position in the box and you immediately will see the corresponding profile in the plot at the left of the gui. You can go trough all the selected profiles with the '<' and '>' buttons. If you have initial guesses (usually plotted in green) then you can try to make a fit by pressing the button. Results will appear in the panels just below the input fields. The output of results will be discussed later. It is possible to exclude data by clicking with the left hand mouse button in the neighbourhood of the point that you want to exclude. In the window there is also an option to exclude points by giving the range of points that you want to include! You can set a clip or range to the Y data in the Clip input field in the same window. Values below the clip (1 value entered) are marked in the plot and get weight zero in the fit. If you entered two values, then all values outside this range are marked and get weight zero in the fit. The button is designed to allow the storage of the fit results on disk. It is activated only after the selection of a proper file name in the menu. The button redraws data and initial estimates. This button starts a menu with the options: 1) 'save plot' This option results in sending your current plot to a printer or a PostScript file. In a sub menu you will be asked to give a file name and to select a 'Plotter'. If the file already exists, you will be warned. Writing is started after pressing the button. The file name that you entered, is displayed in the gui and the GIPSY log. 2) 'save profile data' Save current profile data to a file on disk. This can be handy if you want to import your GIPSY data in a package that reads ASCII files. Be sure that you made a fit before saving the data, otherwise the fitted amplitudes and the residuals are dummies. Here is part of an example output file: ! X GRID X PHYSICAL AMPLITUDE FITTED AMPL RESIDUAL -49 816.783 -0.229097 2.63399e-47 -0.229097 -48 820.927 -0.245038 7.42847e-46 -0.245038 ....... ....... 3) 'Fit output file (short)' Open a file to which the fit results can be sent. If you entered a valid name, the button in the main gui is enabled. Pressing this button results in writing the current fit data to the ASCII file until the program is finished. The layout of the file is such that all output of one fit is written onto one line. A header will inform you about the meaning of the numbers. If a fit could not be made, a message is given. 4) 'Fit output file (long)' Same as above, but the results are written in a so called log format. 5) 'Events input file' If you want to process a couple of profiles, it is sometimes easier to store events in a file. There are two methods. The common method is using the Ggi event player which is started by clicking the right mouse button in the background of the main gui. However, this application has also it's own event player and this one works as follows: Create an ASCII file on disk and write the following lines in that file: setname=real profile=freq box=-10 -10 10 10 GUI_FILE=2 PROMPTER_OK=y PROMPTER_OK=y fit=y fstore=y next=y fit=y fstore=y Give it the name events.txt and enter this name in the input field under , 'Events input file'. The text will be read and events are generated immediately. The events are just keywords with a value. The problem is to find which keyword corresponds to a button. If you are running XGAUPROF, type on the Hermes command line: XGAUPROF MONITOR=Y Press a button and you will see a log of the keyword(s) that is/are involved. However you have to distinguish yourself which events are generated by you and which are generated by the system as a reaction on your action. Note that this mechanism is different from the so called GIPSY defaults file. Pressing will pop up a window for importing data from OTHER sources than a GIPSY image. Both reading from text file and reading from a GIPSY table are supported. The syntax for the rows is displayed in the help text of the 'Rows' input field. Example: Read all x and y values from file 'perfil.dat': Filename: perfil.dat Column X: 1 Y: 2 Errors : Rows : : (Pressing in the rows field has the same effect as entering the colon.) Note that in order to get initial guesses for the least squares fit, your Y values must not be too small (> 1e-8). For the least squares fit routine it doesn't matter how the X-values in your data are distributed, however for the initial estimates routine we need equally spaced X-values. Therefore XGAUPROF will eliminate duplicates in X and uses a spline interpolation to create a dummy profile with equidistant X-values and without blanks. If a spline interpolation could not be made (e.g. all Y values blank except one) then no initial estimates could be found. Note that the fit routine works with the raw data. The use of these inputs is straightforward. If you fix the amplitude range, then each profile in the box will get the same range in Y-values (=amplitudes). The button clears the plot and redraws data and estimates. 'amp. %' Set the percentage of the maximum at which the user wants the width. This option works only for standard gaussians because only for this function there is a simple relation between fitted dispersion and the width at arbitrary height. Turn valleys into peaks with this button. The original data is multiplied by -1. In order to restore the original data you have to push the button again. 'tolerance' This window contains two parts. The first part contains parameters which control the least squares fit routine. The tolerance is a stop criterion. its default is zero, which implies that the program determines the minimum value. If the fit routine cannot improve the reduced chi-squared by 'tolerance', it stops and returns the current fit parameters. 'lambda' Lambda is the mixing parameter. it determines the initial weight of steepest descent method relative to the Taylor method. Lambda should be a small value (i.e. 0.01). It can only be zero when the partial derivatives are independent of the parameters which is not the case in this program. In practice, this values will almost never be changed. 'max. Its' The least squares routine fails to produce results if it needs more than this number of iterations to achieve convergence. 'Include range' Enter the range in physical coordinates for which you want to include the data. Other points are masked and get zero weight in the least squares fit. If you enter an empty input then the entire data range is taken into account and the masked points are reset. 'Clip' You can set a clip or range to the Y data. Values below the clip (1 value entered) are marked in the plot and get weight zero in the fit. If you entered two values, then all values outside this range are marked and get weight zero in the fit. The clip is applied to all the profiles not only the current. 'Min. ampl.' Initial estimates with an amplitude smaller than this value are discarded. See also 'Q'. 'Min. disp.' Initial estimates with a dispersion smaller than this value are discarded. The input has the same units as the X-values in the plot. Usually these are physical units and not grids. See also 'Q'. 'Q' This is called the 'smoothing factor'. The initial estimate routine uses 2Q+1 points around a maximum to fit a second order polynomial. The coefficient of the second-order term of the polynomial is an approximation of the second derivative of the observed profile. With a moments method the center and dispersion are stimated. The maximum amplitude is derived from the observed profile. The default is a series of Q's (1 2 3 4 & 5). If you enter more than one Q then the Q with the best (i.e. lowest chi- squared) initial estimates is used for the least squares fit. The estimate filters 'Min. ampl.' and Min. disp.' apply only to this 'best' Q. Note that there is one special case build in the program. If a profile contains blanks or zero-level values except for one, then the estimate routine fails to return an estimate. One can treat this peak as a peak with known amplitude and center, and a full width at half maximum of one pixel (which is automatically converted to a dispersion). However,... the least squares fit routine will not fit such profiles because any dispersion results in the same reduced chi-squared and therefore there is no convergence. 'Median filter' If a profile contains spurious elements, you can use a three pixel median filter to smooth the data a bit. In many cases this is less work than excluding these points by hand. There are three options in the weights menu. First, the default, is uniform weight, i.e. each non blank point in the profile gets weight 1. The reduced chi-square that is calculated is corrected for the current rms of that profile afterwards. Second option is an automatic weighting. The program calculates a second profile that is a three pixel median filter of the original profile. The distance in Y between these two profiles is a measure for a weight of each point. With this method, spurious elements get less weight in the fit, which improves the fit results for some profiles. The third option allows you to use errors from an ASCII file on disk as specified in the menu. If the uncertainty of a data point is given by E then the corresponding weight is equal to 1/(E*E) normalized to the average of all weighting factors (Bevington, Data Reduction and Error Analysis, Ch 10). This window contains input parmaters and fit results. Gaussian parameters: The param window allows you to control the initial estimates as proposed by the program. If the program cannot find any estimates, perhaps you can and this is the place to do that. This is also the window where you tell the fit routine to set gaussian parameters to fixed or free. The more you fix, the better the fit. Background: The background is taken into account in the fit. It has three parameters. A constant, a linear and a quadratic term. In the default situation only the constant term is a free parameter. Slider: Each parameter has a button. If you push this button, you activate the slider for this parameter. While moving the slider, you can watch the parameter value changing in the value field and in the plot. rms: Initial estimates can only be determined if a reasonable value for the rms of the noise is given. If you have a value in mind, push the button to the right of the rms input field and enter your value. Otherwise, a default is calculated. The profile data is sorted and the distance between the quartiles is taken a a measure of rms of the noise. Usually this value is higher than expected, but it works well to find initial estimates. Below right in this window there are two buttons and that have the same function as the buttons with the same labels in the main gui. A menu is presented with four options: 1) A standard gaussian. Parameters are amplitude, center, and dispersion, 2) Gauss-Hermite polynomial (h3). Parameters are a,b,c, (which are NOT the same as amplitude, center and dispersion) and h3. The parameters are translated to the more familiar amplitude, center, dispersion and skewness. 3) Gauss-Hermite polynomial (h3, h4). Same as above, but an extra parameter h4 is included. The value for this h4 parameter is translated to kurtosis. 4) Voigt. Fit Lorentzian and Doppler widths together with center and amplitude. The mathematics behind all this is worked out in a separate document. Select the number of gaussians that you want to fit. Show this document. Push the button again to get rid of the help window. OUTPUT ------ The different functions generate different output, but there are some quantities that are the same for all. First there is output in the status area; the number of iteration to reach convergence, the total flux, reduced chi-square, the standard deviation of the residuals and the current value of the fit tolerance are displayed. In the so called output panel, the fit parameters are listed. These are the parameters that characterize a distribution like the line strength (or area) of the profile, the mean X value (which is only the center if the distribution is symmetric like the standard Gauss and the Voigt profile), Further the maximum amplitude is given and also the values for the background are listed. For the gauss-Hermite series, the skewness and kurtosis are listed. The real function parameters of the Gauss-Hermite functions are not listed in the gui log, but can be found in the GIPSY log. These 'real' function parameters are transformed into properties of a distributions. An extra parameter gives the position of the maximum which is, as mentioned before, for the GH-series not the same as the value of the mean X. A description of the functions and their parameters can be found in a separate document. Here is a summary for the GH-series. The GH-series with h3 and h4: f = a.Exp[-1/2.g^2] * { 1 + h3[c1.g+c3.g^3] + h4[c5+c2.g^2+c4.g^4] } and: g = (x-b)/c c1..c5 have fixed values. The output in short format to file disk has a table header and columns with results. The first line can look like: ! GAUSS fit parameter table created: Tue Jul 6 16:25:04 1999 !RCHI2 |SIG.RES |AREA | err | MEAN | err | DISP. | err | ...... 8.97364 1.04838 1030.28 141.687 2237.89 ...... 8.97364 1.04838 1030.28 141.687 2237.89 ..... ..... etc. etc. The meaning of the titles: RCHI2: reduced chi-squared. This value is scaled by the rms of the noise that you entered. If you did not enter a rms, the value 1.0 is substituted. SIG.RES: Standard deviation of the residuals. AREA: Total flux under fitted function. err: Error on this parameter. MEAN: Mean of gaussian of first component. DISP: Dispersion FW at 50.0 pct: Full width at certain percentage of maximum. AMPL. Amplitude Z0/1/2 Background terms constant, linear and quadratic. RA PARAM Position of profile in grids. These are the axis names of the box. INITIAL ESTIMATES ----------------- In a separate document we describe the method and give a reference to the article of U.J. Schwarz on which the estimate method is based. LEAST SQUARES FIT ----------------- The least squares fit minimizes the reduced chi-square which relates a model profile to the real data. If the decrement in the value of the reduced chi-squared is less than a given tolerance, then the iteration stops and the fitted parameters are returned. Maximum accuracy can be obtained by setting the tolerance to 0 which is also the default value. BLANKS ------ Blanks in the profile get weight zero in the fit. For the estimate routine, blanks are replaced by the current zero level. In the plot window the blanks are indicated by character 'b'. #end section buttons #begin section keywords KEYWORDS ======== CHANGE THE SIZE OF YOUR GUI AT STARTUP -------------------------------------- keyword: GUISIZE= arguments: size in X in pixels size in Y in pixels height of residual window in pixels For building events files you need to know the keywords used by the application. Here is a (incomplete) list with possible keywords with some default values. You can use the Ggi monitor function to find the keywords corresponding to certain events. The monitor tool is started when you click the right mouse button in the background of the main gui. AMPPERC= 50 CHARSIZE= 1.10 COLCON= 8 COLEST= 3 COLFIT= 2 COLFRAME= 1 COLMED= 4 CONNECT= NO CRITAMP= 0 CRITDISP= 0 DRAWZERO= YES ESTFILTER= NO FITFILE= profilefit.dat FUNCTION= 0 GGIOPT= FIXCOL GUI_FILE= HCFILENAME= profilefit.ps HELP= LAB= 0.01 LWIDTH= 2 LWIDTHCUR= 2 MAXITS= 200 NCOMP= 1 PARAM= PLOTDENS= 200 PLOTEST= YES PLOTMED= NO PLOTRMS= NO PLOTTER= 0 PROFFILE= profiledata.dat Q= 1 2 3 4 5 RMS= 1.0 SUBDIV= 1 1 TOLERANCE= 0 VIEWPORT= 0.15 0.1 0.95 0.9 WEIGHTS= 0 Z0_EST= 0 Z1_EST= 0 Z1_FIX= YES Z2_EST= 0 Z2_FIX= YES #end section keywords #begin section about Purpose: Fit parameters of a standard Gauss function, a Gauss-Hermite series or a Voigt lineshape to data from GIPSY sets or ASCII files Category: ANALYSIS, PROFILES, PLOTTING File: xgauprof.src Author: M.G.R. Vogelaar (GUI: J.P. Terlouw) Keywords: See section keywords Version: 2.1 Example: Notes: Related document: xgauprof.ps Updates: Feb 22, 1999: VOG, Document created. Jun 15, 2001: VOG, Help updated, connect lines added. Apr 15, 2009: VOG, Removed unused macro NINT and inserted string.h in import.c and param.c Jun 28, 2010: JPT, Changed log font. #end section about