Chapter 4. How to use GRASP

4.1. How to run the code

After the installation of GRASP software, the grasp command will be available. As it is explained in a previous section (see Section 3.6, “Running the code”) just typing grasp will print some general information about how the software was compiled.

First arguments that grasp executable expects is a path to the settings file. The settings file describes the inversion strategy and general behaviour of the process: where is the input data, in which format is the input data, where to store the output results... Therefore, a deep knowledge of the settings parameters is the base to understand GRASP.

4.1.1. Settings file

Settings file is written in YAML format and that brings many benefits: easy to write, clear to read, self-explanatory names, flexible and powerful. The concepts are organized in blocks that are translated to YAML thanks to the fixed indentation (we suggest 4 white-spaces). So, for example, the first level defines the different modules:

                    

input:
    # Here settings related with input module
    segment: # It defines a description of the input segment
        x: 2 # It defines the maximum size of x dimension of the segment to 2.
output:
    # Settings linked with the output
    # ...
retrieval:
    # Definition of inversion strategy
    # ...
...
            

When GRASP is called, the first step is to read settings and the second is to prepare the environment for the settings defined in the main structures. If the settings are not valid, an explicit error message will be printed. Please read the first line of it carefully to understand the error.

The settings parameters can also be defined by the command line. In this case, after the first argument (settings file name), extra settings parameters can be defined with the syntax key=value where key is the parameter name in "dot syntax". For example, in GRASP, it is equivalent to be called with the argument input.driver=sdata or to have defined in the settings file the following content:

input:
    driver: sdata
                

It is important to define clearly how the relative paths have to be defined in the settings files. Relative paths are always relative to the file that defines it. In the next section, the reader will learn about how to include other settings file inside of a settings file, but this rule will stay valid: Relative paths are defined from settings file that define it. In the case of usage of command line, relative paths are relative to the current working directory. In the case of absolute paths, all this complexity disappears but the results are less portables.

4.1.1.1. HELP argument

The list of available parameters for GRASP is long. The "help" argument will help you to know the available parameters or to look for something specific. When help argument is present, GRASP is not executed normally, but instead, the help information will appear in the screen print. Additionally, help can be followed by a search string to filter the results. For example, help=input will print only the settings which contain "input" string in its definition and help=input.segment will print only the settings which will be under the block "segment" inside of the block "input" (because "dot" symbol defines block separator).

Table 4.1. List of GRASP options

Field NameField Content
help Shows this help information
version Shows GRASP code version information and stops
resources_path Path to framework resources folder
retrieval. general. path_to_internal_files Path to internal data files
retrieval. mode Define how the code is run. Valid values 'inversion' or 'forward' (inversion=full inversion; forward=only forward model)
retrieval. inversion. convergence. minimization_convention Minimization in absolute or logarithm space
retrieval. inversion. convergence. threshold_for_stopping_Q_iterations Threshold for stopping q - iterations: once the change in the residual smaller than this parameter the iterations are stopped
retrieval. inversion. convergence. scale_for_finite_difference Defines DL for calculations of derivatives: (f(x+DL)-f(x))/DL
retrieval. inversion. convergence. threshold_for_stopping Threshold for stopping p - iterations: once the change in the residual smaller than this parameter the iterations stops
retrieval. inversion. convergence. normal_system_solver Defines the method to solve normal system
retrieval. inversion. convergence. maximum_iterations_of_Levenberg-Marquardt The maximum number of ip iterations where Levenberg-Marquardt correction is applied
retrieval. inversion. convergence. maximum_iterations_for_stopping Maximum number of iterations performed in the retrieval before stop
retrieval. inversion. convergence. shift_for_applying_logarithm_to_negative_values This value (usually 1) is added to negative parameters in order to be able to apply logarithmic transformations
retrieval. inversion. regime_of_multipixel_constraints. inversion_regime Flag for single multi-pixel retrieval (VALUES)
retrieval. inversion. regime_of_measurement_fitting. polarization 1. -We fit: I, Q, U; 2. - We fit: I, Q/I, U/I; 3. - We fit: I, (Q^2+U^2)^(1/2) (Polarized measurements can be defined in input data as Q and U or P); 4- We fit: I, [(Q^2+U^2)^(1/2)]/I (Polarized measurements can be defined in input data as Q and U or P); 5- We fit: P/I (User has to provide P/I in input data)
retrieval. inversion. regime_of_measurement_fitting. scattering_angle_for_normalized_p11 For P11 retrieval. If P11 is given as absolute value, this option has to be -180. 0 (default value). If P11 is relative, indicate with this option the value of the angle that has to be used to divide all p11 values.
retrieval. inversion. noises. noise[]. standard_deviation_synthetic Standard deviation of synthetic random noise added to the corresponding inverted data, if value is 0 no synthetic random noise is added
retrieval. inversion. noises. noise[]. error_type Error type used for definition of covariance matrices used in measurement fitting and in modeling synthetic noise
retrieval. inversion. noises. noise[]. standard_deviation Standard deviation of the random noise expected (this defines the covariance matrix used in fitting)
retrieval. inversion. noises. noise[]. measurement_type[]. type Type of relevant measurements for applying the noise assumptions
retrieval. inversion. noises. noise[]. measurement_type[]. index_of_wavelength_involved List of indices of relevant wavelengths for applying the noise assumptions
retrieval. forward_model. phase_matrix. size_binning_method_for_triangle_bins Determining the scale used for the binning of size distribution
retrieval. forward_model. phase_matrix. number_of_elements Number of phase matrix elements used in the calculations and retrieval
retrieval. forward_model. phase_matrix. kernels_folder Path to kernels when we retrieve size distribution for triangle bins or lognormal size distribution
retrieval. forward_model. phase_matrix. radius. mode[]. bins Number of triangle bins set for size distribution representation
retrieval. forward_model. phase_matrix. radius. mode[]. min Minimum value of the radius used in the triangle bins
retrieval. forward_model. phase_matrix. radius. mode[]. max Maximum value of the radius used in the triangle bins
retrieval. forward_model. phase_matrix. ratio. mode[]. value Values of bins using in retrieved aspect or axis distributions
retrieval. forward_model. radiative_transfer. number_of_layers Maximum number of vertical layers resolved in radiative transfer (minimum and maximum value. If only one number is assigned maximum is by default KNT-1)
retrieval. forward_model. radiative_transfer. molecular_profile_vertical_type It defines which model will be used to describe vertical profile of molecular (Rayleigh) scattering. Values: 'Exponential' describes molecular density profile at altitude h as exp(h/8)/8, 'Standard_atmosphere' uses standard atmosphere model (pressure and temperature profiles) to calculate molecular density at each altitude.
retrieval. forward_model. radiative_transfer. aerosol_profile_vertical_type It defines which model will be used to describe vertical profile of aerosol distribution. Values: 'Exponential' describes aerosol concentration profile at altitude h as exp(h/HM)/HM, 'gaussian' uses normal distribution exp(h-HM)/(sqrt(rpi)*sqrt(2. )*sigma) with sigma and HM parameter taken from parameters.
retrieval. forward_model. radiative_transfer. phase_matrix_truncation Switch on/off the truncation: the technique to calculate scattering effects from the sharp forward peak of the phase function separately from those of the rest of the phase function. It doesn't effect considerably the accuracy while provides much faster calculations
retrieval. forward_model. radiative_transfer. absolute_error_rt_calculations Absolute value of truncation threshold of Fourier and order-of-scattering series expansions in radiative transfer calculations
retrieval. forward_model. radiative_transfer. reference_plane_for_polarization Reference plane for polarization calculations
retrieval. forward_model. radiative_transfer. BOA_reflectance Perform atmospheric correction calculation of surface reflectance and surface BRDF at the Bottom-Of-Atmosphere
retrieval. forward_model. radiative_transfer. rt_kernels. mode LUT Mode for radiative transfer (disable by default)
retrieval. forward_model. radiative_transfer. rt_kernels. folder Folder for aerosol look-up-table
retrieval. forward_model. radiative_transfer. simulating_observation. order_of_scattering Regime of scattering used for modeling diffuse radiation observations
retrieval. forward_model. radiative_transfer. simulating_observation. number_of_gaussian_quadratures_for_expansion_coefficients Number of Gaussian quadratures for calculating expansion coefficients used for multiple scattering simulations
retrieval. forward_model. radiative_transfer. simulating_observation. number_of_guassian_quadratures_for_fourier_expansion_coefficients Number of Gaussian quadratures for calculating Fourier expansion coefficients used for multiple scattering simulations
retrieval. forward_model. radiative_transfer. simulating_observation. number_of_fourier_expansion_coefficients Number of Fourier expansion coefficients used in multiple scattering simulations
retrieval. forward_model. radiative_transfer. simulating_derivatives. order_of_scattering Regime of scattering used in calculation of radiance derivatives
retrieval. forward_model. radiative_transfer. simulating_derivatives. number_of_gaussian_quadratures_for_expansion_coefficients Number of Gaussian quadratures for calculating expansion coefficients used in calculation of radiance derivatives
retrieval. forward_model. radiative_transfer. simulating_derivatives. number_of_guassian_quadratures_for_fourier_expansion_coefficients Number of Gaussian quadratures used for calculating Fourier expansion coefficients used in calculation of radiance derivatives
retrieval. forward_model. radiative_transfer. simulating_derivatives. number_of_fourier_expansion_coefficients Number of Fourier expansion coefficients used in calculation of radiance derivatives
retrieval. product_configuration. wavelength_indices_for_angstrom Indices of wavelengths which will be used to calculate Angstrom exponent
retrieval. product_configuration. aerosol_particulate_matter_diameter Diameters of aerosol particles in microns which will be used to calculate Particulate Matter (PM)
retrieval. product_configuration. wavelength_indices_for_ndvi Indices of wavelengths which will be used to calculate NDVI if it is calculated
retrieval. product_configuration. wavelenght_indices_for_aod_error_estimation Indices of wavelengths which will be used to estimate aod error
retrieval. product_configuration. wavelenght_indices_for_ssa_error_estimation Indices of wavelengths which will be used to estimate ssa error
retrieval. product_configuration. wavelenght_indices_for_lidar_error_estimation Indices of wavelengths which will be used to estimate lidar error
retrieval. products. aerosol. chemistry Retrieve aerosol chemical composition (if retrieved)
retrieval. products. aerosol. lidar Retrieve columnar lidar ratios (e. g. , if lidar data are inverted)
retrieval. products. aerosol. optical_properties Provide aerosol optical properties
retrieval. products. aerosol. phase_matrix Obtain aerosol phase matrix
retrieval. products. aerosol. refractive_index Provide aerosol refractive index
retrieval. products. aerosol. theoretical_bimodal_extinction Provide estimated aerosol extinction for fine and coarse modes
retrieval. products. aerosol. theoretical_bimodal_parameters Provide estimated aerosol microphysical characteristics for fine and coarse modes
retrieval. products. aerosol. particulate_matter Obtain aerosol particulate matter estimation at given particle diameters
retrieval. products. aerosol. type Obtain aerosol type index i. e. 0-Complex mixture,1-Background,2-Maritime,3-Urbn. Polluted,4-Mixed,5-Urbn. clean,6-Smoke flam. ,7-Smoke. sold. ,8-Dust
retrieval. products. error_estimation. aerosol. lidar Implement error estimation for aerosol lidar products
retrieval. products. error_estimation. aerosol. optical_properties Implement error estimation for optical properties of aerosol products
retrieval. products. error_estimation. parameters Implement error estimation for retrieved parameters
retrieval. products. forcing. broadband_flux Provide upward and downward fluxes integrated over the solar spectrum at the user-defined levels for the scenarios with and without aerosols
retrieval. products. forcing. forcing Provide aerosol radiative forcing values at the user-defined levels
retrieval. products. retrieval. fitting Provide obtained measurement fitting
retrieval. products. retrieval. parameters Provide retrieved parameters
retrieval. products. retrieval. residual Provide values of obtained residuals
retrieval. products. surface Provide surface reflectance products
retrieval. edges_size. x Size of edges width in pixels (it have to be lower than KIEDGE compilation constant)
retrieval. edges_size. y Size of edges height in pixels (it have to be lower than KIEDGE compilation constant)
retrieval. edges_size. t Size of temporal dimension of edges (it have to be lower than KIEDGE compilation constant)
retrieval. debug. verbose Retrieval prints progress information while it is performing the inversion
retrieval. debug. additional_information Print some additional information
retrieval. debug. simulated_sdata_file Filename where simulated observation data are to be printed. This option force retrieval. general. stop_before_performing_retrieval=true
retrieval. debug. path_to_extra_files Path to folder with retrieval extra files: debug information, extra resources like image. dat files. . .
retrieval. debug. use_internal_initial_guess Test option which allows to load different initial guess for each pixel (loading them from image. dat files)
retrieval. constraints. characteristic[]. type Type of characteristic
retrieval. constraints. characteristic[]. retrieved Specify if this characteristic will be retrieved or only used in forward model
retrieval. constraints. characteristic[]. mode[]. initial_guess. value Initial values for a specific (determined by type) characteristic
retrieval. constraints. characteristic[]. mode[]. initial_guess. min Minimum value for the specific characteristic
retrieval. constraints. characteristic[]. mode[]. initial_guess. max Maximum value for the specific characteristic
retrieval. constraints. characteristic[]. mode[]. initial_guess. index_of_wavelength_involved Indices of Wavelengths associated to the specific characteristic
retrieval. constraints. characteristic[]. mode[]. initial_guess. estimate_error Flag to retrieve the error of specific parameter if retrieval. products. error_estimation. parameters is true
retrieval. constraints. characteristic[]. mode[]. single_pixel. a_priori_estimates. lagrange_multiplier Value of the Lagrange multiplier associated to a priori estimate of the retrieved characteristics (applied in each single pixel)
retrieval. constraints. characteristic[]. mode[]. single_pixel. smoothness_constraints. difference_order Order of the derivatives/differences used for applying a priori smoothness constrains for the retrieved characteristics
retrieval. constraints. characteristic[]. mode[]. single_pixel. smoothness_constraints. lagrange_multiplier Value of the Lagrange multiplayer used for applying a priori smoothness constraints for the retrieved characteristics
retrieval. constraints. characteristic[]. mode[]. multi_pixel. smoothness_constraints. derivative_order_of_X_variability Order of derivatives/differences used for applying a priori smoothness constraints on parameter inter-pixel X-variability in multi-pixel retrieval regime
retrieval. constraints. characteristic[]. mode[]. multi_pixel. smoothness_constraints. lagrange_multiplier_of_X_variability Value of the Lagrange multiplier used for applying a priori smoothness constraints on parameter inter-pixel X-variability in multi-pixel retrieval regime
retrieval. constraints. characteristic[]. mode[]. multi_pixel. smoothness_constraints. derivative_order_of_Y_variability Order of derivatives/differences used for applying a priori smoothness constraints on parameter inter-pixel Y-variability in multi-pixel retrieval regime
retrieval. constraints. characteristic[]. mode[]. multi_pixel. smoothness_constraints. lagrange_multiplier_of_Y_variability Value of the Lagrange multiplier used for applying a priori smoothness constraints on parameter inter-pixel Y-variability in multi-pixel retrieval regime
retrieval. constraints. characteristic[]. mode[]. multi_pixel. smoothness_constraints. derivative_order_of_T_variability Order of derivatives/differences used for applying a priori smoothness constraints on parameter inter-pixel T-variability in multi-pixel retrieval regime
retrieval. constraints. characteristic[]. mode[]. multi_pixel. smoothness_constraints. lagrange_multiplier_of_T_variability Value of the Lagrange multiplier used for applying a priori smoothness constraints on parameter inter-pixel T-variability in multi-pixel retrieval regime
settings. debug Shows settings read to run this program
settings. strict Force GRASP to continue when there were parse or validation errors in the settings file
settings. dump Stream to dump read settings in short format(experimental)
settings. long_dump Stream to dump read settings in long format (experimental)
input. driver The driver that will be called for inverting data
input. file Name of file(s) which contain input observation
input. center. latitude Latitude of the center of tile to invert
input. center. longitude Latitude of the center of tile to invert
input. corner. row The number of the row for input driver in the native coordinate system of the sensor
input. corner. column The number of the column for input driver in the native coordinate system of the sensor
input. grid_offset. row Information of the first row in the input grid that will be use for normalizing the output. If 0 is used the output coordinate reference will be the same than the input. This offset is for forcing the output row to start at 0 (e. g. if you put 1, output_row == input_row - 1)
input. grid_offset. column Information of the first column in the input grid that will be used for normalizing the output. If 0 is used the output coordinate reference will be the same than the input. This offset is for forcing the output column to start at 0 (e. g. if you put 1, output_column == input_column - 1)
input. area. width The width of the covered area in pixels. It has to be divisible by intput. segment. x value
input. area. height The height of the covered area in pixels. It has to be divisible by intput. segment. y value
input. time. from Initial date and time for data processing
input. time. to Final date and time for data processing
input. segment. x Size of segment width in pixels (it have to be lower than KIX compilation constant)
input. segment. y Size of segment height in pixels (it have to be lower than KIY compilation constant)
input. segment. t Size of segment temporal dimension (it have to be lower than KITIME compilation constant)
input. transformer Name of input data transformer functions to be used after load data
input. debug. raw_segment Stream to print raw segment data loaded
input. debug. clean_segment Stream to print segment information after clean NaN values
input. debug. used_files Stream to print names of the files that have the pixels for inverting
input. sdata. dump Stream where to dump sdata information
input. sdata. dump_original Stream where to dump sdata information just after being generated by the driver. Some transformers can modify sdata information, this setting is thought for debugging purposes where the user is interested in knowing sdata generated by the driver instead of data driving inside the inversion.
input. imagedat. dump Stream where to dump initial guess information (image. dat format)
input. preload_segment. x This parameter specifies how many segments in X dimension will be preloaded in each block
input. preload_segment. y This parameter specifies how many segments in Y dimension will be preloaded in each block
input. preload_segment. t This parameter specifies how many segments in T dimension will be preloaded in each block
input. driver_settings. sdata. debug Print debug information from sdata reader subsystem
input. transformer_settings. segment_imagedat. file File which contains initial guess information in classic input. dat format
output. segment. function Driver to process output for every single retrieval (show information in screen, perform a map, plotting, . . . )
output. segment. stream Stream to dump segment output data
output. iteration. function Driver to process output for every single retrieval (show information in screen, perform a map, plotting, . . . )
output. iteration. stream Stream to dump segment output data
output. tile. function Driver to process output after processing complete tile (show information in screen, perform a map, plotting, . . . )
output. tile. stream Stream to dump tile output data
output. current. function Driver to process output after each retrieval (show information in screen, perform a map, plotting, . . . )
output. current. stream Stream to dump current progress information about the retrieval conducted
output. sdata. simulated_file Stream where to dump sdata fitting information. Fitting product has to be enabled, otherwise this is ignored. Fitting is dumped after retrieval process.
output. segment. function_settings. csv. delimiter Separator between fields
output. segment. function_settings. csv. compression If true output is compressed in GZ format (gz extension is automatically added)
output. segment. function_settings. csv. show_timing If true 'time per pixel' information is added in the output (default). Hide this information is useful to compare results with diff command
output. tile. function_settings. csv. chemical_concentration Calculate and print chemical concentration
output. tile. function_settings. csv. delimiter Separator between fields
output. tile. function_settings. csv. compression If true output is compressed in GZ format (gz extension is automatically added)
output. tile. function_settings. csv. show_timing If true 'time per pixel' information is added in the output (default). Hide this information is useful to compare results with diff command
controller. segment_range This parameter allows specifying a range of segments that will be inverted/processed. If it is a single number a specific retrieval will be processed. Use -1 as undefined. For example [15, -1] will process all retrievals starting by the segment #15
controller. debug. perform_retrieval Allowing controller to call retrieval. If this parameter is false the framework will work without inverting the data and it will force to use only none output functions (results will not be printed). This parameter is useful for debug the framework or prepare input data
controller. debug. compilation_information Print compilation information at the beginning of the process
controller. debug. tracking_memory_stream Stream where printing all information about memory allocated during the execution (debugging information)
controller. stream Stream to dump controller information
controller. mpi. maximum_job_time Maximum number of seconds that the master node will wait for a working node for obtaining results. After this time (specified in seconds) if the job is not finished the controller will kill the task and the segment will be skipped
controller. mpi. polling_time Number of seconds that the master node wait after checking the workers

4.1.1.2. Extending settings: command line, import and template statements

The first argument of GRASP has to be the settings file but this file can be modified by another mechanism proposed by the settings module. The main way is by the command line, which allows to replace every settings parameter with "dot" syntax (replacing indentation by "dot" symbol and colon symbol by equal). All parameters that have been defined before and being replaced in the command line will cause a "note" information during the execution of GRASP. That sentence is just to inform the user that command line arguments always have higher priority than parameters in settings files. The value from command line will be the value that will be used to run the code. The command line is a very powerful feature to be used in the production scripts.

But the command line is not the only way to modify GRASP settings files. Settings files accept "import" and "template" statement. These statements could look similar but theirs behaviour is a bit different. Both of them allow defining other settings files that are read before the current one, but in case of import, the settings can not be overwritten. The template statement allows loading other settings files and then, modifying some settings to customize the loaded file. It is necessary to take into account that these statements can be used in cascade, creating problems to debug the code. So please use these statements carefully.

4.1.1.3. Streams

Some of the settings parameters are defined as "streams". GRASP output streams allow users to create dynamic names avoiding overwriting files or having to change the filenames each time they execute GRASP. When a description of a parameter is defined as "output stream", the user can set up a regular output path, for example ./folder/file.extension or use the "magic" behind the output streams by using a wildcard that will be replaced by dynamic values. For instance:

                    

output:
    segment:
        function: hdf
        stream: "GRASP_Banizoumbou_20080101_20080331_2x2+3286+1376.hdf"
    tile:
        function: [ ascii, hdf ] 
        stream: [ "GRASP_Banizoumbou_20080101_20080331.txt",
                "GRASP_Banizoumbou_20080101_20080331.hdf" ]
                

That definition is ok for many cases but if many tiles or segments are going to be processed, the fixed names will produce name collisions (the content of some files will be overwritten during the process). It would be tedious (and not always possible) for the user to change by himself the dates or other numeric substrings in the file names. For this reason, the configuration system provides some wildcards that will be automatically replaced with the given values depending on the state of the processing. These wildcards are marked with curly braces and their names are quite self-explanatory. The previous example can be rewritten in a more generic way using the stream wildcards:

                    

output:
    segment:
        function: hdf
        stream: "GRASP_Banizoumbou_{tile_from(%Y%m%d)}_{tile_to(%Y%m%d)}_{segment_nx}x{segment_ny}//
        +{segment_corner_column(4)}+{segment_corner_row(4)}.hdf"
    tile:
        function: [ ascii, hdf ] 
        stream: [ "GRASP_Banizoumbou_{tile_from(%Y%m%d)}_{tile_to(%Y%m%d)}.txt",
                "GRASP_Banizoumbou_{tile_from(%Y%m%d)}_{tile_to(%Y%m%d)}.hdf" ]
                

In addition, wildcards will provide the user the capability to set some system streams. If you use the values "true", "screen", "t" or "1", the information will be printed in the terminal (stdout). If the stream is set to "false", "none", "null", "f" or 0, nothing will be printed (like redirect to /dev/null). Finally, using the value "stderr" or "-1", the output will be redirected to the standard error output.

The following list shows all available wildcards that can be used for creating dynamic output filenames:

  • auto(N): itime x icol x irow with N zeros at the left
  • icol(N): current column number with N zeros at the left
  • irow(N): current row number with N zeros at the left
  • itime(N): current time number with N zeros at the left
  • iinversion(N): current inversion id with N zeros at the left
  • segment_nx(N): number of X elements per segment with N zeros at the left
  • segment_ny(N): number of Y elements per segment with N zeros at the left
  • segment_nt(N): number of T elements per segment with N zeros at the left
  • tile_from(FORMAT): start tile date in FORMAT. By default FORMAT is %FT%H_%M_%SZ
  • tile_to(FORMAT): final tile date in FORMAT. By default FORMAT is %FT%H_%M_%SZ
  • tile_corner_column(N): number of the corner (column) of the segment defined in settings file. Requirement: Input data have to be defined using input.corner instead of input.center
  • tile_corner_row(N): number of the corner of (row) the segment defined in settings file. Requirement: Input data have to be defined using input.corner instead of input.center
  • tile_center_longitude(FORMAT): longitude of the center of the tile defined in settings file. Requirement: Input data have to be defined using input.center instead of input.corner
  • tile_center_latitude(FORMAT): latitude of the center of the tile defined in settings file. Requirement: Input data have to be defined using input.center instead of input.corner
  • tile_coordinate_x(I): x input reference of center of the tile defined in settings file. It can be defined by corner or latitude. If is N in case it was defined by corner or 0.I in case it was defined like center
  • tile_coordinate_y(I): y input reference of center of the tile defined in settings file. It can be defined by corner or latitude. If is N in case it was defined by corner or 0.I in case it was defined like center
  • tile_width(N): Number of X elements in tile with N zeros at the left
  • tile_height(N): Number of Y elements in tile with N zeros at the left
  • segment_corner_column(N): number of column of the segment corner with N zeros at the left. Requirement: Input data have to be defined using input.corner instead of input.center
  • segment_corner_row(N): number of row of the segment corner with N zeros at the left. Requirement: Input data have to be defined using input.corner instead of input.center
  • segment_first_date(FORMAT): date of first pixel inside the segment in FORMAT. By default FORMAT is %FT%H_%M_%SZ
  • segment_last_date(FORMAT): date of last pixel inside the segment in FORMAT. By default FORMAT is %FT%H_%M_%SZ
  • iteration(N): Number of iterations with N zeros at the left. (Note. In case of single pixel it returns the number of iterations of first pixel)
  • settings_filename: the name of settings file used to run the retrieval
  • version: version of grasp if it is compiled with saving this information
  • branch: git branch of grasp if it is compiled with saving this information
  • commit: reference of git commit of grasp if it is compiled with saving this information
  • constants_set: constants set used in compilation time
  • pwd: this is replaced by current folder and it is only valid at the beginning of the stream definition
  • yml: this is replaced by current folder of main configuration file and it is only valid at the beginning of the stream definition

4.1.2. Retrieved characteristics

The ensemble of characteristics available to be retrieved or simulated by GRASP is open to user selection in the settings file inside retrieval.constraints.characteristic section. An example of the general structure of them is showed below:

                
retrieval:
    constraints:
        characteristic[1]:
            type: characteristic_name
            retrieved: true
            mode[1]:
                initial_guess:
                    value:                          [0.0, 0.0, 0.0, 0.0, ...]
                    min:                            [0.0, 0.0, 0.0, 0.0, ...]
                    max:                            [0.0, 0.0, 0.0, 0.0, ...]
                    index_of_wavelength_involved:   [0.0, 0.0, 0.0, 0.0, ...]
                single_pixel:
                    smoothness_constraints:
                        difference_order: 0.0
                        lagrange_multiplier: 0.0
                multi_pixel:
                    smoothness_constraints:
                        derivative_order_of_X_variability:    0.0
                        lagrange_multiplier_of_X_variability: 0.0
                        derivative_order_of_Y_variability:    0.0
                        lagrange_multiplier_of_Y_variability: 0.0
                        derivative_order_of_T_variability:    0.0
                        lagrange_multiplier_of_T_variability: 0.0
                    
            mode[2]:
            
                ...
            
            mode[3]:
            
                ...
                
            ...
                
                

The number assigned at each characteristic has no relevance at all while coherence is maintained. The retrieved field can be set to 'true' if the corresponding characteristic is going to be a retrieved parameter in the inversion; or to 'false' if it is just a fixed value for the forward simulation. The number of modes included in each characteristic depends on the nature of each one. For characteristics representing optical or mycrophisical aerosol properties (size distribution or refractive index for example), the number of modes corresponds to the number of aerosol modes selected for the retrieval/simulation (tipically one or two if fine/coarse distinction is made). For characteristics representing surface properties three modes will be needed to describe the associated model; except for polarization that only one is nedeed. The fields present in initial_guess part contain one element for each wavelength following the structure provided in the SDATA file in the case of optical wavelength dependent characteristics; one element for each bin for size distribution related characteristics; and one element for other mycrophysical magnitudes or non-wavelength dependent optical characteristics. In the former case, index_of_wavelength_involved should be filled with zeros for all the corresponding bins. The rest of the elements included in single_pixel or multi_pixel parts are always formed by one single element.

The list of available characteristics can be found below:

Table 4.2. Available GRASP characteristics

Characteristic NameDescription
size_distribution_triangle_bins Normalized Size Distribution dV / dlnr at "triangle" bins
size_distribution_precalculated_lognormal

Normalized Size Distribution dV / dlnr for precalculated lognormal bins

size_distribution_lognormalParameters of bi - modal Lognormal Size Distribution dV / dlnr
aerosol_model_concentrationAerosol model concentration
real_part_of_refractive_index_spectral_dependentSpectral dependent Real part of complex refractive index
real_part_of_refractive_index_constantComplex Refractive Index Real part is spectrally constant
particle_component_volume_fractions_linear_mixtureReal part of complex refractive index is mixture
particle_component_fractions_chemical_mixtureChemistry, fraction of: water, fslbl, finslbl, soot, iron
imaginary_part_of_refractive_index_spectral_dependentSpectral dependent Imaginary part of complex refractive index
imaginary_part_of_refractive_index_constantComplex Refractive Index Imaginary part is spectrally constant
sphere_fractionFraction of spherical particles
aspect_ratio_distributionAxis Ratio Distribution
vertical_profile_parameter_heightScaling factor in case of exponential profile, mean height in case of gaussian distribution
vertical_profile_normalizedAerosol normalized vertical profile
aerosol_concentrationAerosol concentration
lidar_calibration_coefficientCalibration coefficient for lidar
vertical_profile_parameter_standard_deviationStandard deviation for aerosol vertical profile
surface_land_brdf_ross_liBRDF Land normalized parameters according to Ross and Li model
surface_land_brdf_rpvBRDF normalized parameters according to RPV model
surface_land_litvinovBRDF Land normalized parameters according to Litvinov model
surface_land_litvinov_fastBRDF Land normalized parameters according to Litvinov fast model
surface_land_polarized_maignan_breonBRDF Land normalized parameters according to Maignan and Breon model
surface_land_polarized_litvinovBPDF Land normalized parameters according to Litvinov model
surface_water_cox_munk_isoBRDF Water normalized parameters according to Cox and Munk model
surface_water_cox_munk_aniBRDF normalized parameters according to Maignan and Breon model
surface_water_litvinovBRDF Water normalized parameters according to Litvinov model


4.1.2.1. Initial guess through the algorithm

Understanding GRASP inversion procedure means understanding how the algorithm is starting from an initial guess and obtains a results array. This is an iterative procedure explained in the literature and introduced in Section 2.3, “GRASP Scientific Core algorithm”. The purpose of this section is to explain how initial guess is represented inside the code as an array which evolves in each iteration until getting the result array. The following diagram shows how initial guess is read from settings file and translated into an internal array in the code. This detail could look very technical and related with the development but understanding of some internal concepts of the code helps to understand how it works. The following diagram shows this transformation:

Figure 4.1. Translation of settings file into initial guess array

Translation of settings file into initial guess array


Since the initial guess is defined in the settings, it is defined once for the entire segment (one should remember that GRASP can work in a multi-pixel approach). Sometimes, it is desirable to have a different initial guess for each pixel. In that case the already described mechanisms to establish the initial guess are not enough, and GRASP proposes another useful tool to modify the initial guess of each segment. This information is considered as part of the input because it contains pixel dependent information and will be defined in detail in the Section 4.2.3, “Input information for characteristics”.

Once the initial guess is loaded, it is set as first array of characteristics to be retrieved. Retrieval process iterates over it until it gets the results. The results of GRASP retrieval is the array with the same shape as the initial guess, but containing the results of retrieval of these parameters to match them with SDATA file. Then, the GRASP forward model is called using this results array, GRASP obtains the rest of the derived results it provides. So, we can talk about two kinds of results: basic results and derived products. The following diagram shows this process:

Figure 4.2. Evolution of retrieved characteristics during GRASP processing

Evolution of retrieved characteristics during GRASP processing


To know all the details about the products obtained by GRASP, please see Section 4.3.1, “The list of GRASP output parameters”.

4.1.3. Noise simulation

In order to harmonize all possibilities to add random noise to the measurements that are going to be retrieved in GRASP a new setting called add_random_noise in the inversion.noises group has been created. There are three possible values that this setting can take: disable, measurement_fitting (default) and sdata. These options will produce a different effect depending on which retrieval mode (forward or inversion) has been selected. Here it is a small description of the meaning of the options:

  • disable: If this option is selected there is no addition of any noise even if standard_deviation_synthetic greater than 0. All the rest of the settings related to it will be ignored.
  • measurement_fitting: This is the default value and this option maintains the same behaviour of GRASP random noise as before: for both forward and inversion retrieval modes random noise is added to FS measurement vector to be fitted.
    • If the forward option is selected in retrieval.mode any noise addition will be made.
    • If inversion is selected, the noise is added in the FS vector and it is totally transparent to the user. There is no possibility to see what are the noisy measurements that are really being used in the retrieval. This option will always work no matter if the measurements in the sdata do really correspond or not with the ones that are being inverted (this is especially referred to the polarization measurements).
  • sdata: the sdata mode is used to enable the user to see what are the noises added to the measurements. This option will only work if the measurements in the sdata and the measurements that are really retrieved are exactly the same. In other case, meaning that the code is fitting different measurements than those in the input data, an error will arise and the code will not be executed.
    • If the forward option is selected in retrieval.mode, the added noise will be done over the FS vector. It will be seen in the column fit_X of the output file, or if output.sdata.dump setting is present in the yml file.
    • If inversion is selected, the noise will be added to the original sdata provided by the user independently of the forward model configuration. The noise measurements will be available for the user in the meas_X column of output file or via input.sdata.dump.

                
 mode: inversion
    
    inversion:
        regime: single_pixel
        
        convergence:
            minimization_convention: logarithm
            maximum_iterations_for_stopping: 35
            maximum_iterations_of_Levenberg-Marquardt: 35 
            threshold_for_stopping: 1.0e-3
            threshold_for_stopping_Q_iterations: 1e-12
            scale_for_finite_difference: 1.0e-5  
            normal_system_solver: sparse_matrix_solver 

        measurement_fitting:
            polarization: degree_of_polarization                

        noises:
            add_random_noise: measurement_fitting # | disable | sdata
            noise[1]:
                standard_deviation_synthetic: 0.1
                error_type:  relative
                standard_deviation:  0.05
                measurement_type[1]:
                    type: I
                    index_of_wavelength_involved: [ 1, 2, 3, 4 ]
            noise[2]:
                standard_deviation_synthetic: 0.01
                error_type:  absolute
                standard_deviation:  0.005
                measurement_type[1]:
                    type: tod
                    index_of_wavelength_involved: [ 1, 2, 3, 4 ]
...
                
                

4.2. Input module

The input module is responsible for reading the data and for setting the internal input data structures. GRASP provides a very flexible way to inject the input data offering developers the capability to create the input drivers. An input driver is an extension of GRASP that is added during compilation and is selected in the settings file. These generic drivers allow the developers to create a custom way to read a specific database and load it into the GRASP scientific algorithm. Therefore, GRASP can read infinite kinds of input databases, as much as drivers exist. Additionally, drivers can take the responsibility of performing some pre-processing actions (such as calibration corrections) before retrieving the data. Finally, when GRASP is used for massive data processing, drivers are extremely important since they help to prepare input data and load it into scientific module without the use of any intermediate files, everything is performed in memory. For general use, GRASP proposes a generic driver called the SData driver, which reads files in SData format (Sensor-Data format). These files are not standard, they have a specific format proposed by GRASP to start working with the code. This format is described in Section 4.2.1, “The SDATA format”.

Transformers are the second kind of extensions that the input module offers to users (and developers). They are called after getting input data for a segment, it means after calling the driver. The purpose of transformers is to add the capability to modify the segment after obtaining the data. An example of transformation is to load a climatology database and to modify the initial guess of each pixel to optimize the number of iterations needed to retrieve the data. Performing this action after getting the data allows reusing it between different drivers. Default installation of GRASP does not offer any transformers. All of them are considered as optional extensions and have to be externally added manually or using the grasp-manager.

4.2.1. The SDATA format

The SDATA (sensor data) format is the original input data format of the GRASP code. It is a simple text format designed by the science team at the early stages of the development of the scientific code and it's the easiest way to create test data.

In the context of the GRASP project, the SDATA format has a number of pros:

  • Designed by the scientific team at the origin of the GRASP project, it is well adapted to its needs.
  • It is very simple to describe.
  • It is a text format, and therefore portable, quite easy to check (for the accustomed eye!) and to edit and to make some quick experiments.
  • It is a piece of cake to read in Fortran :-)

It has also a number of cons:

  • It is not standard (a by-product of the "not designed by a committee" approach). No off-the-shelf library is available to parse it and validate it (meaning outside of the GRASP project).
  • It lacks of flexibility: the order and the number of values are fixed for one version and any change in the format to meet new requirements is likely to break the compatibility with the former versions.
  • It is fragile: a malformed file may easily make the code crash or produce mysterious bugs. The design of the format, while simple, makes it hard to develop a really reliable validator.
  • Comments are expected only after values, on the same line (after a colon sign)
  • While the format uses a text representation for the data, it contains lots of numeric values, with limited accuracy and with no comment. Large files are tedious to read and it's easy to make shift errors while reading and editing even for the experienced user.
  • Being a text format, it becomes very inefficient for large data volumes. While it can be compressed for archiving, it must be uncompressed for processing, and only sequential access is possible. It is still possible to perform regional processings (several dozens of thousands of pixels that cover more than a few hundreds of kilometres in both directions) with this format (it was actually done for the sake of necessity), but it stresses the computing system a lot and can't be scaled up to the global processing.

Whatever the number and seriousness of the cons, one of the design objectives is to keep the code simple and flexible allowing the scientific community to play with the code. It does not make sense to implement a driver for a single user who wants to do some tests with GRASP. That's why this easy format is maintained by the developer team.

In the following description, the elements in fixed-width font are the snippets of content. The numeric values in these snippets (e.g. in 2 2 2 : NX NY NT) are given only as examples.

Figure 4.3. An example of SDATA file

        

SDATA version 2.0
2   2   2  : NX NY NT

4  2008-01-04T13:15:00Z 70000.0  0  0 : NPIXELS TIMESTAMP HOBS_km NSURF IFGAS
1  1  1	 3286  1377   2.599	13.528   252.0	 100.0  6  0.443  0.490  0.565 ...
2  1  1	 3287  1377   2.657	13.528   242.0	 100.0  6  0.443  0.490  0.565 ...
1  2  1	 3286  1376   2.601	13.583   241.0	 100.0  6  0.443  0.490  0.565 ...
2  2  1	 3287  1376   2.658	13.583   239.0	 100.0  6  0.443  0.490  0.565 ...


4  2008-01-06T13:02:41Z 70000.0  0  0 : NPIXELS TIMESTAMP HOBS_km NSURF IFGAS
1  1  1	 3286  1377   2.599	13.528   252.0	 100.0  6  0.443  0.490  0.565 ...
2  1  1	 3287  1377   2.657	13.528   242.0	 100.0  6  0.443  0.490  0.565 ...
1  2  1	 3286  1376   2.601	13.583   241.0	 100.0  6  0.443  0.490  0.565 ...
2  2  1	 3287  1376   2.658	13.583   239.0	 100.0  6  0.443  0.490  0.565 ...
                

SDATA files have a simple structure:

  1. The first line is the FILE HEADER. It contains a magic identifier SDATA followed by a version number.
  2. The second line is the SEGMENT HEADER. It contains three numbers, NX, NY and NT, the spatial and temporal dimensions of the segment that this SDATA file represents. You can notice a colon and names of fields after the values. This is the way how comments are written in the SDATA files. Everything starting from the colon will be ignored by the SDATA parser.
  3. An empty line follows the SEGMENT HEADER.
  4. Then comes the first CELL (group of neighbouring pixels) of the SEGMENT. Each CELL has at most NX*NY pixels (but it may have less, for various reasons: cloudy pixels that have been filtered, missing pixels, etc.). The number of CELLs in the SEGMENT is given by the NT number provided in the SEGMENT HEADER.

Table 4.3. The SDATA main structure

Field NameField Content
FILE HEADER SDATA version 2.0
SEGMENT HEADER

2 2 2 : NX NY NT

empty line 
CELL 1cell content, look for CELL structure
empty line 
CELL 2cell content, look for CELL structure
empty line 
... 
CELL itcell content, look for CELL structure
empty line 
... 
CELL NTcell content, look for CELL structure
empty line 


A CELL is a set of neighbouring pixels, that form the base of a SEGMENT. Each CELL has a HEADER and a number of PIXELs, supposedly acquired at the same time.

The CELL HEADER contains:

  1. The number of pixels in the CELL (NPIXELS). It may not be larger than NX*NY
  2. The timestamp of acquisition of the pixels, in the ISO8601 time format.
  3. A "height" of observation, in metres. The value here is a bit weird (70000), and doesn't correspond to the satellite altitude (that is at least 10 times larger). Actually, the value doesn't really matter as long as it is large. Historically, the scientific team has used this value of 70000 in many SDATA files.
  4. Two values for the number of surface and gas parameters. These two values are currently not documented and can be set to 0 for the moment.
  5. Comments starting with a colon.

Table 4.4. The CELL structure

Field NameField Content
CELL HEADER
4 2008-01-04T13:15:00Z 70000.0 0 0 : NPIXELS TIMESTAMP HOBS NSURF IFGAS
PIXEL 1a line of values, look for PIXEL structure
PIXEL 2a line of values, look for PIXEL structure
... 
PIXEL NPIXELSa line of values, look for PIXEL structure


Each line of data after the CELL HEADER represents exactly one pixel, with all its fields. The Table 4.5, “The PIXEL structure” describes the order and type of these fields. For the types, the Fortran notation is used: array types are described with the dimensions of arrays between parentheses, and the ordering is such that the first index increases faster. Indices start from 1, not from 0 like in C. For instance, when one reads real(nwl) for wavelengths, that means that one has to read a list of nwl real values that represent wavelengths.

Table 4.5. The PIXEL structure

Field TypeVariable Name (in source code)Field Content
integerpixel[ipix].ixcoordinate x in the current cell, starting at 1 (in the direction EW)
integerpixel[ipix].iycoordinate y in the current cell, starting at 1 (in the direction NS)
integerpixel[ipix].icloudycloud flag: 0 = cloud, 1 = clear [a]
integerpixel[ipix].icolcolumn of the pixel in its original grid or database (can be set to 0 when not relevant)[b]
integerpixel[ipix].irowline of the pixel in its original grid or database (can be set to 0 when not relevant)[b]
realpixel[ipix].x

longitude of the pixel, in decimal degrees, in the range [-180..180]

0: Greenwich meridian, east of Greenwich: positive, west of Greenwich: negative

realpixel[ipix].y

latitude of the pixel, in decimal degrees, in the range [-90..90]

realpixel[ipix].MASLaltitude of the ground, in metres (MASL: metres above sea level)
realpixel[ipix].land_percentpercentage of land, in the range [0 (sea) .. 100 (land)]. Intermediate values correspond to coastal pixels
integerpixel[ipix].nwlnumber of available wavelengths (nwl)
real(nwl)pixel[ipix].meas[nwl].wllist of wavelengths, in micrometers
integer(nwl)pixel[ipix].meas[nwl].nipnumber of types of measurements for each wavelength (nip)
integer(nip, nwl)pixel[ipix].meas[nwl].meas_type[nip]

list of types of measurements meas_type (see Table 4.6, “Types of measurements”)

The ordering is as follows: meas_type(1, wl1) meas_type(2, wl1) ... meas_type(nip, wl1) meas_type(1, wl2) ... meas_type(nip, wln)

integer(nip, nwl)pixel[ipix].meas[nwl].nbvm[nip]

number of valid measurements (nbvm), for each type of measurement and for each wavelength

The ordering is as follows: nbvm(1, wl1) nbvm(2, wl1) ... nbvm(nip, wl1) nbvm(1, wl2) ... meas_type(nip, wln)

real(nwl)pixel[ipix].meas[nwl].szasolar/sounding zenith angle (sza or θs) in decimal degrees ([0..90]), for each wavelength
real(nbvm, nip, nwl)pixel[ipix].meas[nwl].thetav[nip][nbvm]

viewing zenith angle (vza or θv) in decimal degrees ([0..90] for PARASOL, TBD for PHOTOMETERS)

The ordering is as follows: θv(1, 1, wl1) vza(2, 1, wl1) ... θv(nbvm, 1, wl1) θv(1, 2, wl1) ... θv(nbvm, nip, wln)

In case of lidar or vertical observations this field contains altitudes or ranges of the observations in meters.

real(nbvm, nip, nwl)pixel[ipix].meas[nwl].phi[nip][nbvm]

relative azimuth angle (raa or Δϕ) in decimal degrees ([-180..180] or [0..360])

The ordering is as follows: Δϕ(1, 1, wl1) Δϕ(2, 1, wl1) ... Δϕ(nbvm, 1, wl1) Δϕ(1, 2, wl1) ... Δϕ(nbvm, nip, wln)

real(nbvm, nip, nwl)

pixel[ipix].meas[nwl]

.tau[nbvm];...;.P[nbvm]

measurements (depending on meas_type), for each wavelength

The ordering is as follows: meas(1, 1, wl1) meas(2, 1, wl1) ... meas(nbvm, 1, wl1) meas(1, 2, wl1) ... meas(nbvm, nip, wln)

Where meas can be: tau ... P

real(nsurf, nwl)pixel[ipix].meas[nwl].groundpar[nsurf]

ground parameters

This part can be ignored for now. The ordering is as follows: groundpar(1, wl1) groundpar(2, wl1) ... groundpar(nsurf, wl1) groundpar(1, wl2) ... groundpar(nsurf, wln)

real(nwl)pixel[ipix].meas[nwl].gaspar

gas absorption (tau gases) or molecular depolarization ratio. This parameter has to be provided only if the setting IFGAS (in the CELL HEADER) is set to 1.

The ordering is as follows: gaspar(wl1) gaspar(wl2) ... gaspar(wln)

integer(nip, nwl)pixel[ipix].meas[nwl].ifcov[nip]

ifcov (1 if a covariance matrix is available, 0 otherwise)

The ordering is as follows: ifcov(1, wl1) ifcov(2, wl1) ... ifcov(nip, wl1) ifcov(1, wl2) ... ifcov(nip, wln)

real(nbvm, nip, nwl)pixel[ipix].meas[nwl].cmtrx[nip][nbvm]

cmtrx (diagonal of covariance matrix, also known as Ω). These values have to be skipped if ifcov=0

The ordering is as follows: cmtrx(1, 1, wl1) cmtrx(2, 1, wl1) ... cmtrx(nbvm[c], 1, wl1) cmtrx(1, 2, wl2) ... cmtrx(nbvm[c], nip, wln)

integer(nip, nwl)pixel[ipix].meas[nwl].ifmp[nip]

ifmp (1 if a vertical profile (mprof) is available, 0 otherwise)

The ordering is as follows: ifmp(1, wl1) ifmp(2, wl1) ... ifmp(nip, wl1) ifmp(1, wl2) ... ifmp(nip, wln)

real(nbvm, nip, nwl)pixel[ipix].meas[nwl].mprof[nip][nbvm]

mprof (vertical profile of Rayleigh backscattering). These values have to be skipped if ifmp=0

The ordering is as follows: mprof(1, 1, wl1) mprof(2, 1, wl1) ... mprof(nbvm[d], 1, wl1) mprof(1, 2, wl2) ... mprof(nbvm[d], nip, wln)

[a] This fairly counter-intuitive coding has a reason: the cloud flag was at first intended to be a general processing flag (0 = pixel not to be processed, 1 = to be processed), cloud contamination is only one particular case. Now the flag is limited to cloud screening, but unfortunately the coding couldn't be changed right away. Since the framework is still in development, it is planned to correct this unnatural feature in the near future.

[b] These fields are actually not used by the processing and therefore the SDATA implementer is free to put whatever he or she likes here (e.g. 0 for non-gridded data). They are intended mainly for documentation and debugging. For satellite data, they make it possible to retrieve the pixel original information in the original database.

[c] nbvm is actually to be multiplied by ifcov(ip, iwl). If this last number equals 0, the array reduces to an empty set and no value is to be read.

[d] nbvm is actually to be multiplied by ifmp(ip, iwl). If this last number equals 0, the array reduces to an empty set and no value is to be read.


The field pixel[ipix].meas[nwl].meas_type[nip] of pixel structure is a special code which defines the type of measure. The following table describes the valid codes and their interpretation:

Table 4.6. Types of measurements

Constant Name (used in source code)Value (SDATA 2.0)Meaning
MEAS_TYPE_UNKNOWN 0 The measurement type is invalid or not yet implemented
MEAS_TYPE_TOD 11 Total Optical Depth
MEAS_TYPE_AOD 12 Aerosol Optical Depth
MEAS_TYPE_ABS 13 Aerosol absorption optical depth
MEAS_TYPE_P11 21 Phase Matrix Element P11
MEAS_TYPE_P12 22 Phase Matrix Element P12
MEAS_TYPE_P22 23 Phase Matrix Element P22
MEAS_TYPE_P33 24 Phase Matrix Element P33
MEAS_TYPE_P34 25 Phase Matrix Element P34
MEAS_TYPE_P44 26 Phase Matrix Element P44
MEAS_TYPE_P11_rel_ang 27 p11/p11(given_angle) phase matrix element
MEAS_TYPE_P12_rel 28 -p12/p11 phase matrix element
MEAS_TYPE_LS 31 Lidar Signal
MEAS_TYPE_RL 32 Raman Lidar Signal
MEAS_TYPE_DP 35 Volume Depolarization Ratio
MEAS_TYPE_VEXT 36 Vertical Extinction profile
MEAS_TYPE_VBS 39 Vertical Backscatter profile
MEAS_TYPE_I 41 Normalized Radiance I [a]
MEAS_TYPE_Q 42 Polarized radiance Q[a]
MEAS_TYPE_U 43 Polarized Radiance U[a]
MEAS_TYPE_P 44 Polarization Rate: P = sqrt(Q*Q + U*U)/I
MEAS_TYPE_P 44 Polarization Rate: P = sqrt(Q*Q + U*U)
MEAS_TYPE_I_rel_sum 45 Relative Stokes parameter I/sum(I(1:NBVM)), where NBVM is the total number of provided angles
MEAS_TYPE_P_rel 46 Linear polarization sqrt(Q*Q+U*U)/I

[a] All the Stokes Parameters are to be expressed as reduced quantities, without dimension

Equation 4.1. Conversion from absolute radiances to normalized, reduced radiances

I = radiance * π / E0


where radiance is the radiance of the instrument, and E0 the solar spectral flux, that may be both in mW / (m^2 * sr * nm) or equivalent units


4.2.2. Angle definition

This chapter main goal is to describe how the angles should be defined to be used inside of GRASP code. The universal spirit of GRASP, where many different instruments coexist (from satellite to ground based measurements), creates challenges to define a homogeneous way of defining the angles, keeping a unique geometry. As it is shown in the figure Figure 4.4, “Definition of GRASP geometry” GRASP angles are defined to be considered as "normal" for satellite reference.

Figure 4.4. Definition of GRASP geometry

Definition of GRASP geometry


Since GRASP angles are defined using a satellite reference, it provokes some problems to define what we could call as an intuitive "ground based" reference system. That is why we are going to put special emphasis on the definition of the angles to these less intuitive applications. The intuitive reference for "ground based" measurements, in spherical geometry, is given as follows:

  • θgb zenith angle: with the zero established in the zenith
  • ϕgb azimuth angle: with the zero considered in the sun position

where the sub index "gb" makes reference to "ground based".

The conversion to the GRASP geometry is done as follows:

Equation 4.2. Conversion from θgb (ground based) to θG (GRASP)

θG = 180° - θgb


Equation 4.3. Conversion from ϕgb (ground based) to ϕG (GRASP)

ϕG = 180° + ϕgb


Here we propose some examples for better understanding of the process. Before defining the measurement angles introduced in the code, both "intuitive" and "GRASP", we need first to consider the instrument viewing angle for each scenario (θv, ϕv). The following figure and table will provide some examples of angles defined for the ground based applications.

Figure 4.5. Ground based angles definition example

Ground based angles definition example


Table 4.7. Specific examples in ground based angle definition example

ExampleθsθvϕvθgbϕgbθGϕG
125°25°155°180°
225°12.5°12.5°167.5°180°
325°25°180°180°
425°37.5°12.5°180°167.5°
525°50°25°180°155°
625°90°65°180°115°
725°30°25°30°155°210°
825°90°25°90°155°270°

Since sunphptometers are widely used with GRASP, the following table provides information specifically about these instruments, considering the instrument viewing angle for each scenario (θv, ϕv). They can be understood as the "movements of the motors". The process will be as follow:

instrument viewing angle -> angle in (intuitive) ground based -> angle in GRASP

Table 4.8. Sunphotometer angle description

Measure typeAngleInst. View.RangeGr. BasedRangeGRASPRange
Direct sunθθv = 0°[0°]θgb = θs[0° -- 90°][a]θG = 180° - θs[180° -- 90°][a]
ϕϕv = 0°[0°]ϕgb = 0°[0°]ϕG = 180°[180°]
Almucantarθθv = 0°[0°]θgb = θss]θG = 180° - θs[180° -- θs]
ϕϕv = 3°, 3.5°, 4°, 5°, ..., 90°, ... 180°[0° - 180°]ϕgb = 3°, 3.5°, 4°, 5°, ..., 90°, ... 180°[0° -- 180°]ϕG = 183°, 183.5°, 184°, 185°, ... 270°, ... 360°[180° -- 360°]
Principal plane measurement: Before the zenith θθv = -6°, ... , -3°, 3°, θmax < θs [-6° -- (θmax < θs)]θgb = θs + 6° ... θs + 3°, θs - 3°, ... θs - 6°, θs - θmax[(θs + 6°) -- 0][b]θG = 180° - θs - 6°, ... 180° - θs - 3°, 180° - θs + 3°, ... 180° - θs + θmax [(180° - θs - 6°) -- 180°]
ϕϕv = 0°[0°]ϕgb = 0°[0°]ϕG = 180°[180°]
Principal plane measurement: After the zenith θθv = θmin < θs , ... 140° or θmax - θs > 90°[(θmin < θs) -- 140°]θgb = θmin - θs ... 140° - θs or θmax - θs[0° -- 90°][b]θG = 180° + θs - θmin, ... 180° + θs - 140° or 180° + θs - θmax [180° -- 90°][b]
ϕϕv = 0°[180°]ϕgb = 180°[180°]ϕG = 0°[0°]
Polarized principal plane measurement: before the zenith[d]θ  θgb = 85°, 80°, 75°, ... 10°, 5°, 0° θG = 95°, 100°, 105°, ... 170°, 175°, 180° 
ϕ  ϕgb = 0° ϕG = 180° 
Polarized principal plane measurement: after the zenith[d]θ  θgb = 0°, 5°, 10°, ... 75°, 80°, 85° θG = 180°, 175°, 170°, ... 105°, 100°, 95° 
ϕ  ϕgb = 180° ϕG = 0° 

[a] θs refers to the solar zenith angle (for different measurements)

[b] decreasing values

[c] decreasing values

[d] The data of the polarized principal plane correspond always to fixed points in the sky and it is given for the instrument in so-called ground based coordinates.


In the case of nephelometer data angle definition, it is a bit different since only θ angle has to be defined, the rest of the angles will be ignored. To provide nephelometer data (phase matrix), the conversion to the GRASP geometry is done as follows:

Equation 4.4. Conversion from θn (nephelometer scattering angle) to θG (GRASP)

θG = 180° - θn


4.2.3. Input information for characteristics

While in many cases the input data is just represented by SDATA information, sometimes it is necessary to provide extra information for each pixel. The way to inject extra information to the algorithm is to do it through the initial guess. Sometimes, this is just because the user wants to provide different initial guesses for each pixel, or sometimes it is just input information. In the first case, this can be easily done by the imagedat files that will be described below. On the other hand, when a charactacteristic is defined as “retrieved=false”, it is taken as input information.

The “segment_imagedat” tool enables the user to provide an ASCII file containing the complete set of initial guesses for all the pixels. Following illustration is an example of how the imagedat files looks like:

Figure 4.6. Example of a possible use of imagedat

Example of a possible use of imagedat


The columns in this ASCII file correspond to the pixels, and each row is associated with a different characteristic. The first column is just the enumeration of characteristics starting by 1. The total number of both columns and rows has to be consistent with the input file. The value “-999” can be assigned to the initial guesses which the user does not need to modify, in this case the value defined in the characteristics section of the settings file will be taken by default. This mechanism provides a very versatile tool to inject pixel-dependent data to the retrieval code, for adjusting the initial guess or just for providing external input information from climatologies, models, etc.

4.2.4. How to prepare the photometer data

Sunphotometers are widely used with GRASP. They take measurements of sky radiance and direct sun. Many inversion strategies can be used to retrieve sunphotometer data, but this section will explain how to define input data. The information that runs inside of GRASP has to be pre-processed in order to screen clouds, calibrate and normalize the data.

Following Section 4.2.1, “The SDATA format”, the direct sun measurements can be described as AOD or TOD defined in Table 4.6, “Types of measurements” as measurements of type 11 or 12. At this point it is needed to take into account that if "ifgas" field is defined as 1 in the case of AOD, no gaseous absorption optical depth will be accounted, but in the case of TOD they will be subtracted. The gases also affect the radiance measurements but in lower magnitude. In the case of TOD + radiances with ifgas=1, the same model will be applied to all measurements. If AOD is used, some (minor) incongruences could come from the use of different models to calculate gases for AOD and for radiances.

Radiance measurements are defined with the constant MEAS_TYPE_I(41). Polarized measurements can be defined as Q,U (42, 43) or as polarization rate (44). It is also important to check how polarized data is going to be manipulated in the retrieval code based on inversion strategy defined in the settings file.

4.2.5. How to prepare the lidar data

Note that all processing will be considering range corrected profile for one wavelength. The procedure for other profiles from different wavelengths are exactly the same. Range corrected profile implies that at least the background noise was subtracted and altitude corrections were applied to the raw signal, but if it's possible to consider all other corrections (electrical noise and overlap correction, dead time correction, gluing analog and photon-counting signals and all others that your system may have), you should apply them.

Step 1. Background noise subtraction and range correction.

Let B' be the estimation of the background noise. Usually B' is estimated as P(ZB), accumulated and averaged around selected altitude ZBZ, much higher than the maximum altitude of lidar extraction in step 2 (Zmax). For example with maximum altitude Zmax=15km, it is averaged around 30 km and accumulated for the whole period of lidar observation. The noise and range corrected signal will be:

S(Zi) = (P(Zi)-B')*Zi2

Step 2. Altitude range selection and signal cropping.

The minimum Zmin and maximum Zmax altitudes are selected and the signal is cropped, so Zmin<Zi<Zmax. The minimum altitude should be selected as low as possible, preferably in the region where the overlap correction could be correctly applied. The Zmax should be selected from the following considerations: maximum altitude where the noise levels of lidar measurement are acceptable and the amount of atmospheric aerosols is still noticeable.

In case if the lidar used is inclined The sounding zenith angle should be provided for the corresponding wavelenghts. The angle value should be placed in the position of the SZA corresponding to this wavelength. All SZA's for all vertical profile measurements should be the same. Vertically pointed lidar should have this value set to 0. Note: Keep in mind that retrieved aerosol vertical profiles will be retrieved for vertically projected altitudes, i.e. if your lidar is inclined, the maximum altitude of the profile will be equivalent to Zmax*cos(SZA).

Step 3. Backscatter of molecular profile

Backscatter vertical profile is calculated inside of GRASP based on standard atmosphere model for each lidar wavelength[1]:

βmol(Zi,λ) = N * σ(λ) * (P(Zi)/PSA) * (TSA/T(Zi))

        where:

  • N molecular number density
  • σ(λ) total Rayleigh cross section per molecule, which analytical formula can be written like σ(λ) = A * λ -B-C*λ-D/λ
  • Psa pressure of standard atmosphere model
  • Tsa temperature of standard atmosphere model
  • P(Zi) pressure profile of atmosphere
  • T(Zi) temperature profile of atmosphere

Step 4. Reducing the number of points in profiles

GRASP/GARRLiC can use arbitrary altitude/range scale. However, to keep number of retrieved parameters reasonable and to fight higher noise contamination of lidar signals at higher altitudes, it is recommended to use a logarithmical altitude/range scale with NZ points to represent aerosol profiles in the atmosphere. For that, we have to present all vectors (altitude/range vector and profiles of lidar signals) in logarithmically equidistant manner.

1. Move to logarithmic scale and find altitude/range step:

logarithmic scale: Zi lg = lg(Zi)

step in log scale: ΔZ = (Zmax lg - Zminlg)/NZ

logarithmic altitude ranges (from hk to hk+1 ) for averaging data in logarithmically equidistant manner, k = 1 .. Nz: hk = Z0 lg + (k - 1) * ΔZ

2. Average the data profiles

Ak = (Σj=1n Aj(hk, hk+1) ) / n

        where:

  • A altitude vector or profile of lidar signal or molecular backscatter
  • n number of points inside logarithmic altitude ranges

After such procedure, the number of points in three main vectors reduced to NZ points in logarithmically equidistant manner. These values should be placed in places corresponding to the zenith wieving (vza or θv) angles for the vavelenths corresponding to lidar or vertical measurements.

Hint: the altitude/range vectors for measurements at all wavelengths have to be the same.

Step 5. Profile normalization

Values of lidar signals (except for volume depolarization profiles) vary from instrument to instrument, from detector to detector, that is why GRASP/GARRLiC requires normalized lidar signal. For consistency with the molecular optical depth, the profile of the molecular backscatter inside the code has to be normalized as well. Normalized lidar and backscatter profiles:

A'k = Ak / ∫ZminZmax Ak dZ

        where A represents profile of lidar signal or molecular backscatter.

Caution: integration have to be done using meters in altitudes.

At the end, for each wavelength you have to have normalized lidar profiles and altitude vector. Volume depolarization profiles don't need to be normalized, the only requirement for such observations is to be presented in the percentage range i.e. [1.0e-9, 100]

References: Anthony Bucholtz, Rayleigh-scattering calculations for the terrestrial atmosphere, Optical Society of America, 1995.

4.2.6. How to prepare nephelometer data

GRASP is able to retrieve and simulate nephelometer measurements. This section is devoted to the definition of the convention to properly account the units and measurement geometry of these instruments within GRASP standards.

The zero reference of scattering angles in GRASP is 180 °. Thus, in order to introduce the nephelometer geometry in GRASP sdata, only θ angle (zenith angle) has to be defined, the rest of the angles will be ignored. The conversion to the GRASP geometry is done as follows:

θG = 180° - θn

Where θG corresponds to sdata angles and θn to the scattering angles of the nephelometer.

GRASP measurement input corresponds to the normalized Phase Matrix element P1,1, and the relation with scattering function is described below.

GRASP scattering phase function units are:

[F1,1(λ, Θ)] = 1/μm

GRASP scattering phase function units are:

[σ] = 1/μm

The relationship between both is represented by the following expression:

σ = (1/2) * ∫π0 F1,1(λ, Θ) * sin( Θ) d Θ

Which means that the normalized Phase Matrix element P1,1 is defined as:

(1/2) * ∫π0 P1,1(λ, Θ) * sin( Θ) d Θ = 1 P1,1(λ, Θ) = F1,1(λ, Θ)/σ

If polar neph measurements are used, an additional factor of 4* π is necessary to convert into the units in GRASP: μm-1Str-1 to μm-1.

4.3. Output module

The output module is responsible for managing results for each segment or entire tile and driving them to the correspond destination (file or screen). There are two kinds of main output structure into GRASP: tile and segment. Segment output structure represents the output results obtained from the retrieval library. Then, core unit compacts it and stores it in a tile output structure, which contains the results of the entire process.

Output module can be extended in the same way as input module. In the case of the output there are three kinds of extensions:

  • output segment function: it is called after retrieving a segment and is called with the results of that single segment.
  • output current function: it is called after retrieving each segment but it is called with the partial tile processed until that moment. In each call to this function, it gets a partial tile nearest to completion. Last call to that function will send the entire tile results.
  • output tile function: at the end of the process, a function is called sending to it the entire tile results. It can print a complete map of the process.

As the user can see, GRASP is very flexible in the way to work with the output. In a high optimized process, it can be adapted directly to the format of the output database. By default, GRASP comes with some ASCII output functions, which allows to print the results in a readable way (ASCII or specific GRASP format), into a file (using stream library) or on the screen. Additional extensions can be added to get the output in different formats such us HDF, NetCDF, png plots...

4.3.1. The list of GRASP output parameters

The output from GRASP is quite complex and strongly dependent on the settings used. As it was discussed in Section 4.1.2, “Retrieved characteristics”, there are two kinds of output products: direct and derived. Direct results are the values directly inverted in the initial guess array, then after obtaining them, forward model is called one more time to obtain the derived products. The list of direct and derived products obtained depends on the data and the inversion strategy selected (defined in the settings file). This chapter contains a complete list of products that can be obtained with GRASP, but it does not mean that all of them can be obtained with all inversion strategies.

In the internal GRASP output structures there are some products that are repeated between direct products structure (array with the same shape as initial guess) and derived products. The reason is that for some applications they are direct information, for others they are derived results. This depends on input data and inversion strategy, defined in the settings file.

Before defining the output, we are going to define the size of some arrays. This information is needed to understand the list of the products. For example, when a product such as AOD is defined as many times as wavelengths, it is because the output will be wavelength dependent (a value for each wavelength):

  • NW: Number of wavelengths
  • NSD: Number of aerosol components
  • NKNOISE: Number of noises defined
  • NPARS: Number of parameters to be retrieved

The following list includes all products that can be obtained using GRASP:

  • Number of iterations (niter)
  • Total final measurement fitting residual for multi-pixel retrieval (rest)
  • Detailed (i.e., separated by the type of observation) final absolute measurement fitting residuals for segment (resat)
  • Detailed (i.e., separated by the type of observation) final relative measurement residuals for segment (resrt)
  • If SD (size distribution) retrieved in form of binned SD, the number of grid radii for SD (used to print output) (radius)
  • If SD retrieved in form of binned SD, grid radii
  • If the pre-calculated lognormal bins where used, the function describing each lognormal SD bin (used to print output) (SDL)
  • Main output for each single pixel (for both single- and multiple-pixel retrievals):
    • Single pixel total residual (meas. + smoothness constraints) (res)
    • Detailed absolute measurement residuals (resa[NKNOISE])
    • Detailed relative measurement residuals (resr[NKNOISE])
    • Retrieved aerosol and surface reflectance parameters (par[NPARS])
    • Angstrom exponent (Aexp)
    • For each wavelength:
      • Spectral total aerosol extinction (extt)
      • Spectral total aerosol single scattering albedo (ssat)
      • Spectral total aerosol absorption extinction (aext)
      If retrieved aerosol consist of several components:
      • Spectral extinction for each component (ext[NSD])
      • Spectral single scattering albedo for each component (ssa[NSD])
      • Real part of refractive index for each aerosol component (mreal[NSD])
      • Imaginary part of refractive index for each component(mimag[NSD])
  • If SD retrieved in form of binned SD, the optical properties can be calculated for fine and coarse modes, separated using chosen inflection radii:
    • Volume median radius (rv)
    • Standard deviation (std)
    • Concentration (cv)
    • Effective radius (reff)
    • Spectral extinction for each wavelength (ext[NW])
  • Phase matrix parameters:
    • Number of scattering angles (nangle)
    • Values of scattering angles (angles)
    • For each pixel and for all the wavelengths:
      • Phase matrix elements for each aerosol component ph11, ph12, ph22, ph33, ph34, ph44
      • Total aerosol phase matrix elements pht11, pht12, pht22, pht33, pht34, pht44
      • Lidar and depolarization ratios for each aerosol component (lr, dlpr)
      • Total aerosol lidar and depolarization ratios (lrt, dlprt)
    • Phase matrix norm. The analytical expresion used to calculate the norm is:

      Norm = (1/2) * ∫π0 P1,1(λ, Θ) * sin( Θ) d Θ

      However, in order to achieve an optimal degree of accuracy it is recomended to perform this integration trought the Simpson rule with 721 bins in the logaritmic space.
    • Asymmetry parameter. The analytical expresion used to calculate this magnitude is:

      < cos( Θ) > = (1/2) * ∫π0 P1,1(λ, Θ) * sin( Θ) * cos( Θ) d Θ

      It is recomended to perform the numerical calculation of the integral analogously to the case of the phase matrix norm.
  • If chemical composition is retrieved, then for each pixel and each aerosol component:
    • Relative humidity (rh(NSD))
    • Fraction of soluble (fslbl(NSD))
    • Insoluble fractions of soot (fsoot(NSD))
    • Insoluble fractions of iron (firon(NSD))
    • Insoluble fractions of “quartz” (fslbl(NSD))
    • Water fraction (fwtr(NSD))
  • Surface reflectance parameters for each pixel and for each wavelength:
    • All parameters of BRDF
    • All parameters of BPRDF
    • Surface Albedo
  • Lidar characteristics for each pixel:
    • Levels for vertical profiles
    • Vertical profiles for retrieved aerosol components (avp(NSD))
  • Lidar optical characteristics for each pixel and for each wavelength:
    • Aerosol extinction profiles (σaer(λ, h)). In order to obtain this magnitude from the variables in the standard GRASP output file, the column aerosol optical depth just needs to be weighted by the normalized aerosol vertical profile for each mode.

      σaer,i(λ, h) = τi * avpi(h)

      If two or more aerosol modes are included, the total extinction profile can be obtained adding up all profiles.

      σaer(λ, h) = Σi=1n τi * avpi(h)

    • Aerosol backscatter profiles (βaer(λ, h)). To obtain this magnitude it is only necessary to divide the aerosol extinction profile of each mode by its corresponding Lidar ratio.

      βaer(λ, h) = Σi=1n σaer,i(λ, h) / Si(λ )

    • Aerosol absorption profiles (σaerabs(λ, h)). In order to obtain this magnitude from the variables in the standard GRASP output file, the column aerosol absorption optical depth (τiabs) just needs to be weighted by the normalized aerosol vertical profile for each mode.

      σaerabs(λ, h) = Σi=1n τiabs * avpi(h)

    • SSA profiles (ω0(λ, h)). Once all the previous described magnitudes has been calculated the calculation of SSA profiles is inmidiate:

      ω0(λ, h) = σaerscat(λ, h) / σaer(λ, h) = (σaer(λ, h) - σaerabs(λ, h)) / σaer(λ, h)

    • Lidar ratio profiles (S(λ, h)):

      Saer(λ, h) = σaer(λ, h) / βaer(λ, h)

    • Phase matrix (Pj, k(λ, Θ, h)):

      Pj,k(λ, Θ, h) = Σi=1n Pj,ki(λ, Θ) σaer,i scat(λ, h) / σaerscat(λ, h)

    • Lidar depolarization profiles (δ(λ, h)):

      δ(λ,h) = (P1,1(λ, 180°, h) - P2,2(λ, 180°, h)) / (P1,1(λ, 180°, h) + P2,2(λ, 180°, h))

    • Retrieved lidar calibration coefficients (for lidar wavelength only)
  • Fit of every measured characteristic for each pixel and for each wavelength
  • Error estimation for each pixel:
    • Standard deviations of the random errors of the retrieved parameter logarithms (~relative errors) (ERRP)
    • Standard deviation of systematic errors of the retrieved parameter logarithms (BIASP)
    • Standard deviations of the random errors of the retrieved extinction for each aerosol component (~relative errors) (ERR_ext)
    • Standard deviations of systematic errors of the retrieved extinction for each aerosol component (BIAS_ext)
    • Standard deviations of the random errors of the retrieved total extinction (~relative errors) (ERR_extt)
    • Standard deviations of systematic errors of retrieved total extinction (BIAS_extt)
    • Standard deviations of the random errors of the retrieved single scattering albedo for each aerosol component (~relative errors) (ERR_ssa)
    • Standard deviations of systematic errors of retrieved single scattering albedo for each aerosol component (BIAS_ssa)
    • Standard deviations of the random errors of the of retrieved total single scattering albedo (~relative errors) (ERR_ssat)
    • Standard deviations of systematic errors of the retrieved total single scattering albedo (BIAS_ssat)
    • Standard deviations of the random errors of the of retrieved lidar ratio for each aerosol component (ERR_lr)
    • Standard deviations of systematic errors of the lidar ratio for each aerosol component (BIAS_lr)
    • Standard deviations of the random errors of the of retrieved depolarization ratio for each aerosol component (ERR_dr)
    • Standard deviations of systematic errors of the depolarization ratio for each aerosol component (BIAS_dr)
  • Radiative forcing for each pixel:
    • The heights for forcing output (HLV)
    • Broad band up-ward flux without aerosol at each height (BBUFX0)
    • Broad band down-ward flux without aerosol at each height (BBDFX0)
    • Broad band up-ward flux with aerosol at each height (BBUFXA)
    • Broad band down-ward flux with aerosol at each height (BBDFXA)
  • Estimations of aerosol particulate matter at the ground level (PM)
  • Aerosol type for each pixel (requires that the optical properties for fine and coarse modes are included in the calculated output)

4.3.2. GRASP classic output description

In this section the GRASP classic output format is going to be described for both if the output.segment.stream setting parameter has been set to "screen", in this case all output information will be printed on terminal; or alternatively if a path to an ascii file is provided. However, note that there are more possibilities of GRASP output formatting which can differ from what is going to be shown here.

The GRASP classic output is divided in three main sections:

  • Information of the residuals.

    This information is place in the head of the classic output. It contains one line per pixel with information about convergence and residuals after the last iteration. The first float number in each of these lines corresponds to the absolute value of the total inversion residual of the corresponding pixel. Then, the absolute and relative residuals for each noise type defined in settings can be found.

    Figure 4.7. An example of the residual information in GRASP classic output

                                    
          noise  abs         rel      noise  abs        rel
    0.358 1:     0.777E-03   0.262 %  2:     0.317E-03  0.115 %   pixel # 1  Residual after iteration # 21
    0.288 1:     0.124E-03   0.425 %  2:     0.341E-03  0.155 %   pixel # 2  Residual after iteration # 17
                                       
                                

    After the residuals, the information of date, time, longitude and latitude of each pixel can be found.

    The next part of the output is a list of all the retrieved parameters as they are calculated in the inversion matrix, each column is associated with a different pixel.

    Figure 4.8. An example of the vector of retrieved parameters in GRASP classic output

                                    
                                    Parameter #, Vector of retrieved parameters
                                    1    0.31000E-03    0.21000E-03
                                    2    0.38939E-02    0.89239E-02
                                    3    0.19772E-01    0.13772E-01
                                    4    0.45246E-01    0.55232E-01
                                    5    0.50813E-01    0.68823E-01
                                    6    0.35267E-01    0.23481E-01
                                    7    0.16086E-01    0.24567E-01
                                    8    0.41738E-02    0.67799E-02
                                    9    0.58516E-03    0.24514E-03
                                    10   0.44740E-04    0.44740E-04
                                    11   0.50000E-05    0.43456E-05
                                    12   0.13710E-04    0.19340E-04
                                    13   0.50331E-04    0.548370-04
                                    14   0.18693E-03    0.20019E-03
                                    15   0.67248E-03    0.23535E-03
                                    16   0.22029E-02    0.11244E-02
                                    17   0.59918E-02    0.60022E-02
                                    18   0.12137E-01    0.23731E-01
                                    19   0.17703E-01    0.19822E-01
                                    20   0.18786E-01    0.12441E-01
                                    21   0.14473E-01    0.29249E-01
                                    22   0.79026E-02    0.08708E-02
                                    23   0.29037E-02    0.33077E-02
                                    24   0.69835E-03    0.54197E-03
                                    25   0.11000E-03    0.21002E-03
                                

  • Detailed parameters.

    In the next part of the output different atmospheric and surface parameters can be found. Which of them are shown and which are not is defined by the settings parameters in retrieval.products part. Each product is identified by a header where the name and some information is provided. The different columns corresponds to the different pixels analougously to the case of the vector of retrieved parameters.

    If there are more than one atmospheric component in the retrieval the information corresponding to each mode is expressed differently depending on if it is an optical or microphysical product. In the case of microphysical products the information of the different components is expressed in different lines which start with an identificative integer.

    Figure 4.9. An example of the aerosol volume concentration in GRASP classic output for two aerosol modes

                                    
    Aerosol volume concentration (um^3/um^2 or um^3/um^3)
        1   0.47855E-01 0.68932E-01
        2   0.22771E-01 0.13793E-01
                                

    However, in the case of optical products, the corresponding mode of the product is indicated in the product header name, and each line is associated with a wavalength which is also indicated in the beginning of it.

    Figure 4.10. An example of the Aerosol Optical Depth in GRASP classic output for two aerosol modes

                                    
    Wavelength (um), AOD_Total (unitless or 1/um)
           0.44000   0.45107E+00   0.354561E+00
           0.67500   0.18630E+00   0.223511E+00
           0.87000   0.99709E-01   0.782688E-01
           1.02000   0.66909E-01   0.456723E-01
    Wavelength (um), AOD_Particle_mode_1 (unitless or 1/um)
           0.44000   0.43880E+00   0.30345E+00
           0.67500   0.17344E+00   0.20230E+00
           0.87000   0.86553E-01   0.701544-01
           1.02000   0.53769E-01   0.405612E-01
    Wavelength (um), AOD_Particle_mode_2 (unitless or 1/um)
           0.44000   0.12274E-01   0.54561E-01
           0.67500   0.12856E-01   0.23462E-01
           0.87000   0.13156E-01   0.82655E-02
           1.02000   0.13140E-01   0.65231E-02
                                

  • Information of the fitting.

    In the final part of the GRASP classic output the information is organized in different blocks. In each of these blocks the measurements associated to each wavelength and pixel of the SDATA and the fitted measurements after the final iteration can be found. Each measurement is separated by a header where the measurement in the sdata is indicated as meas_+"measurement name", and the fitted as fit_+"measurement name". If the measurement is defined by a specific geometry, this geometry is also included here. As for example the Solar Zenith Angle (sza), the zenith angle of the measurement (vis), the azimuth angle of the measurement (fis) or the corresponding scattering angle (sca_ang).

    Figure 4.11. An example of the fitting information in GRASP classic output for one wavelenght of one pixel with TOD and irradiance measurements

                                
    ------------------------------------------------------------------------
    pixel #    1   wavelength #   1     0.440 (um)
    ------------------------------------------------------------------------
           meas_tod        fit_tod
        0.70302E+00    0.70290E+00
       #      sza      vis      fis  sca_ang         meas_I          fit_I
       1    70.00   110.00   183.50     3.29    0.67436E+00    0.67369E+00
       2    70.00   110.00   184.00     3.76    0.62857E+00    0.62765E+00
       3    70.00   110.00   185.00     4.70    0.57672E+00    0.57630E+00
       4    70.00   110.00   186.00     5.64    0.55140E+00    0.55151E+00
       5    70.00   110.00   187.00     6.58    0.53713E+00    0.53745E+00
       6    70.00   110.00   188.00     7.52    0.52757E+00    0.52796E+00
       7    70.00   110.00   190.00     9.40    0.51327E+00    0.51358E+00
       8    70.00   110.00   192.00    11.27    0.50043E+00    0.50059E+00
       9    70.00   110.00   194.00    13.15    0.48706E+00    0.48712E+00
      10    70.00   110.00   196.00    15.03    0.47341E+00    0.47335E+00
      11    70.00   110.00   198.00    16.91    0.45892E+00    0.45876E+00
      12    70.00   110.00   200.00    18.78    0.44384E+00    0.44361E+00
      13    70.00   110.00   205.00    23.47    0.40481E+00    0.40449E+00
      14    70.00   110.00   210.00    28.15    0.36583E+00    0.36555E+00
      15    70.00   110.00   215.00    32.83    0.32882E+00    0.32864E+00
      16    70.00   110.00   220.00    37.49    0.29500E+00    0.29489E+00
      17    70.00   110.00   225.00    42.15    0.26492E+00    0.26483E+00
      18    70.00   110.00   230.00    46.80    0.23864E+00    0.23853E+00
      19    70.00   110.00   240.00    56.05    0.19656E+00    0.19638E+00
      20    70.00   110.00   250.00    65.23    0.16617E+00    0.16595E+00
      21    70.00   110.00   260.00    74.32    0.14471E+00    0.14453E+00
      22    70.00   110.00   270.00    83.28    0.13010E+00    0.13000E+00
      23    70.00   110.00   280.00    92.08    0.12075E+00    0.12071E+00
      24    70.00   110.00   300.00   108.94    0.11330E+00    0.11323E+00
      25    70.00   110.00   320.00   124.02    0.11489E+00    0.11461E+00
      26    70.00   110.00   340.00   135.46    0.11897E+00    0.11849E+00
                            

4.4. Forward model

GRASP has several forward models and each of them are used (or not) depending on the application. For example, to retrieve nephelometer data, just single scattering (particle properties) will be used. For other applications, GRASP has also a multiple scattering module (radiative transfer) and a lidar signal module.

4.4.1. How to use the forward model: Derived products and reprocesing data

As it was explained in Section 4.1.2, “Retrieved characteristics” section, the retrieval algorithm works iteratively over an array of parameters (in its first definition it is called the initial guess), until it represents the best solution. This solution array has same shape as the initial guess (the same parameters and defined in the same position). Once it is obtained, a final call of the forward model with the resulting array provides a complete list of output products.

For some applications, it can be useful to use the forward model without inverting any data. It can be done easily in GRASP with the use of the setting parameter retrieval.mode=forward. When no retrieval is performed, just one call of forward model is performed. If in the initial guess array the user has set an aerosol model, it will be used inside of the forward model obtaining therefore an entire output structure, with information in all fields.

This procedure can be used also to reprocess some data. If output parameters of a retrieval are stored, then they can be set as initial guess and then, running GRASP with the same settings, except for retrieval.mode=forward the entire output can be obtained again. This procedure can be used to reprocess data with many objectives such as saving storage space (just save the output array of grasp and reprocess to obtain the rest, if it is needed) or obtaining extra products in the future.

4.4.2. Synthetic data

The previous procedure can also help to simulate the input data. It is useful because a valid geometry is necessary in the input data to set an SDATA file. In this case, an aerosol model is set as an initial guess and the code works just for the forward run. Then, by using retrieval.debug.simulated_sdata_file parameter, the user can set a path to the simulated files. A SDATA file will be dump to that path where the geometry is the same and the measurements are filled with the output results of the forward run. Then, this SDATA file can be used to self-consistency tests, where synthetic data is retrieved.

4.5. Aerosol modeling in GRASP

Due to the versatility and the high degree of generalisation of GRASP, there are different approaches to model aerosols in order to maximize the possibilities of the different retrieval schemes. Each of the approaches described below presents different advantages depending on the information available in the input measurements, the desired output products or the available computation time.

4.5.1. Kernels

In this approach the retrieved characteristics which determine the optical properties of aerosols are the size distribution, the real and imaginary refractive indexes, and the percentage of spherical particles (the non-spherical part is modeled as spheroids). The calculation of the radiative properties (phase matrix, scattering and absorption cross-sections) from these initial characteristics is made through four-dimensional look-up-tables called "Kernels".

The four dimensions of the Kernels are: the size parameter, sphericity and the real and the imaginary refractive indices. Due to the extensive nature of these kernels, this approach represents the less restricted methodology, because any possible combination of the represented characteristics can be obtained as the solution.

The lack of intrinsic restrictions of this approach (note that as in all GRASP retrievals, apriori constraints can be set for any characteristic) presents obvious advantages. However, the general nature of Kernels methodology normally requires input measurements containing a higher amount of information in comparison with other approaches.

There are three different ways to represent the aerosol size distribution, from the most general to the simplest: triangle bins, lognormal bins and precalculated lognormal bins. The corresponding characteristics names in the GRASP YAML settings file are: “size_distribution_triangle_bins”, “size_distribution_lognormal” and “size_distribution_precalculated_lognormal”. Independently of the selected representation of the aerosol size distribution, the rest of the steps of this methodology to obtain aerosol radiative properties are the same.

In the triangle bins representation the value of each bin is retrieved independently from the others, and it corresponds to the integrated value following the trapezoidal rule between the provided size limits. The lognormal bins approach is a simplified version of the former, where each aerosol size distribution mode is modeled as a gaussian function only represented by three parameters: the center, the standard deviation and the norm. The former value is the aerosol concentration of the corresponding mode. The precalculated lognormal bins size distribution represents a step further in the simplification of this characteristic. In this case each aerosol mode is also represented by a Gaussian function. However, the center and the standard deviation are fixed to a predefined value and only the norm (aerosol concentration) of each mode is retrieved.

4.5.2. Models

The main retrieved aerosol products of this approach are the fractions of the total aerosol concentration of precalculated aerosol models. These precalculated models correspond to the main aerosol types established by all our previous experiences (Ex: smoke, urban, oceanic and dust). Each of these models correspond to a fixed particle size distribution and refractive indices, containing the already calculated phase matrix, and the extinction and absorption cross-sections. The total aerosol retrieved characteristics can be obtained by weighting the characteristics used to calculate each of these models by its corresponding volume concentration.

The significant reduction of retrieved parameters makes this approach very suitable for the retrieval of input measurements with a reduced amount of information. The absence of Kernels in the whole process makes this option by far the fastest of all. One of the main drawbacks of this methodology is the fact that the inversion is intrinsically constrained by the selected models. However, these models have been carefully selected to be suitable to cover almost all atmospheric situations. Moreover, they can be recalculated or extended in order to cover specific situations.

An example of the necessary settings to use this approach can be found below. Where in the phase matrix section the bin of each mode no longer represents the limit radius but the accounted aerosol models. Particle size distribution, refractive indexes and sphericity characteristics are substituted by another characteristic called "aerosol_model_concentration".


forward_model:

      phase_matrix:
          size_binning_method_for_triangle_bins: logarithm
          number_of_elements: 4
          kernels_folder: "models_ang35_wl22_optimized"
          radius:
              mode[1]:
                  bins: [ 1.,  2., 3., 4., 5. ]
                  
  .
  .
  .


  constraints:
      characteristic[1]:
          type: aerosol_model_concentration
          retrieved: true
          mode[1]:
              initial_guess:                       #smoke       #urban    #oceanic   #dust
                  value:                          [0.9,           0.4,      0.4,     0.4.,      0.4]
                  min:                            [0.000005, 0.000005, 0.000005, 0.000005, 0.000005]
                  max:                            [1.0,           1.0,      1.0,      1.0,      1.0]
                  index_of_wavelength_involved:   [0,               0,        0,        0,        0]
              single_pixel:
              .
              .
              .

      characteristic[2]:
          type: aerosol_concentration
          retrieved: true
          mode[1]:
              initial_guess:
                  value:                          [0.05   ]
                  min:                            [0.0001 ]
                  max:                            [5.0    ]
                  index_of_wavelength_involved:   [0      ]
              single_pixel:
              .
              .
              .

        

4.5.3. Chemistry

The chemistry approach can be seen somehow as an intermediate approach between Kernels and Models. In this case size distribution and sphericity percentage are directly retrieved as in the case of Kernels. However, instead of the refractive indexes, here the fractions of the total aerosol concentration corresponding to the different chemical components are retrieved. The refractive indices corresponding to each of these components are predefined. Thus, the weighted refractive indexes of all the fractions in combination with the particle size distribution and the sphericity parameter are used as input for the Kernels look-up-tables to calculate the radiative properties.

In comparison with the Models approach, where the radiative properties can only be linearly weighted, the refractive index presents an extra layer of complexity when they are mixed. Because the weighting methodology used to obtain the total refractive index retains information about the aerosol internal structure. Two mixing possibilities are available now in GRASP: an internal mixture, where a linear mixture of the different components is performed; but also a Maxwell-Garnett mixture is available, where one element is considered as the main host and the rest of the chemical elements are taken as inclusions inside of it.

In this case the necessary settings to use this approach are very similar to Kernels. The size distribution can be represented using the three already described possibilities in the Kernels approach. However, the refractive index characteristics are substituted by "particle_component_volume_fractions_linear_mixture" in the case of the linear mixture, or by "particle_component_fractions_chemical_mixture" in the case of a Maxwell-Garnett mixture. An extra section called "chemistry" has to be added in the Forward model part, where it is provided: the names of each chemical component accounted in the inversion, the soluble component (if necessary) and the path to the directory where these predefined refractive index look-up-tables are located.

    
    forward_model:

        aerosol:
            chemistry:
                folder: "chemistry_refractive_indexes/"
                soluble: "ammnm_ntrt"
                species:
                    mode[1]: [ 'black_carbon', 'mix_dust','iron_oxide','water']
                    
        phase_matrix:
            size_binning_method_for_triangle_bins: logarithm
            number_of_elements: 1
            kernels_folder: "KERNELS_BASE/"
            radius:
                mode[1]:
                    min: 0.05
                    max: 15.0
                    
    .
    .
    .


    constraints:
        characteristic[1]:
            type: size_distribution_triangle_bins
            retrieved: true
            mode[1]:
                initial_guess:
                    value:                        [1.4715e-05, 1.4715e-03, ..., 1.4715e-03, 1.4715e-05]
                    min:                          [0.00001,  0.00001,  ...,  0.00001,  0.00001]
                    max:                          [1.0,      1.0,      ...,      1.0,      1.0]
                    index_of_wavelength_involved: [0,        0,       ...,        0,        0]
            .
            .
            .
                        
        characteristic[2]:
            type: particle_component_volume_fractions_linear_mixture
            retrieved: true
            mode[1]:
                initial_guess:             #1         #2        #3         #4
                    value:                 [0.00001,       0.0001,      0.0001,      0.0001      ]
                    min:                   [0.000001,    0.000001,   0.000001,    0.000001   ]
                    max:                   [0.2,        1.0,       1.0,        1.0      ]
                    index_of_wavelength_involved:   [0,        0,      0,     0      ]
                .
                .
                .

        
            

4.5.4. Transport models

The transport models approach constitutes the most complex methodology to model atmospheric aerosols. This approach has been designed to facilitate the interface between sophisticated aerosol transport models and GRASP code, both for forward modeling and retrieval. Normally, these complex transport models work simultaneously with a high number of different aerosol particles with specific size distributions, shape, optical characteristics, vertical distributions… In the rest of the GRASP approaches for aerosol modeling only one or two aerosol modes are considered which have all these aforementioned characteristics totally independent between each other. However, the GRASP transport model approach is not only an additional interface that enables the possibility to work with multiple independent aerosol modes. But it is also a tool that allows the conversion between the characteristics that define the particles in the aerosol transport models (as mass mixing ratio) to the normal microphysical and optical characteristics that are normally used in GRASP.

GRASP transport model approach consists of 5 main aerosol components similar to MERRA-2 and CAMS aerosol models: sulphate (SU), desert dust (DU), sea salt (SS), organic (OC) and black carbon (BC). Each component can exhibit hydrophobic or hydrophilic properties resulting in 15 aerosol tracers: hydrophilic sulphate (SU), five size bins for hydrophobic dust (DU) and hydrophilic sea salt (SS), hydrophobic and hydrophilic modes of organic (OC) and black carbon (BC) aerosol . The concentration and optical depth of each tracer in the model is refined through the mass mixing ratio. Furthermore, each of these aerosol modes count with its own vertical profile, particle size distribution, sphericity parameter and all the rest of optical characteristics. Note that the different aerosol tracers are externally mixed with vertically dependent mass mixing ratios to obtain the corresponding total aerosol optical properties.

An example of the necessary settings to define the additional aerosol parameters of the transport models can be found below:

                
                forward_model:
                    phase_matrix:            
                        size_binning_method_for_triangle_bins: logarithm
                        number_of_elements: 4 
                        use_transport_model: true
                        number_of_bins_for_lognormal_size_distribution: [ ^repeat(43;15) ]
                        transport_model: 
                            vertical_profile: column_average  #tracer_average
                            tracers:     ['du1', 'du2', 'du3', 'du4', 'du5','bc1','bc2','oc1','oc2','ss1', 'ss2', 'ss3', 'ss4', 'ss5','su1'  ]
                            hydrophilic: [  0,     0,     0,     0,     0,    0,    1,    0,    1,    1,     1,     1,     1,     1,    1   ]
                            #hydrophilic: [  0,     0,     0,     0,     0,    0,    0,    0,    0,    0,     0,     0,     0,     0,    0   ]
                            density:     [2500., 2650., 2650., 2650., 2650.,1000.,1000.,1800.,1800., 2200., 2200., 2200., 2200.,2200., 1700.]                            
                        kernels_folder:  'KERNELS_BASE' # 'merra2_data_base'
                        radius:
                            mode[1]: # du1
                                min: 0.08
                                max: 20.0
                            mode[2]: # du2
                                min: 0.08
                                max: 20.0
                            mode[3]: # du3
                                min: 0.08
                                max: 20.0
                            mode[4]: # du4
                                min: 0.08
                                max: 20.0
                            mode[5]: # du5
                                min: 0.08
                                max: 20.0
                            mode[6]: # bc1
                                min: 0.0001 #0.002
                                max: 5.0
                            mode[7]: # bc2
                                min: 0.0001 #0.002
                                max: 5.0
                            .
                            .
                            . 
            

Where “forward_model.phase_matrix.transport_model.tracers” refers to the 15 aerosol tracers. “forward_model.phase_matrix.transport_model.density” and “forward_model.phase_matrix.transport_model.hyfrophilic” mark correspondingly the density of dry tracers and their hygroscopic properties (“0” corresponds to hydrophobic and “1” corresponds to hydrophilic aerosol).

The rest of the aerosol characteristics are defined following the standard GRASP settings convention for different aerosol modes:

                
                constraints:
                    characteristic[1]: #1
                        type: size_distribution_lognormal
                        retrieved: true
                        mode[1]: 
                            initial_guess: 
                                value:                        [0.6576, 0.2828]
                                min:                          [0.1,   0.1]
                                max:                          [3.0,   0.9]
                                index_of_wavelength_involved: [^repeat(0;2)]
                            single_pixel:
                                a_priori_estimates:
                                    lagrange_multiplier:      [1.0e-03,  1.0e-03]
                                smoothness_constraints:
                                    difference_order: 0
                                    lagrange_multiplier: 0.0
                        mode[2]: 
                            initial_guess: 
                                value:                        [1.2436, 0.25]
                                min:                          [0.1,   0.1]
                                max:                          [3.0,   0.9]
                                index_of_wavelength_involved: [^repeat(0;2)]
                            single_pixel:
                                a_priori_estimates:
                                    lagrange_multiplier:      [1.0e-03,  1.0e-03]
                                smoothness_constraints:
                                    difference_order: 0
                                    lagrange_multiplier: 0.0
                        mode[3]: 
                            initial_guess: 
                                value:                        [2.2969, 0.344]
                                min:                          [0.1,   0.1]
                                max:                          [3.0,   0.9]
                                index_of_wavelength_involved: [^repeat(0;2)]
                            single_pixel:
                                a_priori_estimates:
                                    lagrange_multiplier:      [1.0e-03,  1.0e-03]
                                smoothness_constraints:
                                    difference_order: 0
                                    lagrange_multiplier: 0.0
                        mode[4]: 
                            initial_guess: 
                                value:                        [5.3645, 0.707]
                                min:                          [0.1,   0.1]
                                max:                          [6.0,   0.9]
                                index_of_wavelength_involved: [^repeat(0;2)]
                            single_pixel:
                                a_priori_estimates:
                                    lagrange_multiplier:      [1.0e-03,  1.0e-03]
                                smoothness_constraints:
                                    difference_order: 0
                                    lagrange_multiplier: 0.0
                        .
                        .
                        .
                

                    characteristic[2]: 
                                type: tracer_level_concentration
                                retrieved: true
                                mode[1]: 
                                    initial_guess: 
                                        value:                        [^repeat(0.007;72)]
                                        min:                          [^repeat(1e-15;72)]
                                        max:                          [^repeat(0.500;72)]
                                        index_of_wavelength_involved: [^repeat(0;72)]
                                    single_pixel:
                                        a_priori_estimates:
                                            lagrange_multiplier:      [^repeat(0;72)]
                                        smoothness_constraints:
                                            difference_order: 0
                                            lagrange_multiplier: 0.0
                                mode[2]: 
                                    initial_guess: 
                                        value:                        [^repeat(0.007;72)]
                                        min:                          [^repeat(1e-15;72)]
                                        max:                          [^repeat(0.500;72)]
                                        index_of_wavelength_involved: [^repeat(0;72)]
                                    single_pixel:
                                        a_priori_estimates:
                                            lagrange_multiplier:      [^repeat(0;72)]
                                        smoothness_constraints:
                                            difference_order: 0
                                            lagrange_multiplier: 0.0
                                mode[3]: 
                                    initial_guess: 
                                        value:                        [^repeat(0.007;72)]
                                        min:                          [^repeat(1e-15;72)]
                                        max:                          [^repeat(0.500;72)]
                                        index_of_wavelength_involved: [^repeat(0;72)]
                                    single_pixel:
                                        a_priori_estimates:
                                            lagrange_multiplier:      [^repeat(0;72)]
                                        smoothness_constraints:
                                            difference_order: 0
                                            lagrange_multiplier: 0.0
                                mode[4]: 
                                    initial_guess: 
                                        value:                        [^repeat(0.007;72)]
                                        min:                          [^repeat(1e-15;72)]
                                        max:                          [^repeat(0.500;72)]
                                        index_of_wavelength_involved: [^repeat(0;72)]
                                    single_pixel:
                                        a_priori_estimates:
                                            lagrange_multiplier:      [^repeat(0;72)]
                                        smoothness_constraints:
                                            difference_order: 0
                                            lagrange_multiplier: 0.0

                    characteristic[3]: 
                                type: real_part_of_refractive_index_spectral_dependent
                                retrieved: true
                                mode[1]:
                                    initial_guess:                       #1         #2        #3        #4        #5        #6       
                                        value:                        [^repeat(1.53;7)]
                                        min:                          [^repeat(1.33;7)]
                                        max:                          [^repeat(1.6;7)]
                                        index_of_wavelength_involved: [^repeat(0;7)]
                                    single_pixel:
                                        smoothness_constraints:
                                            difference_order: 1
                                            lagrange_multiplier: 1.0e+1
                                    multi_pixel:
                                        smoothness_constraints:
                                            derivative_order_of_X_variability:    1
                                            lagrange_multiplier_of_X_variability: 1.0e-1                        
                                            derivative_order_of_Y_variability:    1
                                            lagrange_multiplier_of_Y_variability: 1.0e-1                     
                                            derivative_order_of_T_variability:    1
                                            lagrange_multiplier_of_T_variability: 2.0e-2                                               
                                mode[2]:
                                    initial_guess:                       #1         #2        #3        #4        #5        #6       
                                        value:                        [^repeat(1.53;7)]
                                        min:                          [^repeat(1.33;7)]
                                        max:                          [^repeat(1.6;7)]
                                        index_of_wavelength_involved: [^repeat(0;7)]
                                    single_pixel:
                                        smoothness_constraints:
                                            difference_order: 1
                                            lagrange_multiplier: 1.0e+1
                                    multi_pixel:
                                        smoothness_constraints:
                                            derivative_order_of_X_variability:    1
                                            lagrange_multiplier_of_X_variability: 1.0e-1                        
                                            derivative_order_of_Y_variability:    1
                                            lagrange_multiplier_of_Y_variability: 1.0e-1                     
                                            derivative_order_of_T_variability:    1
                                            lagrange_multiplier_of_T_variability: 2.0e-2                                               
                                mode[3]:
                                    initial_guess:                       #1         #2        #3        #4        #5        #6       
                                        value:                        [^repeat(1.53;7)]
                                        min:                          [^repeat(1.33;7)]
                                        max:                          [^repeat(1.6;7)]
                                        index_of_wavelength_involved: [^repeat(0;7)]
                                    single_pixel:
                                        smoothness_constraints:
                                            difference_order: 1
                                            lagrange_multiplier: 1.0e+1
                                    multi_pixel:
                                        smoothness_constraints:
                                            derivative_order_of_X_variability:    1
                                            lagrange_multiplier_of_X_variability: 1.0e-1                        
                                            derivative_order_of_Y_variability:    1
                                            lagrange_multiplier_of_Y_variability: 1.0e-1                     
                                            derivative_order_of_T_variability:    1
                                            lagrange_multiplier_of_T_variability: 2.0e-2       

                                .
                                .
                                .
            

4.6. Error estimation

GRASP provides rigorous dynamic error estimates of the retrieved characteristics, but also these calculations are available for some of the derived parameters that GRASP calculates internally.

In order to perform the calculations some extra information for each measurement has to be provided in the GRASP settings file. Two different kinds of bias can be assumed: “bias_measurements_synthetic” and “bias_equation”. The units of these bias are the same as the corresponding measurement and the “error_type” option establishes if they are taken as an absolute or relative value to the measurement. An example of the definition of these settings can be found below:

                
                noises:
                    noise[1]:
                        standard_deviation_synthetic: 0.05
                        bias_measurements_synthetic: 0.05
                        bias_equation: 0.05
                        error_type:  relative
                        standard_deviation:  0.03
                        measurement_type[1]:
                            type: I
                            index_of_wavelength_involved: [ 1, 2, 3, 4 ]
                    noise[2]:
                        standard_deviation_synthetic: 0.01
                        bias_measurements_synthetic: 0.01
                        bias_equation: 0.01
                        error_type:  absolute
                        standard_deviation:  0.01
                        measurement_type[1]:
                            type: aod
                            index_of_wavelength_involved: [ 1, 2, 3, 4 ]

        

In a similar way to other products provided by GRASP, in the “products” section the “error_estimation” group it is possible to configure how these error estimates are made with the option “using_Levenberg-Marquardt”. In order to select to what magnitudes are applied, the setting “retrieved” establishes to provide the error estimates for all retrieved characteristics, and in the “derived” group the user can select to what of the derived products (ex.: AOD, Angstrom Exponent, SSA…) these calculations will be applied:

                
                products:  
                    error_estimation:
                        using_Levenberg-Marquardt: true
                        derived:
                                aerosol:
                                    lidar: true
                                    optical_properties: true
                        retrieved: true 
        

The error estimation of the selected magnitudes, retrieved and derived, can be found in the end of the classic GRASP output file separated in three: “Total standard deviations”, “Standard deviations” and “BIAS - Standard deviation”. An example of an output of the error estimates can be seen below:

                
                -------------------------------------------------------------------------
                Total standard deviations of retrieved parameter logarithms (~relative errors) :
                -------------------------------------------------------------------------
                Date:             2014-08-22
                Time:               14:58:12
                            1   0.84496E+00
                            2   0.59295E+00
                            3   0.32039E+00
                            4   0.12886E+00
                            5   0.10403E+00
                    .
                    .
                    .
                -------------------------------------------------------------------------
                Standard deviations of retrieved parameter logarithms (~relative errors) :
                -------------------------------------------------------------------------
                Date:             2014-08-22
                Time:               14:58:12
                            1   0.84343E+00
                            2   0.35544E+00
                            3   0.15784E+00
                            4   0.11046E+00
                            5   0.84152E-01
                        .
                    .
                    .


                ----------------------------------------------------------------------------------
                BIAS - Standard deviation of systematic errors of retrieved parameter logarithms :
                ----------------------------------------------------------------------------------
                Date:             2014-08-22
                Time:               14:58:12
                            1   0.50815E-01
                            2   0.47461E+00
                            3   0.27880E+00
                            4   0.66360E-01
                            5  -0.61162E-01
                                .
                    .
                    .

                INVSING = 0
                --------------------------------------------------------------------------------------------
                Total standard deviations of retrieved optical characteristic logarithms (~relative errors) :
                --------------------------------------------------------------------------------------------
                Date:             2014-08-22
                Time:               14:58:12
                Wavelength (um), Aerosol Optical Depth (Random) for Particle component 1
                        0.4400   0.17618E-01
                        0.6750   0.68020E-01
                        0.8700   0.86484E-01
                        1.0200   0.82026E-01
                Wavelength (um), Single Scattering Albedo (Random) for Particle component 1
                        0.4400   0.59979E-01
                        0.6750   0.49679E-01
                        0.8700   0.25789E-01
                        1.0200   0.22196E-01
                --------------------------------------------------------------------------------------
                Standard deviations of retrieved optical characteristic logarithms (~relative errors) :
                --------------------------------------------------------------------------------------
                Date:             2014-08-22
                Time:               14:58:12
                Wavelength (um), Aerosol Optical Depth (Random) for Particle component   1
                        0.4400   0.65241E-02
                        0.6750   0.13940E-01
                        0.8700   0.19455E-01
                        1.0200   0.25821E-01
                Wavelength (um), Single Scattering Albedo (Random) for Particle component   1
                        0.4400   0.12061E-01
                        0.6750   0.99402E-02
                        0.8700   0.13855E-01
                        1.0200   0.18953E-01
                -----------------------------------------------------------------------------------------------
                BIAS - Standard deviations of systematic errors of retrieved optical characteristic logarithms :
                -----------------------------------------------------------------------------------------------
                Date:             2014-08-22
                Time:               14:58:12
                Wavelength (um), Aerosol Optical Depth (Bias) for Particle mode 1
                        0.4400   0.16370E-01
                        0.6750   0.66575E-01
                        0.8700   0.84366E-01
                        1.0200   0.77988E-01
                Wavelength (um), Single Scattering Albedo (Bias) for Particle mode 1
                        0.4400   0.58718E-01
                        0.6750   0.48694E-01
                        0.8700   0.21934E-01
                        1.0200   0.13280E-01
        

In order to make a proper interpretation of the error estimates provided by GRASP two considerations have to be considered. First, the random (standard deviation) and the bias components for each parameter have to be added quadratically to obtain the total error:

Secondly, GRASP operates in the logarithmic space, which includes de error estimates. Thus, some calculations are needed in order to represent together the value of the different magnitudes and the corresponding error estimates:

ln(a*) = ln(a) ± σa

Thus:

a exp( σa) = a*high

a exp( - σa) = a*low