CR2RE Pipeline Reference Manual 1.6.10
Modules | Functions
Resample 2D/3D images/cubes

Algorithms to resample 2D images and 3D cubes onto a common grid based on the MUSE code. More...

Modules

 HDRL parameter - interpolation method
 HDRL parameter controlling the interpolation method.
 
 HDRL parameter - output grid
 Hdrl parameter defining the final output grid.
 

Functions

hdrl_resample_result * hdrl_resample_compute (const cpl_table *ResTable, hdrl_parameter *method, hdrl_parameter *outputgrid, const cpl_wcs *wcs)
 High level resampling function.
 
void hdrl_resample_result_delete (hdrl_resample_result *aCube)
 Deallocates the memory associated to a hdrl_resample_result object.
 
cpl_table * hdrl_resample_image_to_table (const hdrl_image *hima, const cpl_wcs *wcs)
 Convert a hdrl image into a cpl table that can be given as input to hdrl_resample_compute()
 
cpl_table * hdrl_resample_imagelist_to_table (const hdrl_imagelist *himlist, const cpl_wcs *wcs)
 Convert a hdrl imagelist into a cpl table that can be given as input to hdrl_resample_compute().
 

Detailed Description

Algorithms to resample 2D images and 3D cubes onto a common grid based on the MUSE code.

The routine currently implements the following algorithms:

The 2D and 3D interpolation is done in 2 dimension and 3 dimension, respectively. Moreover, additional error based weights can be taken into account.

The calculation is performed by calling the top-level function hdrl_resample_compute(). The latter does not directly work on images but on a cpl table. The table is created from a 2D image by calling the function hdrl_resample_image_to_table() or from a 3D cube by calling the function hdrl_resample_imagelist_to_table(). The advantage of this is that the user can combine many images/cubes into a single table and perform the interpolation based on all information in one step.

Beside the world coordinate system (wcs), there are two hdrl_parameter passed to the function, one controlling the interpolation method and one defining the output grid.

The parameter controlling the interpolation method can be created by

that sets the interpolation method one would like to use.

The parameter controlling the output grid can be created by

depending on a 2D or 3D interpolation and if the user wants to specify the final output grid on a high granularity or the hdrl routine should derive it from the input data.

Function Documentation

◆ hdrl_resample_compute()

hdrl_resample_result * hdrl_resample_compute ( const cpl_table *  ResTable,
hdrl_parameter *  method,
hdrl_parameter *  outputgrid,
const cpl_wcs *  wcs 
)

High level resampling function.

Parameters
ResTableThe cpl table to be resampled. Should be derived by the hdrl_resample_imagelist_to_table() or by the hdrl_resample_image_to_table() function.
methodThe hdrl parameter defining the resampling method
outputgridThe hdrl parameter defining the output grid
wcsThe cpl wcs parameter do derive scalings/normalisations
Returns
The hdrl_resample_result structure containing all informations of the resampled output (data, error, bpm, wcs encoded in a cpl_propertylist). The returned hdrl_resample_result structure must be deleted with hdrl_resample_result_delete()

The routine currently implements the following algorithms:

  • Nearest neighbor resampling
  • Weighted resampling using Renka weighting function
  • Weighted resampling using inverse distance weighting function
  • Weighted resampling using quadratic inverse distance weighting function
  • Weighted resampling using a drizzle-like weighting scheme
  • Weighted resampling using a lanczos-like restricted sinc for weighting

The 2D and 3D interpolation is done in 2 dimension and 3 dimension, respectively. Moreover, additional error based weights can be taken into account.

As can be seen from the interface, this function does not directly work on images but on the cpl table ResTable. The table is created from a 2D image by calling the function hdrl_resample_image_to_table() or from a 3D cube by calling the function hdrl_resample_imagelist_to_table(). In case you create the table on your own without the above mentioned functions, make sure to flag all(!) non normal pixels as bad (e.g. by using the isfinite() function).

The parameter controlling the interpolation method (method) is created with

depending on the interpolation method one would like to use.

The parameter controlling the output grid (outputgrid) is created with

depending on a 2D or 3D interpolation and if the user wants to specify the final output grid on a high granularity or the hdrl routine should derive it from the input data.

Definition at line 735 of file hdrl_resample.c.

References hdrl_imagelist_get_size(), hdrl_resample_parameter_method_verify(), and hdrl_resample_parameter_outgrid_verify().

◆ hdrl_resample_image_to_table()

cpl_table * hdrl_resample_image_to_table ( const hdrl_image *  hima,
const cpl_wcs *  wcs 
)

Convert a hdrl image into a cpl table that can be given as input to hdrl_resample_compute()

Parameters
himahdrl_image containing the data/error/bpm
wcsThe world coordinate system used for the conversion
Returns
cpl_table that is used by the hdrl_resample_compute() function

Definition at line 2582 of file hdrl_resample.c.

References hdrl_imagelist_delete(), hdrl_imagelist_new(), hdrl_imagelist_set(), hdrl_imagelist_unset(), and hdrl_resample_imagelist_to_table().

◆ hdrl_resample_imagelist_to_table()

cpl_table * hdrl_resample_imagelist_to_table ( const hdrl_imagelist *  himlist,
const cpl_wcs *  wcs 
)

Convert a hdrl imagelist into a cpl table that can be given as input to hdrl_resample_compute().

Parameters
himlisthdrl_imagelist containing the data/error/bpm
wcsThe world coordinate system used for the conversion
Returns
cpl_table that is used by the hdrl_resample_compute() function

Definition at line 2612 of file hdrl_resample.c.

References hdrl_image_get_error_const(), hdrl_image_get_image_const(), hdrl_image_get_mask_const(), hdrl_imagelist_get_const(), hdrl_imagelist_get_size(), hdrl_imagelist_get_size_x(), and hdrl_imagelist_get_size_y().

Referenced by hdrl_resample_image_to_table().

◆ hdrl_resample_result_delete()

void hdrl_resample_result_delete ( hdrl_resample_result *  aCube)

Deallocates the memory associated to a hdrl_resample_result object.

Parameters
aCubeinput hdrl_resample_result object
Returns
void

Definition at line 849 of file hdrl_resample.c.

References hdrl_imagelist_delete().