Overview

In this section you will find documentation how to obtain and use the commercial program matlab resp. its open source (clone) version octave in combination with McArtim result files to create inversion algorithms. The code examples will be in the octave syntax, since this is the authors preferred software.

How to obtain octave or matlab

matlab is a commercial program that can be obtained from mathworks.de. There are both windows an linux versions. If you don't want to pay the licence fee or if you support the open source movement you can alternatively use the open source (clone) version octave or any other software with the same functionality. The syntax of matlab and octave is very similar but there are some minor differences. A short introduction to matlab is kindly provided by Iske and Behrens from the TU München. octave can be obtained for linux and windows platforms from the octave homepage. Documentation can be found there too or at

Loading McArtim results or other text files

Usually the output of McArtimQuestion of McArtimScript is in table form. Such a table can be loaded by the octave command

Result=load('tablefile.txt');
The content of the table file is therewith stored into the matrix Result. Please be aware that one maybe has to remove all comment lines in the McArtim output file. Assume that a single line of a question result is of the form:
uuuxexexe....xe,
see for example the output of absorber_derivative.conf). Single characters in the string correspond to TAB-separated numbers. The 'u's are unteresting values, 'x' are the quantities we are interested in and 'e' is the respective error (Monte Carlo noise). Thus, in the example above we are interested in the columns 4,6,8,... and separately 5,7,9,... containing the errors.
The following octave script loadMcArtimResultMatrix.octave loads a McArtim output file filename and extracts two matrices M and Me with each n columns beginning with the column first_column: the file loadMcArtimResultMatrix.octave must be located in the same folder as the calling script. The calling script line could look like:
[M Me]=loadMcArtimResultMatrix('c:\McArtim\BoxSensitivities.txt', 4, 10);

Loading grid files and linearely interpolated profile files

A grid file (see section on grids) is simply a list of height values in rows. The command

h=load('grid/height_grid.txt');
thus stores this list into the column vector h (in fact this is a N x 1 matrix). This grid vector h can then be used to create (script: loadHeightGridVectors.octave) a vector hm containing layer middles:

With the layer middles one can write a octave script loadLIProfile.octave to obtain a profile vector by evaluating a linearely interpolated profile at the layer middle heights:

Building essential quantities

For optimal estimation or other inversion schemes (see thesis inversion chapter), some essential vectors and matrices (auch as the state and measurement vector or the covariance matrices) have to be built. In octave (script file: buildDiagonalCovarianceMatrix.octave) this can be done in very few lines:

Some people also use non diagonal a priori covariance matrix elements, which can be built (script file buildDampingStateCovarianceMatrix.octave) as follows:

The next example script buildRegularisationMatrix.octave shows how to build a regularisation matrix used for a first derivative constrained inversion:

DSCDs and DBoxAMFs

Many DOAS evaluations produce differential SCDs (DSCDs). The components of the measurement vector are then to be understood as

yi=log(Ir/Ii)/σ=log(Ir)/σ-log(Ii)/σ
where σ is the cross section of the absorber of interest. In order to use McArtim absorber derivatives (familiar with BoxAMFs), one needs differences of logarithmic derivatives. The following script buildDifferentialJacobiMatrix.octave achieves this: Practically the Jacobians for the reference geometry have to be in the first row of the output file. Using McArtimScript, the reference geometry thus must be the first to be processed.

A helper function: write matrices to files

(script file: writeMatrixToFile.octave)

Linear inversion example

The "vanilla" linear inversion optionally with a priori and/or regularisation (due to insufficient information content, too few geometries) is as follows.

cost function

χ2=[F(x)-y]TSε-1[F(x)-y] +α[x-xa]TSa-1[x-xa] +βxTRx
The regularisation matrix can be obtained with the script buildRegularisationMatrix.octave. α should be set either to 0 or 1. This controls whether the a priori profile is considered or not. β is a parameter for the regularisation that can be found by trial and error (depening on the accuracy of the retrieved profile). The feature can be skipped by setting β=0. The forward model is linearized by Taylor expansion:
F(x)≈F(xa)+K(x-xa)

Jacobian or gradient

The cost function is to be minimised so we need the gradient:

½∇χ2 =KTSε-1[F(xa)+K(x-xa)-y] +αSa-1[x-xa] +βRx

Solution

Setting the gradient to zero yields:

x=[ KTSε-1K +αSa-1 +βR ]-1[ KTSε-1 [y-F(xa)+Kxa] +αSa-1xa ]

here is the script (linear_inversion_vanilla.octave):