Skip to content

Commit

Permalink
smurf: streamlined parameter file handling and activated clip option.
Browse files Browse the repository at this point in the history
(cherry picked from commit 4b40efe)
  • Loading branch information
Remo Tilanus authored and MalcolmCurrie committed Aug 24, 2012
1 parent 9ed097d commit 921dd83
Show file tree
Hide file tree
Showing 4 changed files with 279 additions and 133 deletions.
33 changes: 19 additions & 14 deletions applications/smurf/defaults/smurf_fit1d.def
Expand Up @@ -27,31 +27,36 @@
#
# IMPORTANT: for gausshermites and voigts the first three parameters
# (referred to as the amplitude, position, and width) do *NOT*
# correspond to the maximum, centre, and FWHM as for a basic Gaussian!
# correspond to the maximum, centre, and FWHM as for a basic Gaussian!
# Read the "smurfhelp fit1d" section on "Fitting_Functions" to see the
# relations.
# relations.

# Axis the fit along (nr starting at 1). '0' translates fit along
# highest dimension i.e. Vlsr in a Ra, Dec, Vlsr cube.

axis = 0

# Coordinate range along axis to find and fit profiles. The format is
# (x1, x2) including the (). E.g. Vlsr -20 to 35 km/s: (-20,35).
# Coordinate range along axis to find and fit profiles. The format is
# (x1, x2) including the (). E.g. Vlsr -20 to 35 km/s: (-20,35).
# Default is to use the full extent of the axis.

range = <undef>

# Values in the input profiles outside the specified clip-range will
# be not be used in the fit. The format is (min, max) including the ().

clip = <undef>

# Function to fit. See documentation: currently implemented are gaussian,
# gausshermite1, gausshermite2, voigt.

function = gaussian
function = gausshermite2

# Maximum number of 'component' functions to fit to each profile, e.g.
# a multi-component spectrum of max. 3 gaussians. Note: The complete fit
# of the gaussians is done concurrently, not iteratively starting e.g.
# with the strongest component. The routine will try to find and fit ncomp
# functions along each profile, but may settle for less.
# of the gaussians is done concurrently, not iteratively starting e.g.
# with the strongest component. The routine will try to find and fit ncomp
# functions along each profile, but may settle for less.

ncomp = 3

Expand All @@ -63,11 +68,11 @@ ncomp = 3
#
# Sorting can be helpful, but be cautioned that it can also complicate
# things: if there are two components one at -10 km/s and one at 10 km/s
# sorting by amplitude or width can result in the parameter file for
# sorting by amplitude or width can result in the parameter file for
# component 1 to be a mix of the -10 and 10 features depending on which
# one was relatively stronger or wider. Similarly, sorting by position
# can result in low-amplitude fits of noise spikes to be mixed with
# stronger components. For more precise control try to run the routine
# stronger components. For more precise control try to run the routine
# iteratively with e.g. a different restricted velocity range to try pick
# out the different components. Default sorting is by amplitude.

Expand All @@ -84,14 +89,14 @@ sort = amp

minamp = 3

# Minimum value for the width-like parameter (see at top) to accept as
# Minimum value for the width-like parameter (see at top) to accept as
# a genuine fit as FWHM (~2.35*Dispersion) in terms of ==PIXELS== (!).
# This parameter can make a huge difference to the number of spurious or
# rejected fits in the end result.

minwidth = 1.88

# Maximum value for the FHWM of the Loretzian component of a voigt in terms
# Maximum value for the FHWM of the Loretzian component of a voigt in terms
# of ==PIXELS== (!). A negative value implies no constraint.
# Voigt fits only.

Expand All @@ -100,7 +105,7 @@ maxlorz = <undef>
# Limits on the Gausshermite parameters h3 (~skew) and peakiness (~h4).
# High values here typically mean that multi-component fit may be better
# or,in case of h3, large negative features. It is possible for the
# fitting routine to compensate for a negative wing (due to h3) with a
# fitting routine to compensate for a negative wing (due to h3) with a
# secondary positive component. I.e. unless a fit of negative wings is
# desired it is better to restrict h3. Small values for h3 and h4 are
# probably better left fixed at 0. Whenever a hermite1 or hermite2 is
Expand All @@ -118,7 +123,7 @@ absh4max = 0.35
# fixed at 0. This means that the fit resorts to a simpler function:
# gausshermite2 -> gausshermite1 -> gaussian; voigt -> gaussian.
# The result is that there are valid fits for more profiles, but the
# fitted function may vary.
# function actually fitted may vary with position.
# Setting the retry value to 0 prevents this from happening and may
# cause the fit to fail or be (very) poor.

Expand Down
12 changes: 5 additions & 7 deletions applications/smurf/libsmf/smf_fit_profile.c
Expand Up @@ -194,8 +194,7 @@ static int getestimates( const double fdata[], const float weight[], int ndat,
static int dolsqfit( smf_math_function fid, const double pcoord[],
const double fdata[], float *weight, int npts,
double *parlist, double *errlist, const int fixmask[],
int npar, int *ncomp, void *pfcntrl, int *fitopt,
int *status );
int npar, int *ncomp, void *pfcntrl, int *fitopt );

static void adjustestimates( smf_math_function fid, int nfound,
double *parlist, int npar );
Expand Down Expand Up @@ -917,7 +916,7 @@ static void FitProfileThread ( void *job_data_ptr, int *status ) {
#endif
iters = dolsqfit( SMF__MATH_GAUSS, pcoord, fdata, weight, npts,
parlist2, errlist2, fixmask2, npar2,
&nfound2, fcntrl, fitopt, status );
&nfound2, fcntrl, fitopt );
if ( iters >= 0 && nfound2 > 0 ) {
for ( i = 0; (int) i < nfound2-1; i++ ) {
int offset = i*npar;
Expand Down Expand Up @@ -950,7 +949,7 @@ static void FitProfileThread ( void *job_data_ptr, int *status ) {
#endif
int *fixmask = fcntrl->fixmask;
iters = dolsqfit( fid, pcoord, fdata, weight, npts, parlist, errlist,
fixmask, npar, &nfound, fcntrl, fitopt, status );
fixmask, npar, &nfound, fcntrl, fitopt );

if ( iters < 0 ) {
nfound = 0;
Expand Down Expand Up @@ -1121,8 +1120,7 @@ static void FitProfileThread ( void *job_data_ptr, int *status ) {
static int dolsqfit( smf_math_function fid, const double pcoord[],
const double fdata[], float *weight, int npts,
double *parlist, double *errlist, const int fixmask[],
int npar, int *ncomp, void *pfcntrl, int *fitopt,
int *status )
int npar, int *ncomp, void *pfcntrl, int *fitopt )
/*------------------------------------------------------------*/
/* PURPOSE: Fit profiles. Wrapper routine for smf_lsqfit. */
/* Return the number of iteration or a number < 0 for an */
Expand Down Expand Up @@ -1806,7 +1804,7 @@ static int fillfromparndf( void *pfcntrl, const smfArray *pardata, int pbase,
/*------------------------------------------------------------*/
{

int iters, pfound; /* parndf: fit error code and nr fitted */
int pfound; /* parndf: nr fitted */
double *fitval, *fiterr; /* Pointer to supplied data */
size_t index; /* Current index into array */

Expand Down

0 comments on commit 921dd83

Please sign in to comment.