Home > atmlab > scattering > tmatrix.m

tmatrix

PURPOSE ^

TMATRIX T-matrix scattering data, following ARTS

SYNOPSIS ^

function [D,M] = tmatrix( workfolder,f_grid,t_grid,za_grid,aa_grid,N,pshape,ptype,pdiameter,aspect_ratio)

DESCRIPTION ^

 TMATRIX   T-matrix scattering data, following ARTS

    Basically an interface to ARTS for obtaining T-matrix data. 
    For detailed information, such as options for *pshape* see: 
       arts -d scat_meta_arrayAddTmatrix
    
    This function can handle combinations of shapes, types and sizes
    (and then differs a bit from scat_meta_arrayAddTmatrix).
    The variables that can be varied are N, pshape, ptype, pdiameter
    and aspect_ratio. These variables must have length 1 or *np*, 
    where *np* is the number of particles considered. If one or 
    several of these variables have length 1, these data are assumed 
    to be valid for all particles. For example, to compare scattering  
    of ice and water for identical particles:

      % set shape, type, diameter and aratio to scalar values
      N(1) = complex_refr_indexFromFunc(fg,tg,@eps_water_liebe93,@sqrt);
      N(2) = complex_refr_indexFromFunc(fg,tg,@eps_ice_liebe93,@sqrt);
      D    = tmatrix( ...
    
    The grids are common for all particles.

 FORMAT   [D,M] = tmatrix( workfolder,f_grid,t_grid,za_grid,aa_grid,N, ...
                       pshape,ptype,pdiameter,aspect_ratio)
        
 OUT   D            Array of ARTS single scattering data
       M            The "meta data" corresponding to D.
 IN    workfolder   A folder where temporary files can be stored.
                    If set to [], a temporary folder is created.
       f_grid       Frequency grid.
       t_grid       Temperature grid.
       za_grid      Zenith angle grid.
       aa_grid      Azimuth angle grid.
       N            Data of complex refractive index.
       pshape       Particle shape. A string or a cell string array.
       ptype        Particle type. A string or a cell string array.
       pdiameter    Particle diameter.
       aspect_ratio Particle aspect ratio.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

DOWNLOAD ^

tmatrix.m

SOURCE CODE ^

0001 % TMATRIX   T-matrix scattering data, following ARTS
0002 %
0003 %    Basically an interface to ARTS for obtaining T-matrix data.
0004 %    For detailed information, such as options for *pshape* see:
0005 %       arts -d scat_meta_arrayAddTmatrix
0006 %
0007 %    This function can handle combinations of shapes, types and sizes
0008 %    (and then differs a bit from scat_meta_arrayAddTmatrix).
0009 %    The variables that can be varied are N, pshape, ptype, pdiameter
0010 %    and aspect_ratio. These variables must have length 1 or *np*,
0011 %    where *np* is the number of particles considered. If one or
0012 %    several of these variables have length 1, these data are assumed
0013 %    to be valid for all particles. For example, to compare scattering
0014 %    of ice and water for identical particles:
0015 %
0016 %      % set shape, type, diameter and aratio to scalar values
0017 %      N(1) = complex_refr_indexFromFunc(fg,tg,@eps_water_liebe93,@sqrt);
0018 %      N(2) = complex_refr_indexFromFunc(fg,tg,@eps_ice_liebe93,@sqrt);
0019 %      D    = tmatrix( ...
0020 %
0021 %    The grids are common for all particles.
0022 %
0023 % FORMAT   [D,M] = tmatrix( workfolder,f_grid,t_grid,za_grid,aa_grid,N, ...
0024 %                       pshape,ptype,pdiameter,aspect_ratio)
0025 %
0026 % OUT   D            Array of ARTS single scattering data
0027 %       M            The "meta data" corresponding to D.
0028 % IN    workfolder   A folder where temporary files can be stored.
0029 %                    If set to [], a temporary folder is created.
0030 %       f_grid       Frequency grid.
0031 %       t_grid       Temperature grid.
0032 %       za_grid      Zenith angle grid.
0033 %       aa_grid      Azimuth angle grid.
0034 %       N            Data of complex refractive index.
0035 %       pshape       Particle shape. A string or a cell string array.
0036 %       ptype        Particle type. A string or a cell string array.
0037 %       pdiameter    Particle diameter.
0038 %       aspect_ratio Particle aspect ratio.
0039 
0040 % 2013-08-19   Created by Patrick Eriksson.
0041 
0042 function [D,M] = tmatrix( workfolder,f_grid,t_grid,za_grid,aa_grid,N, ...
0043                       pshape,ptype,pdiameter,aspect_ratio)
0044 
0045 
0046 %- Basic checks of input
0047 %
0048 rqre_datatype( f_grid, @istensor1 );
0049 rqre_datatype( t_grid, @istensor1 );
0050 rqre_datatype( za_grid, @istensor1 );
0051 rqre_datatype( aa_grid, @istensor1 );
0052 rqre_datatype( N, @isstruct );
0053 rqre_datatype( pshape, { @ischar, @iscellstr } );
0054 rqre_datatype( ptype, { @ischar, @iscellstr } );
0055 rqre_datatype( pdiameter, @istensor1 );
0056 rqre_datatype( aspect_ratio, @istensor1 );
0057 
0058 %- Convert char input to cellstrs
0059 %
0060 if ischar( pshape ) 
0061   dummy = pshape;
0062   clear pshape;
0063   pshape{1} = dummy; 
0064   clear dummy;
0065 end
0066 if ischar( ptype )
0067   dummy = ptype;
0068   clear ptype;
0069   ptype{1} = dummy;  
0070   clear dummy;
0071 end
0072 
0073 %- How many particle types?
0074 %
0075 np = max( [ length(N), length(pshape), length(ptype), ...
0076             length(pdiameter), length(aspect_ratio) ] );
0077 
0078 %- Check if length 1 or np
0079 %
0080 if ~( length(N)==1 | length(N)==np )
0081   error( 'In this case, length of *N* must be 1 or %d', np );
0082 end
0083 if ~( length(pshape)==1 | length(pshape)==np )
0084   error( 'In this case, length of *pshape* must be 1 or %d', np );
0085 end
0086 if ~( length(ptype)==1 | length(ptype)==np )
0087   error( 'In this case, length of *ptype* must be 1 or %d', np );
0088 end
0089 if ~( length(pdiameter)==1 | length(pdiameter)==np )
0090   error( 'In this case, length of *pdiameter* must be 1 or %d', np );
0091 end
0092 if ~( length(aspect_ratio)==1 | length(aspect_ratio)==np )
0093   error( 'In this case, length of *aspect_ratio* must be 1 or %d', np );
0094 end
0095 
0096   
0097 %- Create workfolder?
0098 %
0099 if isempty( workfolder )
0100   workfolder = create_tmpfolder;
0101   cu = onCleanup( @()delete_tmpfolder( workfolder ) );
0102 end
0103 
0104 
0105 %- Start defining cfile
0106 %
0107 S{1} = 'Arts2{';
0108 
0109 
0110 %- Grids
0111 %
0112 filename = fullfile( workfolder, 'scat_f_grid.xml' ); 
0113 xmlStore( filename, f_grid, 'Vector' );
0114 S{end+1} = 'VectorCreate(scat_f_grid)';
0115 S{end+1} = sprintf( 'ReadXML(scat_f_grid,"%s")', filename );
0116 %
0117 filename = fullfile( workfolder, 'scat_t_grid.xml' ); 
0118 xmlStore( filename, t_grid, 'Vector' );
0119 S{end+1} = 'VectorCreate(scat_t_grid)';
0120 S{end+1} = sprintf( 'ReadXML(scat_t_grid,"%s")', filename );
0121 %
0122 filename = fullfile( workfolder, 'scat_za_grid.xml' ); 
0123 xmlStore( filename, za_grid, 'Vector' );
0124 S{end+1} = sprintf( 'ReadXML(scat_za_grid,"%s")', filename );
0125 %
0126 filename = fullfile( workfolder, 'scat_aa_grid.xml' ); 
0127 xmlStore( filename, aa_grid, 'Vector' );
0128 S{end+1} = sprintf( 'ReadXML(scat_aa_grid,"%s")', filename );
0129   
0130 
0131 %- Create meta data
0132 %
0133 S{end+1} = 'scat_meta_arrayInit';
0134 %
0135 for i = 1 : np
0136   
0137   % Store complex_refr_index
0138   if i==1 | length(N) > 1
0139     filename = fullfile( workfolder, ...
0140                               sprintf('complex_refr_index_%d.xml',i) ); 
0141     xmlStore( filename, N(i), 'GriddedField3' );
0142     S{end+1} = sprintf( 'ReadXML(complex_refr_index,"%s")', filename );
0143   end
0144   
0145   % Expand meta data
0146   S{end+1} = 'scat_meta_arrayAddTmatrixOldVersion(';
0147   S{end+1} = '  complex_refr_index = complex_refr_index,';  
0148   if length(pshape) > 1
0149     S{end+1} = sprintf( '  shape              = "%s",', pshape{i} );
0150     desc{i}  = sprintf( '%s', pshape{i} );
0151   else
0152     S{end+1} = sprintf( '  shape              = "%s",', pshape{1} );
0153     desc{i}  = sprintf( '%s', pshape{1} );
0154   end
0155   if length(ptype) > 1
0156     S{end+1} = sprintf( '  particle_type      = "%s",', ptype{i} );
0157   else
0158     S{end+1} = sprintf( '  particle_type      = "%s",', ptype{1} );
0159   end
0160   if length(aspect_ratio) > 1
0161     S{end+1} = sprintf( '  aspect_ratio       = %.6d,', aspect_ratio(i) );
0162     desc{i}  = sprintf( '%s, aspect ratio %.3f', desc{i}, aspect_ratio(i) );
0163   else
0164     S{end+1} = sprintf( '  aspect_ratio       = %.6d,', aspect_ratio(1) );
0165     desc{i}  = sprintf( '%s, aspect ratio %.3f', desc{i}, aspect_ratio(1) );
0166   end  
0167   if length(pdiameter) > 1
0168     S{end+1} = sprintf( '  diameter_grid      = [%.6e],', pdiameter(i) );
0169     desc{i}  = sprintf( '%s, %.3f um', desc{i}, 1e6*pdiameter(i) );
0170   else
0171     S{end+1} = sprintf( '  diameter_grid      = [%.6e],', pdiameter(1) );
0172     desc{i}  = sprintf( '%s, %.3f um', desc{i}, 1e6*pdiameter(1) );
0173   end
0174   S{end+1} = '  scat_f_grid        = scat_f_grid,';
0175   S{end+1} = '  scat_T_grid        = scat_t_grid )';
0176 end
0177 
0178 
0179 %- T-matrix part and save
0180 %
0181 S{end+1} = 'scat_data_arrayFromMeta(';
0182 S{end+1} = '  za_grid   = scat_za_grid,';
0183 S{end+1} = '  aa_grid   = scat_aa_grid,';
0184 S{end+1} = '  precision = 0.001 )';
0185 S{end+1} = 'WriteXML( "ascii", scat_meta_array )';
0186 S{end+1} = 'WriteXML( "ascii", scat_data_array )';
0187 
0188 
0189 %- Create cfile
0190 %
0191 S{end+1} = '}';
0192 cfile = fullfile( workfolder, 'cfile.arts' );
0193 strs2file( cfile, S );
0194 
0195 
0196 %- Run arts and load result
0197 %
0198 arts( cfile );
0199 %
0200 D = xmlLoad( fullfile( workfolder, 'cfile.scat_data_array.xml' ) );
0201 %
0202 for i = 1 : np
0203   D{i}.description = desc{i};
0204 end  
0205   
0206 
0207 if nargout > 1
0208   M = xmlLoad( fullfile( workfolder, ...
0209                                   'cfile.scat_meta_array.xml' ) );
0210 end

Generated on Mon 15-Sep-2014 13:31:28 by m2html © 2005