ARTS  2.4.0(git:4fb77825)
tmatrix.h File Reference

Declarations for the T-Matrix interface. More...

#include "messages.h"
#include "optproperties.h"

Go to the source code of this file.

Functions

void calcSingleScatteringDataProperties (SingleScatteringData &ssd, ConstMatrixView ref_index_real, ConstMatrixView ref_index_imag, const Numeric equiv_radius=200, const Index np=-1, const Numeric axial_ratio=1.000001, const Numeric precision=0.001, const Index ndgs=2, const Index robust=0, const Index quiet=1)
 Calculate SingleScatteringData properties. More...
 
void tmatrix_ampld_test (const Verbosity &verbosity)
 T-Matrix validation test. More...
 
void tmatrix_tmd_test (const Verbosity &verbosity)
 T-Matrix validation test. More...
 
void calc_ssp_random_test (const Verbosity &verbosity)
 Single scattering properties calculation for randomly oriented particles. More...
 
void calc_ssp_fixed_test (const Verbosity &verbosity)
 Single scattering properties calculation for particles with fixed orientation. More...
 

Detailed Description

Declarations for the T-Matrix interface.

Author
Oliver Lemke
Date
2013-06-25

Definition in file tmatrix.h.

Function Documentation

◆ calc_ssp_fixed_test()

void calc_ssp_fixed_test ( const Verbosity verbosity)

Single scattering properties calculation for particles with fixed orientation.

Two cases are calculated. One with oblate particles which is equivalent to the following PyARTS case:

from PyARTS import arts_types
from PyARTS import constants
params = {'ptype': constants.PTYPE_AZIMUTH_RND,
          'f_grid': [230e9, 240e9],
          'T_grid': [220, 250],
          'za_grid': numpy.arange(0, 181, 10),
          'aa_grid': numpy.arange(0, 181, 10),
          'equiv_radius': 200, # equivalent volume radius
          'NP':-1, # -1 for spheroid, -2 for cylinder, positive for chebyshev
          'phase':'ice',
          'mrr': numpy.array([[1.78031135, 1.78150475], [1.78037238, 1.78147686]]),
          'mri': numpy.array([[0.00278706, 0.00507565], [0.00287245, 0.00523012]]),
          'axial_ratio': 1.5}
s = arts_types.SingleScatteringData(params)
s.calc()
 

And one with prolate particles which is equivalent to the following PyARTS case:

from PyARTS import arts_types
from PyARTS import constants
params = {'ptype': constants.PTYPE_AZIMUTH_RND,
          'f_grid': [230e9, 240e9],
          'T_grid': [220, 250],
          'za_grid': numpy.arange(0, 181, 10),
          'aa_grid': numpy.arange(0, 181, 10),
          'equiv_radius': 200, # equivalent volume radius
          'NP':-1, # -1 for spheroid, -2 for cylinder, positive for chebyshev
          'phase':'ice',
          'mrr': numpy.array([[1.78031135, 1.78150475], [1.78037238, 1.78147686]]),
          'mri': numpy.array([[0.00278706, 0.00507565], [0.00287245, 0.00523012]]),
          'axial_ratio': 0.7}
s = arts_types.SingleScatteringData(params)
s.calc()
 
Author
Oliver Lemke

Definition at line 1507 of file tmatrix.cc.

References SingleScatteringData::aa_grid, calcSingleScatteringDataProperties(), CREATE_OUT0, SingleScatteringData::f_grid, ConstVectorView::nelem(), nlinspace(), SingleScatteringData::ptype, PTYPE_AZIMUTH_RND, SingleScatteringData::T_grid, and SingleScatteringData::za_grid.

Referenced by TMatrixTest().

◆ calc_ssp_random_test()

void calc_ssp_random_test ( const Verbosity verbosity)

Single scattering properties calculation for randomly oriented particles.

Two cases are calculated. One with oblate particles which is equivalent to the following PyARTS case:

from PyARTS import arts_types
params = {'ptype': constants.PTYPE_AZIMUTH_RND,
          'f_grid': [230e9, 240e9],
          'T_grid': [220, 250],
          'za_grid': numpy.arange(0, 181, 10),
          'aa_grid': numpy.arange(0, 181, 10),
          'equiv_radius': 200, # equivalent volume radius
          'NP':-1, # -1 for spheroid, -2 for cylinder, positive for chebyshev
          'phase':'ice',
          'mrr': numpy.array([[1.78031135, 1.78150475], [1.78037238, 1.78147686]]),
          'mri': numpy.array([[0.00278706, 0.00507565], [0.00287245, 0.00523012]]),
          'axial_ratio': 1.5}
s = arts_types.SingleScatteringData(params)
s.calc()
 

And one with prolate particles which is equivalent to the following PyARTS case:

from PyARTS import arts_types
params = {'ptype': constants.PTYPE_AZIMUTH_RND,
          'f_grid': [230e9, 240e9],
          'T_grid': [220, 250],
          'za_grid': numpy.arange(0, 181, 10),
          'aa_grid': numpy.arange(0, 181, 10),
          'equiv_radius': 200, # equivalent volume radius
          'NP':-1, # -1 for spheroid, -2 for cylinder, positive for chebyshev
          'phase':'ice',
          'mrr': numpy.array([[1.78031135, 1.78150475], [1.78037238, 1.78147686]]),
          'mri': numpy.array([[0.00278706, 0.00507565], [0.00287245, 0.00523012]]),
          'axial_ratio': 0.7}
s = arts_types.SingleScatteringData(params)
s.calc()
 
Author
Oliver Lemke

Definition at line 1454 of file tmatrix.cc.

References SingleScatteringData::aa_grid, calcSingleScatteringDataProperties(), CREATE_OUT0, SingleScatteringData::f_grid, ConstVectorView::nelem(), nlinspace(), SingleScatteringData::ptype, PTYPE_TOTAL_RND, SingleScatteringData::T_grid, and SingleScatteringData::za_grid.

Referenced by TMatrixTest().

◆ calcSingleScatteringDataProperties()

void calcSingleScatteringDataProperties ( SingleScatteringData ssd,
ConstMatrixView  ref_index_real,
ConstMatrixView  ref_index_imag,
const Numeric  equiv_radius = 200,
const Index  np = -1,
const Numeric  axial_ratio = 1.000001,
const Numeric  precision = 0.001,
const Index  ndgs = 2,
const Index  robust = 0,
const Index  quiet = 1 
)

Calculate SingleScatteringData properties.

Port of calc_SSP function from PyARTS.

Parameters
[in,out]ssdAs input all grids must be set (f, T, za and aa). These grids are given by ssd are used to calculate pha_mat_data, ext_mat_data and abs_vec_data. Also ssd.ptype must be set as input. The output ssd has remaining fields set.
[in]ref_index_realVector with real parts of refractive index Number of rows must match elements in ssd.f_grid Number of cols must match elements in ssd.T_grid
[in]ref_index_imagVector with imaginary parts of refractive index
[in]equiv_radiusequivalent volume radius [micrometer]
[in]npParticle type (-1 for spheroid, -2 for cylinder)
[in]axial_ratioAxial ratio of particles
[in]precisionAccuracy of the computations
[in]ndgsThe number of division points in computing integrals over the particle surface.
Author
Oliver Lemke

Definition at line 979 of file tmatrix.cc.

References SingleScatteringData::abs_vec_data, SingleScatteringData::ext_mat_data, SingleScatteringData::f_grid, ARTS::Var::f_index(), ConstMatrixView::ncols(), ConstVectorView::nelem(), ConstMatrixView::nrows(), SingleScatteringData::pha_mat_data, PI, precision, SingleScatteringData::ptype, PTYPE_TOTAL_RND, Tensor7::resize(), Tensor5::resize(), SPEED_OF_LIGHT, SingleScatteringData::T_grid, tmatrix_random_orientation(), and SingleScatteringData::za_grid.

Referenced by calc_ssp_fixed_test(), and calc_ssp_random_test().

◆ tmatrix_ampld_test()

void tmatrix_ampld_test ( const Verbosity verbosity)

T-Matrix validation test.

Executes the standard test included with the double precision T-Matrix code for nonspherical particles in a fixed orientation. Should give the same as running the 3rdparty/tmatrix/tmatrix_ampld executable.

Author
Oliver Lemke

Definition at line 1297 of file tmatrix.cc.

References CREATE_OUT0, and tmatrix_().

Referenced by TMatrixTest().

◆ tmatrix_tmd_test()

void tmatrix_tmd_test ( const Verbosity verbosity)

T-Matrix validation test.

Executes the standard test included with the double precision T-Matrix code for randomly oriented nonspherical particles. Should give the same as running the 3rdparty/tmatrix/tmatrix_tmd executable.

Author
Oliver Lemke

Definition at line 1361 of file tmatrix.cc.

References CREATE_OUT0, VectorView::get_c_array(), and tmd_().

Referenced by TMatrixTest().