BVLS_LH | Bounded value least-square minimization (BVLS) |
BVLS_PS | Bounded value least-square minimization (BVLS) |
MREGRESS | Perform a multiple linear regression fit |
SLICES | Slice an image in horizontal/vertical direction |
TICKSLOG | Compute values for the ticks of a logarithmic axis |
ULY_HTML_HELP | Generate HTML documentation from the embedded doc |
ULY_LINE | Define a spectral-line component |
ULY_STAR | Define a star component |
ULY_STAR_INIT | Initialize a single star fit component |
UNDEFINE | Destroy a variable |
XREBIN | Rebin a vector |
NAME: BVLS_LH PURPOSE: Bounded value least-square minimization (BVLS) USAGE: bvls_lh, a, b, bnd, x, STATUS=status, /DOUBLE, /QUIET RNORM=rnorm, NSETP=nsetp, W=w, INDEX=index, ITMAX=itmax, EPS=eps ARGUMENTS: A On entry A contains the M by N matrix, A. On return A contains the product matrix, Q*A, where Q is an M by M orthogonal matrix generated by this subroutine. B On entry B contains the M-vector, B. On return, B contains Q*B. The same Q multiplies A. BND Input, unchanged on return BND[0,J] is the lower bound for X[J]. BND[1,J] is the upper bound for X[J]. Require: BND[0,J] <= BND[1,J]. X On return, X will contain the solution N-vector. KEYWORDS: [Optional input] ITMAX=itmax, Set to 3*N. Maximum number of iterations permitted. This is usually larger than required. EPS=eps, Determines the relative linear dependence of a column vector for a variable moved from its initial value. The default value is EPS=(machar()).eps (single precision). Other values, larger or smaller may be needed for some problems. /DOUBLE, Setting this keyword forces some computations to be performed in double precision. This may be required if: (i) N is large (>10^6) or (ii) if the A or B values span a large amplitude (e.g. 10^20). This slows the computation, and is not required in most cases. /QUIET Verbosity control. By default the procedure prints the messages. KEYWORDS: [Optional output] RNORM=rnorm, Euclidean norm of the residual vector, b - A*X. NSETP=nsetp, Indicates the number of components of the solution vector, X, that are not at their constraint values. W=w, An N-array. On return, W() will contain the dual solution vector. Using the following definitions: W[J] = 0 for all j in Set P (i.e. X not at a limit), W[J] <= 0 for all j in Set Z, such that X[J] is at its lower bound, W[J] >= 0 for all j in Set Z, such that X[J] is at its upper bound. If BND[0,J] = BND[1,J], so the variable X[J] is fixed, then W[J] will have an arbitrary value. INDEX=index, An INTEGER N-array. On exit the contents of this array define the sets P, Z, and F as follows: INDEX[0] through INDEX[NSETP-1] = Set P. INDEX[NSETP] through INDEX[IZ2] = Set Z. INDEX[IZ2+1] through INDEX[N-1] = Set F. Any of these sets may be empty. Set F is those components that are constrained to a unique value by the given constraints. Sets P and Z are those that are allowed a non-zero range of values. Of these, set Z are those whose final value is a constraint value, while set P are those whose final value is not a constraint. The value of IZ2 is computable as the number of bounds constraining a component of X uniquely. STATUS=status, Indicates status on return. 0 Solution completed. 1 M <= 0 or N <= 0 2 B, X, BND, W, or INDEX size or shape violation. 3 Input bounds are inconsistent. 4 Exceed maximum number of iterations. DESCRIPTION: Given an M by N matrix, A, and an M-vector, B, compute an N-vector, X, that solves the least-squares problem A *X = B subject to X[J] satisfying BND[0,J] <= X[J] <= BND[1,J] Values BND[0,J] <= -(machar()).xmax and BND[1,J] >= (machar()).xmax designate that there is no constraint in that direction. Note that this flags are the single precision largest values, therefore double precision bounds values exceeding this range will be taken as no bound. Non finite values in BND are also interpreted as absence of constraint in the corresponding direction. AUTHOR: R. J. Hanson (1995 April) HISTORY: This algorithm is a generalization of NNLS, that solves the least-squares problem, A * X = B, subject to all X[J] >= 0. The subroutine NNLS appeared in 'SOLVING LEAST SQUARES PROBLEMS', by Lawson and Hanson, Prentice-Hall, 1974. Work on BVLS was started by C. L. Lawson and R. J. Hanson at Jet Propulsion Laboratory, 1973 June 12. Many modifications were subsequently made. The Fortran 90 code was completed in April, 1995 by R. J. Hanson. The BVLS package is an additional item for the reprinting of the book by SIAM Publications and the distribution of the code package using netlib and Internet or network facilities. The GDL/IDL transcoding was completed by Philippe Prugniel on 2010 November 10th. A similar transcoding by M. Cappellari was used in earlier releases of the ULySS package. There are some differences in the arguments of both versions, and the present one was found slightly faster. See also the implementation of BVSL by R. Parker and P. Stark. It allows a "warm" start and may be more efficient than the present one.
(See pgm/bvls_lh.pro)
NAME: BVLS_PS PURPOSE: Bounded value least-square minimization (BVLS) USAGE: bvls_ps, a, b, bnd, x, istate, loopa, RNORM=rnorm ITMAX=itmax, EPS=eps STATUS=status, /QUIET, /DOUBLE, ARGUMENTS: A [input] On entry A contains the M by N matrix, A. B [input] On entry B contains the M-vector, B. BND [input] BND[0,J] is the lower bound for X[J]. BND[1,J] is the upper bound for X[J]. Require: BND[0,J] <= BND[1,J]. X [output] Solution N-vector. ISTATE [optional input/output] (N+1)-vector containing the state of each variable. If it is given in input, it indicates a "warm" start, i.e. the program will start from a guess of which variables are active (i.e. strictly within their bounds). In output, it contains the final state. This may be used when a series of related problems is treated (pass the vector computed at a preious call for the next). istate[n] is the total number of components at their bounds. The absolute values of the first nbound entries of istate are the indices of these `bound' components of x. The sign of istate[j] j=0...nbound-1, indicates whether x[|istate[j]|-1] is at its upper or lower bound. istate[j] is positive if the component is at its upper bound, negative if the component is at its lower bound. istate[j], j=nbound,...,n-1 contain the indices of the active components of x. KEYWORDS: [Optional input] LOOPA=loopa number of iterations taken in the main loop, Loop A. ITMAX=itmax, Maximum number of iterations permitted. Defaulted to 3*N, this is usually larger than required. EPS=eps, Determines the relative linear dependence of a column vector for a variable moved from its initial value. The default value is EPS=(machar()).eps (single precision). /DOUBLE, Setting this keyword forces some computations to be performed in double precision. This may be required if: (i) N is large (>10^6) or (ii) if the A or B values span a large amplitude (e.g. 10^20). This slows the computation, and is not required in most cases. /QUIET Verbosity control. By default the procedure prints the messages. KEYWORDS: [Optional output] RNORM=rnorm, Euclidean norm of the residual vector, || b - A*X ||. STATUS=status, Indicates status on return. 0 Solution completed. 1 A is not a 2D array 2 B is not a vector 3 Sizes of A and B and not compatible 4 BND is not a 2D array 5 Size of BND is not consistent with B 6 Some lower bounds are excessively large 7 Some upper bounds are excessively small 8 Inconsistent bounds 9 No free variable 10 Too many active variables in input STATE 11 Too many free variables 12 Failed to converge (reached ITMAX) DESCRIPTION: BVLS_PS solves linear least-squares problems with upper and lower bounds on the variables, using an active set strategy. It solves the problem: min || a.x - b || such that bl <= x <= bu 2 where x is an unknown n-vector a is a given m by n matrix b is a given m-vector bl is a given n-vector of lower bounds on the components of x. bu is a given n-vector of upper bounds on the components of x. The unconstrained least-squares problems for each candidate set of free variables are solved using the QR decomposition. BVLS has a "warm-start" feature permitting some of the variables to be initialized at their upper or lower bounds, which speeds the solution of a sequence of related problems. See the article ``Bounded Variable Least Squares: An Algorithm and Applications'' by P.B. Stark and R.L. Parker, in the journal: Computational Statistics, vol.10(2), (1995) for further description and applications to minimum l-1, l-2 and l-infinity fitting problems, as well as finding bounds on linear functionals subject to bounds on variables and fitting linear data within l-1, l-2 or l-infinity measures of misfit. AUTHOR: Robert L. Parker and Philip B. Stark (Version 3/19/90) Copyright of this software is reserved by the authors; however, this algorithm and subroutine may be used freely for non-commercial purposes, and may be distributed freely for non-commercial purposes. The authors do not warrant this software in any way: use it at your own risk. HISTORY: The original Fortran77 code, can be found in the authors web pages: http://www.stat.berkeley.edu/~stark/ It was converted to Fortran90 by Alan Miller (2000-11-23 Time: 23:41:42) This can be found from: http://jblevins.org/mirror/amiller/ The GDL/IDL transcoding of this latter was completed by Philippe Prugniel on 2011 January 18th. It use LA_LEAST_SQUARES (IDL) to perform the QR factorization. This routine is used in ULySS for replace the Lawson and Hanson version, BVLS_LH. It is about 3 times faster, and using the warm start capability provides another factor 3 gain when a minimization with several components is made.
(See pgm/bvls_ps.pro)
NAME: MREGRESS PURPOSE: Perform a multiple linear regression fit USAGE: Result = MREGRESS(X, Y, $ MEASURE_ERRORS=measure_errors, $ SIGMA=sigma, $ INVERSE=inv, $ STATUS=status) ARGUMENTS: X: (Input) The array of independent variable data. X must be dimensioned as an array of Npoints by Nterms, where there are Nterms coefficients (independent variables) to be found and Npoints of samples. Note that this array is transposed with respect to the input of the IDL routine REGRESS. Y: (Input) The vector of dependent variable points. Y must have Npoints elements. KEYWORDS: MEASURE_ERRORS : Set this keyword to a vector containing standard measurement errors for each point Y[i]. This vector must be the same length as X and Y. For Gaussian errors (e.g. instrumental uncertainties), MEASURE_ERRORS should be set to the standard deviations of each point in Y. For Poisson or statistical weighting MEASURE_ERRORS should be set to sqrt(Y). SIGMA : Set this keyword to a named variable to recieve the errors on the returned coefficients INVERSE : Set this keyword to a named value to receive the structure containing the covariance matrix and other terms that can be re-used to solve the system with the same X and other Y. STATUS : Set this keyword to a named variable to receive the status of the operation. Possible status values are: - 0 for successful completion, - 1 for a singular array (which indicates that the inversion is invalid), and - 2 which is a warning that a small pivot element was used and that significant accuracy was probably lost. If STATUS is not specified then any error messages will be output to the screen. RETURN: MREGRESS returns a column vector of coefficients that has Nterms elements. DESCRIPTION: Adapted from the IDL program REGRESS Perform a multiple linear regression fit: Y[i] = A[0]*X[0,i] + A[1]*X[1,i] + ... + A[Nterms-1]*X[Nterms-1,i] This is a variant of the IDL funtion REGRESS, where the covariance matrix and other intermediate arrays can be saved to solve faster the system with successive values of Y and the same X. An noticeable difference is also that MREGRESS does not include a constant term (if necessary, the constant may be included as one of the terms of X). HISTORY: Written, Ph. Prugniel 2008/01/17 (from the IDL REGRESS)
(See pgm/mregress.pro)
NAME: SLICES PURPOSE: Slice an image in horizontal/vertical direction USAGE: slices, image, xx, yy INPUTS: IMAGE: Input image, which should be already displayed XX: The x coordinates of the image YY: The y corrdinates of the image PROCEDURE: Plot interactively 2d slices of a 3d image in a seperate window. Use the left mouse button to switch between row and column profiles, while if the right mouse button is pressed the program exits. This routine is a modified version of the original IDL routine profiles.pro. The differences are: - The analysed image does not need to be scaled to the window size. - It takes as input the x/y coordinates of the image - It plots the profiles in DATA coordinates. EXAMPLE: Create an chi2 map: rangex = dindgen(10)*(10000.-3000.)/9.+3000 rangey = dindgen(10)/9.-1 map = uly_chimap(uly_root+'/data/VazMiles_z-0.40t07.94.fits', uly_ssp(),[[0,0],[0,1]], /QUIET, FILE_OUT='chimap.fits', XNODE=alog(rangex), YNODE=rangey) Plot the resulting chi2 map: uly_chimap_plot, map, /QUIET Invoke slices: slices, map.map, map.xnode, map.ynode HISTORY: This routine is a modification of the original IDL Library routine Profiles, M. Koleva, 2012/11/08
(See pgm/slices.pro)
NAME: TICKSLOG PURPOSE: Compute values for the ticks of a logarithmic axis USAGE: result = LOGLEVELS([range | MIN=min,MAX=max], [NMIN=nmin]) DESCRIPTION: Compute the position of the ticks for a graphic with a logarithmic scale (using /XLOG in PLOT). Rounded values, approximately logaritmically spaced, are generated. The minimum number of ticks positions to generate may be defined (by default it is 3), and the maximum number is about twice the minimum. INPUTS: RANGE -> A 2-element vector with the minimum and maximum value to be returned. Only levels _within_ this range will be returned. If RANGE contains only one element, this is interpreted as MAX and MIN will be assumed as 3 decades smaller. RANGE superseeds the MIN and MAX keywords. Note that RANGE must be positive definite but can be given in descending order in which case the labels will be reversed. KEYWORD PARAMETERS: MIN, MAX -> alternative way of specifying a RANGE. If only one keyword is given, the other one is computed as 3 decades smaller/larger than the given parameter. RANGE superseeds MIN and MAX. NMIN -> the minimum number of ticks to gnerate OUTPUTS: A vector with "round" logarithmic values within the given range. EXAMPLE: data = 5 + findgen(101)/100. *995 range = [min(data), max(data) ] ticks = tickslog(range) nticks = n_elements(ticks) plot, data, /YLOG, YRANGE=range, YTICKS=nticks-1, YTICKV=ticks HISTORY: 2008/08/29, Philippe Prugniel, created from LOGLEVELS by Martin Schultz mgs, 17 Mar 1999: VERSION 1.00
(See pgm/tickslog.pro)
NAME: ULY_HTML_HELP PURPOSE: Generate HTML documentation from the embedded doc USAGE: ULY_HTML_HELP, <sources>, <outfile> [, VERBOSE=<verbose>][, TITLE=<title>] [, HEADER=<header>] [, STYLE_SHEET=<style_sheet>] [, CROSSREF=<crossref>] [, XREF_FILE=<xref_file>] ARGUMENTS: SOURCES: A string or string array containing the name(s) of the .pro files (or the names of directories containing such files) for which help is desired. If a source file is an IDL procedure, it must include the .PRO file extension. All other source files are assumed to be directories. OUTFILE: The name of the output file which will be generated. KEYWORDS: VERBOSE: Normally, ULY_HTML_HELP does its work silently. Setting this keyword to a non-zero value causes the procedure to issue informational messages that indicate what it is currently doing. !QUIET must be 0 for these messages to appear. TITLE: If present, a string which supplies the name that should appear as the Document and page Title. By default, the title is "Extended IDL help". HEADER: If present, this string or string-array will be included at the beginning of the body element, before the title heading. STYLE_SHEET: If present, a string giving the name of an external style sheet (CSS) that will be linked from the HEAD element of the generated HTML document. CROSSREF: If present, the program will only find the cross-references and append them to the file specified in XREF_FILE XREF_FILE: If given in combination with CROSSREF, name of the output file where the cross-reference are appended. Otherwise if CROSSREF is not set, name of the file used to read the cross-references. Cross-references will make the name of other documentation blocks clickable. To use the documentation with cross-references, the file 'uly.js' must be copied on the same directory as the html files in order to have the purpose of a function displayed when the cursor move over its name (this is a javascript code) KEEPPRO refs to routine_name.pro not to routine_name.html DESCRIPTION: This routine is a modified version of the original IDL function MK_HTML_HELP: The differences are: - Make reference to a style-sheet (CSS) that can be passed through the keyword STYLE_SHEET. Define the following classes for some components of the created page: list : list of routines routine : one DIV for each routine documentation In addition, some elements have particular IDs: page_title : The main title (H1 header) doc_tag : For each tag in the documentation - The list of routines is presented as a table with their name and purpose. - Generate cross references through the different documented modules. - Generate HTMLized files with the source of the program, using the function PRO2HTML. Note that PURPOSE is meant to be a short (one line) description of the function. The detailed description should come in the tag DECRIPTION. See another variant of the same original program: http://astro.berkeley.edu/~marc/idlshare/general/html/ HISTORY: This routine is a modification of the original IDL Library routine MK_HTML_HELP, Ph. Prugniel, 2008/08
(See pgm/uly_html_help.pro)
NAME: ULY_LINE PURPOSE: Define a spectral-line component USAGE: cmp = uly_line(wave, sigma, LSF=lsf_file, NAME=name, LL=lim_wave, SL=lim_sigma, WL=lim_weight) ARGUMENT: WAVE: Input. [Angstroem] Wavelength of the line SIGMA: Input. [km/s] Broadening of the line KEYWORDS: LSF: Name of the file containing a relative LSF to be injected in the template. NAME: A string, name of this component. By default a unique name is generated. LL: Permitted limits for the fitted wavelength (in Angstroem). If given, it must be an array of two values. SL: Permitted limits for the velocity dispersion, sigma (in km/s). If given, it must be an array of two values. WL: Permitted limits for the weight of each model component [min_w, max_w]; given not normalized. If a single value is passed, take it for 'min_w'. Note that this constraints works only in the case of multiple components. DESCRIPTION: Return a component to fit a single line (emission or absorption). The model is either an unresolved line (Dirac, if SIGMA is not given) or a Gaussian of standard deviation SIGMA, in km/s. (note that the line may be convolved with a LOSVD in the fitting process if ULYSS is called with KMOMENT > 1). The free parameters are the wavelength of the line, its gaussian broadening and the weight of the component, determining the height or depth (intensity) of the line. In the future we may modify this procedure to make the Gauss-Hermite coefficients (h3 h4) of the line other free parameters. This cmp can be passed to ULYSS. The functions ULY_LINE_INIT and ULY_LINE_EVAL attached in this file are automatically called by the fitting routines. OUTPUT: cmp structure or array of structures, for more detailed explanation check ULY_FIT. EXAMPLE: The following example shows how to fit the H & K lines in a spectrum. Note that using an additive term (AD=0) is required to fit the continuum. spect = uly_spect_read(uly_root+'/data/m67.fits') cmp1 = uly_line(3969e, 400e) cmp2 = uly_line(3934e, 400e) cmp = [cmp1, cmp2] ulyss, spect, cmp, AD=0, /PLO SEE ALSO: Definition of other components: ULY_SSP, ULY_TGM, ULY_STAR REQUIRED FUNCTION: ULY_CMP_NAME ULY_SPECT_ALLOC ULY_SPECT_FREE ULY_SPECT_LSFCONVOL ULY_SPECT_LOSVDCONVOL HISTORY: Mina Koleva, 2008/05 created
(See pgm/uly_line.pro)
NAME: ULY_STAR PURPOSE: Define a star component USAGE: cmp = uly_star(file_in, LSF=lsf_file, NAME=name, $ WL=lim_weight) ARGUMENT: FILE_IN: Input spectrum file name, spectrum structure (see ULY_SPECT_ALLOC) or array of file names or spectrum structures. In the package we supply a solar spectrum as template, the users can find and use their own template star(s). KEYWORDS: LSF: Name of the file containing a relative LSF to be injected in the template. NAME: A string, name of this component. By default a unique name is generated. WL: Limits for the weight of each model component [min_w, max_w]. The weight is in data units over component units. By default the weight is constrained to be positive, to suppress any constraint, set: WL=[-1,1]*machar(/DOUBLE)).xmax This constraint is ignored when ULySS fits a single component. DESCRIPTION: Return a star fitting component or array of star component to be passed to ULYSS. A star component has no free parameter. It can be a stellar template or any sort of spectrum. It is typically used to fit the LOSVD. Using an array of such components allows to search for the best positive linear combination of spectra matchin an observation: this approach, called 'optimal template matching', is commonly used to analyse the kinematics of galaxies. Note that a TGM component (see ULY_TGM) can be used to determine in the same time the LOSVD and the atmospheric parameters that would best represent the fitted spectrum. OUTPUT: Single star cmp struct, for more detailed explanation please check ULY_FIT. REQUIRED FUNCTION: ULY_CMP_NAME EXAMPLE: cmp = uly_star(uly_root+'/models/sun.fits') HISTORY: Philippe Prugniel, 2008/05 created
(See pgm/uly_star.pro)
NAME: ULY_STAR_INIT PURPOSE: Initialize a single star fit component USAGE: cmp = uly_star_init(cmp, WAVERANGE=waverange, VELSCALE=velscale, /QUIET) ARGUMENTS: CMP: A component defined using ULY_STAR, please check it for more details. KEYWORDS: WAVERANGE: Wavelength range that the ULYSS fitting procedure will use. VELSCALE: Specifies the sampling step to be used by ULYSS, in km/s /QUIET: Verbosity control: suppress the messages SEE ALSO: Definition of other components: ULY_SSP, ULY_TGM, ULY_LINE HISTORY: Philippe Prugniel, 2008/05 created
(See pgm/uly_star.pro)
NAME: UNDEFINE PURPOSE: Destroy a variable USAGE: undefine, varname ARGUMENT: varname : Named variable to delete DESCRIPTION: Deletes a variable so that after n_elements(varname) will be 0. Note that if the variable is a pointer or structure containing pointers, this routine does not free the memory. HISTORY: Philippe Prugniel, 2008, created
(See pgm/undefine.pro)
NAME: XREBIN PURPOSE: Rebin a vector USAGE: yout = xrebin(xin, yin, xout [,/SPLINE] [,/LSQUADRATIC] [,/QUADRATIC], [/SPLINF], [/CUBIC], [/NAN] ARGUMENTS: XIN (N+1) points of the borders of each input bin YIN values (N) of input bins (1 or 2D) XOUT points of the borders of output bins KEYWORDS: SPLINE Spline interpolation (default: linear) LSQUADRATIC Least square quadratic interpolation QUADRATIC Quadratic interpolation SPLINF Other version of spline interpolation 3-4 times faster than SPLINE, the results are slightly different, but seemingly better. CUBIC Use the function INTERPOLATE with CUBIC=-05 Cubic convolution is an interpolation method that closely approximates the theoretically optimum sinc interpolation function using cubic polynomials. NAN To treat the missing values in YIN as zeros when interpolating it. Use with caution. OUTPUT YOUT output 1 or 2D array, rebined at the borders requested by XOUT. YOUT is in double precision (it shall be converted in float if appropriate. DESCRIPTION: This function rebins a vector or a n-dimensional array into different bins along the first axis. The rebining is made as integrating the signal, interpolating and differentiating the result. The integrated signal is preserved. Different interpolating algorithms are proposed, and only one of the keywords SPLINE LSQUADRATIC QUADRATIC SPLINF CUBIC must be given. If not of them is specified, the routine INTERPOL will make a linear interpolation. This latter routin also handle the SPLINE LSQUADRATIC QUADRATIC cases. The difference between SPLINF and SPLINE is that the spline representation is computed on the whole vector instead of on 4-points blocks (this is faster and often our preferred option). CUBIC is using the function INTERPOLATE with CUBIC=-0.5 (sinc). No check is performed to prevent extrapolation (beware to give proper input). The computations are systematically made in double precision. The reason is that the algorithm results in loss of relative precision equal to the number of pixels. In float, this is often a unacceptable. AUTHOR: Mina Koleva (2007/03/01)
(See pgm/xrebin.pro)
uly_html_help
.Contact: ulyss at obs.univ-lyon1.fr |
Last modified: Wed Feb 24 10:42:37 2016. |