This code behind this document is an application of the Hilbert Class Library. So the first thing you need to do is to make HCL available. If HCL is already installed on your system, you do this by adding a line to your .cshrc file or appropriate component file. Otherwise you must install HCL first. The easiest way to do this is to download it from the TRIP web page:
http://www.trip.caam.rice.eduand follow the installation instructions.
The code also depends on the SU/SEGY vector class package sVector, which therefore must also be installed. It will be part of the next HCL release.
NB: At Rice/TRIP and Stanford/SEP, no installation is necessary: the packages are already installed. Simply add the following lines to your shell environment files:
setenv HCLROOT /import/masc39c/symes/hclr0.9 setenv KBDAROOT /import/masc39c/symes/kbda setenv QPROOT <path to the root directory of this package>
setenv HCLROOT /jon/symes/hclr0.9 setenv KBDAROOT /jon/symes/kbda setenv QPROOT <path to the root directory of this package>
The ``root directory of this package'' referenced in these instructions is the directory you create by downloading the tar file containing this report. In so doing, you create a directory tree with root named qcqm. This is the ``root'' in question. All pathnames in the following discussion are relative to this root.
HCL includes a set of make rules which evolved from the SEP rule set as it was about two years ago. I am sure that SEP's rules have also evolved, and differently. The makefiles for this report are all output of the HCL makefile autowriter, maw, and use HCL's rules. If you get as far as modifying makefiles by hand, bear the possible HCL/SEP incompatibility in mind.
You can rebuild the entire package simply by entering make build in the root directory. You will construct all of the executables and final results (.ps files in the Fig) directories, including this postscript version of this report. You can make individual results in the usual way.
Note that all figures in this report except the first two will vary from build to build, as you choose a new pseudorandom seed each time you execute the commands.
The command for deconvolution is sfilter/decon.x. You execute it by following it with a parameter file name: decon.x par - it reads all of its parameters from the file. Regrettably the ``getpar'' device used in curren HCL programs is not flexible enough to permit specification of parameters on the command line. Probably I should just steal SEP's getpar!
You can use the executable sfilter/decon.x with other data by altering its input parameters, and so explore the capabilities of the algorithm using data other than that supplied with this report. Parameter control requires you to edit the parameter files manually. HCL parameter files are keyword=value lists; the values can be integers, floating point numbers with any size of mantissa, and strings. Parameters to be read only by one part of the program (typically a class constructor) get an identifying string prepended, with a double colon. Thus the parameter Tol for the conjugate gradient algorithm becomes CG::Tol. That is, the parameter file can specify many variables named Tol, so ong as they have been equipped with qualifiers which allow each program unit to choose a unique value.
Files should be in either SU (stripped SEGY) or SEGY formats. sVector currently supports only native binary floating point representation, so if you port data from Linux to SGI etc. you will have to byte-swap it.
Here is the parameter file structure for the deconvolution example:
Sigma=0.5 Lambda=0.0001 Wavelet="rick15.su" DataTimeSeries="fnd.su" CG::Tol=1.e-4 CG::MaxItn=50 CG::DispFlag=3 QCQM::Tol=0.01 QCQM::MaxItn=10 QCQM::DispFlag=1
The parameters are
http://www.trip.caam.rice.eduand follow the links to the HCL reference manual, page on the conjugate gradient algorithms). Level 3 is max output, including a summary of the progress at each iteration. Note that the residual reported is the normal residual for this application!
Amongst the diagnostics printed out at the end of a run when QCQM::DispFlag is set you will find ``Lagrange cosine''. This is the cosine of the angle between the constraint gradient and the objective gradient. It diagnoses the success of the constrained optimization: if it is very close to 1, then the two gradients are parallel and the first order necessary condition has been satisfied. This occurs in the deconvolution examples in all cases except that depicted in Figure 8, in which nothing converges and you can't fit the data. Apparently failure to get the Lagrange cosine close to 1 in a reasonable number of CG and Moré-Hebden iterations implies that the noise level has been set too small and you are trying to match data components associated with very small singular values.