ARTS  2.2.66
m_optproperties.cc
Go to the documentation of this file.
1 /* Copyright (C) 2002-2012
2  Sreerekha T.R. <rekha@uni-bremen.de>
3  Claudia Emde <claudia.emde@dlr.de>
4  Cory Davies <cory@met.ed.ac.uk>
5 
6  This program is free software; you can redistribute it and/or
7  modify it under the terms of the GNU General Public License as
8  published by the Free Software Foundation; either version 2 of the
9  License, or (at your option) any later version.
10 
11  This program is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with this program; if not, write to the Free Software
18  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
19  USA. */
20 
35 /*===========================================================================
36  === External declarations
37  ===========================================================================*/
38 
39 #include <cmath>
40 #include "arts.h"
41 #include "exceptions.h"
42 #include "array.h"
43 #include "matpackIII.h"
44 #include "matpackVII.h"
45 #include "logic.h"
46 #include "interpolation.h"
47 #include "messages.h"
48 #include "xml_io.h"
49 #include "optproperties.h"
50 #include "math_funcs.h"
51 #include "sorting.h"
52 #include "check_input.h"
53 #include "auto_md.h"
54 
55 extern const Numeric PI;
56 extern const Numeric DEG2RAD;
57 extern const Numeric RAD2DEG;
58 
59 #define PART_TYPE scat_data_array[i_pt].particle_type
60 #define F_DATAGRID scat_data_array[i_pt].f_grid
61 #define T_DATAGRID scat_data_array[i_pt].T_grid
62 #define ZA_DATAGRID scat_data_array[i_pt].za_grid
63 #define AA_DATAGRID scat_data_array[i_pt].aa_grid
64 #define PHA_MAT_DATA_RAW scat_data_array[i_pt].pha_mat_data //CPD: changed from pha_mat_data
65 #define EXT_MAT_DATA_RAW scat_data_array[i_pt].ext_mat_data //which wouldn't let me play with
66 #define ABS_VEC_DATA_RAW scat_data_array[i_pt].abs_vec_data //scat_data_array_mono.
67 #define PND_LIMIT 1e-12 // If particle number density is below this value,
68  // no transformations will be performed
69 
70 
71 /* Workspace method: Doxygen documentation will be auto-generated */
72 void pha_mat_sptFromData( // Output:
73  Tensor5& pha_mat_spt,
74  // Input:
75  const ArrayOfSingleScatteringData& scat_data_array,
76  const Vector& scat_za_grid,
77  const Vector& scat_aa_grid,
78  const Index& scat_za_index, // propagation directions
79  const Index& scat_aa_index,
80  const Index& f_index,
81  const Vector& f_grid,
82  const Numeric& rtp_temperature,
83  const Tensor4& pnd_field,
84  const Index& scat_p_index,
85  const Index& scat_lat_index,
86  const Index& scat_lon_index,
87  const Verbosity& verbosity
88  )
89 {
91 
92  out3 << "Calculate *pha_mat_spt* from database\n";
93 
94  const Index N_pt = scat_data_array.nelem();
95  const Index stokes_dim = pha_mat_spt.ncols();
96 
97  if (stokes_dim > 4 || stokes_dim < 1){
98  throw runtime_error("The dimension of the stokes vector \n"
99  "must be 1,2,3 or 4");
100  }
101 
102  assert( pha_mat_spt.nshelves() == N_pt );
103 
104  // Phase matrix in laboratory coordinate system. Dimensions:
105  // [frequency, za_inc, aa_inc, stokes_dim, stokes_dim]
106  Tensor5 pha_mat_data_int;
107 
108 
109  // Loop over the included particle_types
110  for (Index i_pt = 0; i_pt < N_pt; i_pt++)
111  {
112  // If the particle number density at a specific point in the atmosphere for
113  // the i_pt particle type is zero, we don't need to do the transfromation!
114  if (pnd_field(i_pt, scat_p_index, scat_lat_index, scat_lon_index) > PND_LIMIT)
115  {
116 
117  // First we have to transform the data from the coordinate system
118  // used in the database (depending on the kind of particle type
119  // specified by *particle_type*) to the laboratory coordinate sytem.
120 
121  // Frequency interpolation:
122 
123  // The data is interpolated on one frequency.
124  pha_mat_data_int.resize(PHA_MAT_DATA_RAW.nshelves(), PHA_MAT_DATA_RAW.nbooks(),
125  PHA_MAT_DATA_RAW.npages(), PHA_MAT_DATA_RAW.nrows(),
126  PHA_MAT_DATA_RAW.ncols());
127 
128 
129  // Gridpositions:
130  GridPos freq_gp;
131  gridpos(freq_gp, F_DATAGRID, f_grid[f_index]);
132 
133  GridPos t_gp;
134  gridpos(t_gp, T_DATAGRID, rtp_temperature);
135 
136  // Interpolationweights:
137  Vector itw(4);
138  interpweights(itw, freq_gp, t_gp);
139 
140  for (Index i_za_sca = 0; i_za_sca < PHA_MAT_DATA_RAW.nshelves() ; i_za_sca++)
141  {
142  for (Index i_aa_sca = 0; i_aa_sca < PHA_MAT_DATA_RAW.nbooks(); i_aa_sca++)
143  {
144  for (Index i_za_inc = 0; i_za_inc < PHA_MAT_DATA_RAW.npages();
145  i_za_inc++)
146  {
147  for (Index i_aa_inc = 0; i_aa_inc < PHA_MAT_DATA_RAW.nrows();
148  i_aa_inc++)
149  {
150  for (Index i = 0; i < PHA_MAT_DATA_RAW.ncols(); i++)
151  {
152  pha_mat_data_int(i_za_sca,
153  i_aa_sca, i_za_inc,
154  i_aa_inc, i) =
155  interp(itw,
156  PHA_MAT_DATA_RAW(joker, joker, i_za_sca,
157  i_aa_sca, i_za_inc,
158  i_aa_inc, i),
159  freq_gp, t_gp);
160  }
161  }
162  }
163  }
164  }
165 
166  // Do the transformation into the laboratory coordinate system.
167  for (Index za_inc_idx = 0; za_inc_idx < scat_za_grid.nelem();
168  za_inc_idx ++)
169  {
170  for (Index aa_inc_idx = 0; aa_inc_idx < scat_aa_grid.nelem();
171  aa_inc_idx ++)
172  {
173  pha_matTransform(pha_mat_spt
174  (i_pt, za_inc_idx, aa_inc_idx, joker, joker),
175  pha_mat_data_int,
177  PART_TYPE, scat_za_index, scat_aa_index,
178  za_inc_idx,
179  aa_inc_idx, scat_za_grid, scat_aa_grid,
180  verbosity);
181  }
182  }
183  }
184  }
185 }
186 
187 
188 /* Workspace method: Doxygen documentation will be auto-generated */
190  Tensor5& pha_mat_spt,
191  // Input:
192  const ArrayOfTensor7& pha_mat_sptDOITOpt,
193  const ArrayOfSingleScatteringData& scat_data_array_mono,
194  const Index& doit_za_grid_size,
195  const Vector& scat_aa_grid,
196  const Index& scat_za_index, // propagation directions
197  const Index& scat_aa_index,
198  const Numeric& rtp_temperature,
199  const Tensor4& pnd_field,
200  const Index& scat_p_index,
201  const Index& scat_lat_index,
202  const Index& scat_lon_index,
203  const Verbosity&)
204 {
205  // atmosphere_dim = 3
206  if (pnd_field.ncols() > 1)
207  {
208  assert(pha_mat_sptDOITOpt.nelem() == scat_data_array_mono.nelem());
209  // I assume that if the size is o.k. for one particle type is will
210  // also be o.k. for more particle types.
211  assert(pha_mat_sptDOITOpt[0].nlibraries() == scat_data_array_mono[0].T_grid.nelem());
212  assert(pha_mat_sptDOITOpt[0].nvitrines() == doit_za_grid_size);
213  assert(pha_mat_sptDOITOpt[0].nshelves() == scat_aa_grid.nelem() );
214  assert(pha_mat_sptDOITOpt[0].nbooks() == doit_za_grid_size);
215  assert(pha_mat_sptDOITOpt[0].npages() == scat_aa_grid.nelem());
216  }
217 
218  // atmosphere_dim = 1, only zenith angle required for scattered directions.
219  else if ( pnd_field.ncols() == 1 )
220  {
221  //assert(is_size(scat_theta, doit_za_grid_size, 1,
222  // doit_za_grid_size, scat_aa_grid.nelem()));
223 
224  assert(pha_mat_sptDOITOpt.nelem() == scat_data_array_mono.nelem());
225  // I assume that if the size is o.k. for one particle type is will
226  // also be o.k. for more particle types.
227  assert(pha_mat_sptDOITOpt[0].nlibraries() == scat_data_array_mono[0].T_grid.nelem());
228  assert(pha_mat_sptDOITOpt[0].nvitrines() == doit_za_grid_size);
229  assert(pha_mat_sptDOITOpt[0].nshelves() == 1);
230  assert(pha_mat_sptDOITOpt[0].nbooks() == doit_za_grid_size);
231  assert(pha_mat_sptDOITOpt[0].npages() == scat_aa_grid.nelem());
232  }
233 
234  assert(doit_za_grid_size > 0);
235 
236  // Create equidistant zenith angle grid
237  Vector za_grid;
238  nlinspace(za_grid, 0, 180, doit_za_grid_size);
239 
240  const Index N_pt = scat_data_array_mono.nelem();
241  const Index stokes_dim = pha_mat_spt.ncols();
242 
243  if (stokes_dim > 4 || stokes_dim < 1){
244  throw runtime_error("The dimension of the stokes vector \n"
245  "must be 1,2,3 or 4");
246  }
247 
248  assert( pha_mat_spt.nshelves() == N_pt );
249 
250  GridPos T_gp;
251  Vector itw(2);
252 
253  // Initialisation
254  pha_mat_spt = 0.;
255 
256  // Do the transformation into the laboratory coordinate system.
257  for (Index i_pt = 0; i_pt < N_pt; i_pt ++)
258  {
259  // If the particle number density at a specific point in the atmosphere
260  // for the i_pt particle type is zero, we don't need to do the
261  // transfromation!
262  if (pnd_field(i_pt, scat_p_index, scat_lat_index, scat_lon_index) > PND_LIMIT) //TRS
263  {
264  if( scat_data_array_mono[i_pt].T_grid.nelem() > 1)
265  {
266  ostringstream os;
267  os << "The temperature grid of the scattering data does not cover the \n"
268  "atmospheric temperature at cloud location. The data should \n"
269  "include the value T="<< rtp_temperature << " K. \n";
270  chk_interpolation_grids(os.str(), scat_data_array_mono[i_pt].T_grid, rtp_temperature);
271 
272  // Gridpositions:
273  gridpos(T_gp, scat_data_array_mono[i_pt].T_grid, rtp_temperature);
274  // Interpolationweights:
275  interpweights(itw, T_gp);
276  }
277 
278 
279 
280  for (Index za_inc_idx = 0; za_inc_idx < doit_za_grid_size;
281  za_inc_idx ++)
282  {
283  for (Index aa_inc_idx = 0; aa_inc_idx < scat_aa_grid.nelem();
284  aa_inc_idx ++)
285  {
286  if( scat_data_array_mono[i_pt].T_grid.nelem() == 1)
287  {
288  pha_mat_spt(i_pt, za_inc_idx, aa_inc_idx, joker, joker) =
289  pha_mat_sptDOITOpt[i_pt](0, scat_za_index,
290  scat_aa_index, za_inc_idx,
291  aa_inc_idx, joker, joker);
292  }
293 
294  // Temperature interpolation
295  else
296  {
297  for (Index i = 0; i< stokes_dim; i++)
298  {
299  for (Index j = 0; j< stokes_dim; j++)
300  {
301  pha_mat_spt(i_pt, za_inc_idx, aa_inc_idx, i, j)=
302  interp(itw,pha_mat_sptDOITOpt[i_pt]
303  (joker, scat_za_index,
304  scat_aa_index, za_inc_idx,
305  aa_inc_idx, i, j) , T_gp);
306  }
307  }
308  }
309  }
310  }
311  }// TRS
312  }
313 }
314 
315 
316 /* Workspace method: Doxygen documentation will be auto-generated */
317 void opt_prop_sptFromData(// Output and Input:
318  Tensor3& ext_mat_spt,
319  Matrix& abs_vec_spt,
320  // Input:
321  const ArrayOfSingleScatteringData& scat_data_array,
322  const Vector& scat_za_grid,
323  const Vector& scat_aa_grid,
324  const Index& scat_za_index, // propagation directions
325  const Index& scat_aa_index,
326  const Index& f_index,
327  const Vector& f_grid,
328  const Numeric& rtp_temperature,
329  const Tensor4& pnd_field,
330  const Index& scat_p_index,
331  const Index& scat_lat_index,
332  const Index& scat_lon_index,
333  const Verbosity& verbosity)
334 {
335 
336  const Index N_pt = scat_data_array.nelem();
337  const Index stokes_dim = ext_mat_spt.ncols();
338  const Numeric za_sca = scat_za_grid[scat_za_index];
339  const Numeric aa_sca = scat_aa_grid[scat_aa_index];
340 
341  if (stokes_dim > 4 || stokes_dim < 1){
342  throw runtime_error("The dimension of the stokes vector \n"
343  "must be 1,2,3 or 4");
344  }
345 
346  assert( ext_mat_spt.npages() == N_pt );
347  assert( abs_vec_spt.nrows() == N_pt );
348 
349  // Phase matrix in laboratory coordinate system. Dimensions:
350  // [frequency, za_inc, aa_inc, stokes_dim, stokes_dim]
351  Tensor3 ext_mat_data_int;
352  Tensor3 abs_vec_data_int;
353 
354  // Initialisation
355  ext_mat_spt = 0.;
356  abs_vec_spt = 0.;
357 
358 
359  // Loop over the included particle_types
360  for (Index i_pt = 0; i_pt < N_pt; i_pt++)
361  {
362  // If the particle number density at a specific point in the atmosphere for
363  // the i_pt particle type is zero, we don't need to do the transfromation
364 
365  if (pnd_field(i_pt, scat_p_index, scat_lat_index, scat_lon_index) > PND_LIMIT)
366  {
367  // First we have to transform the data from the coordinate system
368  // used in the database (depending on the kind of particle type
369  // specified by *particle_type*) to the laboratory coordinate sytem.
370 
371  // Frequency interpolation:
372 
373  // The data is interpolated on one frequency.
374  //
375  // Resize the variables for the interpolated data:
376  //
377  ext_mat_data_int.resize(EXT_MAT_DATA_RAW.npages(),
378  EXT_MAT_DATA_RAW.nrows(),
379  EXT_MAT_DATA_RAW.ncols());
380  //
381  abs_vec_data_int.resize(ABS_VEC_DATA_RAW.npages(),
382  ABS_VEC_DATA_RAW.nrows(),
383  ABS_VEC_DATA_RAW.ncols());
384 
385 
386  // Gridpositions:
387  GridPos freq_gp;
388  gridpos(freq_gp, F_DATAGRID, f_grid[f_index]);
389  GridPos t_gp;
390  Vector itw;
391 
392  if ( T_DATAGRID.nelem() > 1)
393  {
394  ostringstream os;
395  os << "The temperature grid of the scattering data does not cover the \n"
396  "atmospheric temperature at cloud location. The data should \n"
397  "include the value T="<< rtp_temperature << " K. \n";
398  chk_interpolation_grids(os.str(), T_DATAGRID, rtp_temperature);
399 
400  gridpos(t_gp, T_DATAGRID, rtp_temperature);
401 
402  // Interpolationweights:
403  itw.resize(4);
404  interpweights(itw, freq_gp, t_gp);
405 
406  for (Index i_za_sca = 0; i_za_sca < EXT_MAT_DATA_RAW.npages();
407  i_za_sca++)
408  {
409  for(Index i_aa_sca = 0; i_aa_sca < EXT_MAT_DATA_RAW.nrows();
410  i_aa_sca++)
411  {
412  //
413  // Interpolation of extinction matrix:
414  //
415  for (Index i = 0; i < EXT_MAT_DATA_RAW.ncols(); i++)
416  {
417  ext_mat_data_int(i_za_sca, i_aa_sca, i) =
419  i_za_sca, i_aa_sca, i),
420  freq_gp, t_gp);
421  }
422  }
423  }
424 
425  for (Index i_za_sca = 0; i_za_sca < ABS_VEC_DATA_RAW.npages();
426  i_za_sca++)
427  {
428  for(Index i_aa_sca = 0; i_aa_sca < ABS_VEC_DATA_RAW.nrows();
429  i_aa_sca++)
430  {
431  //
432  // Interpolation of absorption vector:
433  //
434  for (Index i = 0; i < ABS_VEC_DATA_RAW.ncols(); i++)
435  {
436  abs_vec_data_int(i_za_sca, i_aa_sca, i) =
437  interp(itw, ABS_VEC_DATA_RAW(joker, joker, i_za_sca,
438  i_aa_sca, i),
439  freq_gp, t_gp);
440  }
441  }
442  }
443  }
444  else
445  {
446  // Interpolationweights:
447  itw.resize(2);
448  interpweights(itw, freq_gp);
449 
450  for (Index i_za_sca = 0; i_za_sca < EXT_MAT_DATA_RAW.npages();
451  i_za_sca++)
452  {
453  for(Index i_aa_sca = 0; i_aa_sca < EXT_MAT_DATA_RAW.nrows();
454  i_aa_sca++)
455  {
456  //
457  // Interpolation of extinction matrix:
458  //
459  for (Index i = 0; i < EXT_MAT_DATA_RAW.ncols(); i++)
460  {
461  ext_mat_data_int(i_za_sca, i_aa_sca, i) =
462  interp(itw, EXT_MAT_DATA_RAW(joker, 0,
463  i_za_sca, i_aa_sca, i),
464  freq_gp);
465  }
466  }
467  }
468 
469  for (Index i_za_sca = 0; i_za_sca < ABS_VEC_DATA_RAW.npages();
470  i_za_sca++)
471  {
472  for(Index i_aa_sca = 0; i_aa_sca < ABS_VEC_DATA_RAW.nrows();
473  i_aa_sca++)
474  {
475  //
476  // Interpolation of absorption vector:
477  //
478  for (Index i = 0; i < ABS_VEC_DATA_RAW.ncols(); i++)
479  {
480  abs_vec_data_int(i_za_sca, i_aa_sca, i) =
481  interp(itw, ABS_VEC_DATA_RAW(joker, 0, i_za_sca,
482  i_aa_sca, i),
483  freq_gp);
484  }
485  }
486  }
487  }
488 
489 
490  //
491  // Do the transformation into the laboratory coordinate system.
492  //
493  // Extinction matrix:
494  //
495 
496 
497  ext_matTransform(ext_mat_spt(i_pt, joker, joker),
498  ext_mat_data_int,
500  za_sca, aa_sca,
501  verbosity);
502  //
503  // Absorption vector:
504  //
505  abs_vecTransform(abs_vec_spt(i_pt, joker),
506  abs_vec_data_int,
508  za_sca, aa_sca, verbosity);
509  }
510 
511  }
512 }
513 
514 
515 /* Workspace method: Doxygen documentation will be auto-generated */
516 void ext_matAddPart(Tensor3& ext_mat,
517  const Tensor3& ext_mat_spt,
518  const Tensor4& pnd_field,
519  const Index& atmosphere_dim,
520  const Index& scat_p_index,
521  const Index& scat_lat_index,
522  const Index& scat_lon_index,
523  const Verbosity&)
524 
525 {
526  Index N_pt = ext_mat_spt.npages();
527  Index stokes_dim = ext_mat_spt.nrows();
528 
529  Matrix ext_mat_part(stokes_dim, stokes_dim, 0.0);
530 
531 
532  if (stokes_dim > 4 || stokes_dim < 1){
533  throw runtime_error(
534  "The dimension of stokes vector can be "
535  "only 1,2,3, or 4");
536  }
537  if ( ext_mat_spt.ncols() != stokes_dim){
538 
539  throw runtime_error(" The columns of ext_mat_spt should "
540  "agree to stokes_dim");
541  }
542 
543  if (atmosphere_dim == 1)
544  {
545  // this is a loop over the different particle types
546  for (Index l = 0; l < N_pt; l++)
547  {
548 
549  // now the last two loops over the stokes dimension.
550  for (Index m = 0; m < stokes_dim; m++)
551  {
552  for (Index n = 0; n < stokes_dim; n++)
553  //summation of the product of pnd_field and
554  //ext_mat_spt.
555  ext_mat_part(m, n) +=
556  (ext_mat_spt(l, m, n) * pnd_field(l, scat_p_index, 0, 0));
557  }
558  }
559 
560  //Add particle extinction matrix to *ext_mat*.
561  ext_mat(0, Range(joker), Range(joker)) += ext_mat_part;
562  }
563 
564  if (atmosphere_dim == 3)
565  {
566 
567  // this is a loop over the different particle types
568  for (Index l = 0; l < N_pt; l++)
569  {
570 
571  // now the last two loops over the stokes dimension.
572  for (Index m = 0; m < stokes_dim; m++)
573  {
574  for (Index n = 0; n < stokes_dim; n++)
575  //summation of the product of pnd_field and
576  //ext_mat_spt.
577  ext_mat_part(m, n) += (ext_mat_spt(l, m, n) *
578  pnd_field(l, scat_p_index,
579  scat_lat_index,
580  scat_lon_index));
581 
582  }
583  }
584 
585  //Add particle extinction matrix to *ext_mat*.
586  ext_mat(0, Range(joker), Range(joker)) += ext_mat_part;
587 
588  }
589 
590 }
591 
592 
593 /* Workspace method: Doxygen documentation will be auto-generated */
594 void abs_vecAddPart(Matrix& abs_vec,
595  const Matrix& abs_vec_spt,
596  const Tensor4& pnd_field,
597  const Index& atmosphere_dim,
598  const Index& scat_p_index,
599  const Index& scat_lat_index,
600  const Index& scat_lon_index,
601  const Verbosity&)
602 
603 {
604  Index N_pt = abs_vec_spt.nrows();
605  Index stokes_dim = abs_vec_spt.ncols();
606 
607  Vector abs_vec_part(stokes_dim, 0.0);
608 
609  if ((stokes_dim > 4) || (stokes_dim <1)){
610  throw runtime_error("The dimension of stokes vector "
611  "can be only 1,2,3, or 4");
612  }
613 
614  if (atmosphere_dim == 1)
615  {
616  // this is a loop over the different particle types
617  for (Index l = 0; l < N_pt ; ++l)
618  {
619  // now the loop over the stokes dimension.
620  //(CE:) in the middle was l instead of m
621  for (Index m = 0; m < stokes_dim; ++m)
622  //summation of the product of pnd_field and
623  //abs_vec_spt.
624  abs_vec_part[m] +=
625  (abs_vec_spt(l, m) * pnd_field(l, scat_p_index, 0, 0));
626 
627  }
628  //Add the particle absorption
629  abs_vec(0, Range(joker)) += abs_vec_part;
630  }
631 
632  if (atmosphere_dim == 3)
633  {
634  // this is a loop over the different particle types
635  for (Index l = 0; l < N_pt ; ++l)
636  {
637 
638  // now the loop over the stokes dimension.
639  for (Index m = 0; m < stokes_dim; ++m)
640  //summation of the product of pnd_field and
641  //abs_vec_spt.
642  abs_vec_part[m] += (abs_vec_spt(l, m) *
643  pnd_field(l, scat_p_index,
644  scat_lat_index,
645  scat_lon_index));
646 
647  }
648  //Add the particle absorption
649  abs_vec(0,Range(joker)) += abs_vec_part;
650  }
651 }
652 
653 
654 /* Workspace method: Doxygen documentation will be auto-generated */
655 void ext_matInit(Tensor3& ext_mat,
656  const Vector& f_grid,
657  const Index& stokes_dim,
658  const Index& f_index,
659  const Verbosity& verbosity)
660 {
661  CREATE_OUT2;
662 
663  Index freq_dim;
664 
665  if( f_index < 0 )
666  freq_dim = f_grid.nelem();
667  else
668  freq_dim = 1;
669 
670  ext_mat.resize( freq_dim,
671  stokes_dim,
672  stokes_dim );
673  ext_mat = 0; // Initialize to zero!
674 
675  out2 << "Set dimensions of ext_mat as ["
676  << freq_dim << ","
677  << stokes_dim << ","
678  << stokes_dim << "] and initialized to 0.\n";
679 }
680 
681 
682 /* Workspace method: Doxygen documentation will be auto-generated */
683 void ext_matAddGas(Tensor3& ext_mat,
684  const Tensor4& propmat_clearsky,
685  const Verbosity&)
686 {
687  // Number of Stokes parameters:
688  const Index stokes_dim = ext_mat.ncols();
689 
690  // The second dimension of ext_mat must also match the number of
691  // Stokes parameters:
692  if ( stokes_dim != ext_mat.nrows() )
693  throw runtime_error("Row dimension of ext_mat inconsistent with "
694  "column dimension.");
695  if ( stokes_dim != propmat_clearsky.ncols() )
696  throw runtime_error("Col dimension of propmat_clearsky "
697  "inconsistent with col dimension in ext_mat.");
698 
699  // Number of frequencies:
700  const Index f_dim = ext_mat.npages();
701 
702  // This must be consistent with the second dimension of
703  // propmat_clearsky. Check this:
704  if ( f_dim != propmat_clearsky.npages() )
705  throw runtime_error("Frequency dimension of ext_mat and propmat_clearsky\n"
706  "are inconsistent in ext_matAddGas.");
707 
708  // Sum up absorption over all species.
709  // This gives us an absorption vector for all frequencies. Of course
710  // this includes the special case that there is only one frequency.
711  Tensor3 abs_total(f_dim,stokes_dim,stokes_dim);
712  abs_total = 0;
713 
714 // for ( Index i=0; i<f_dim; ++i )
715 // abs_total[i] = abs_scalar_gas(i,joker).sum();
716  for ( Index iv=0; iv<f_dim; ++iv )
717  for ( Index is1=0; is1<stokes_dim; ++is1 )
718  for ( Index is2=0; is2<stokes_dim; ++is2 )
719  abs_total(iv,is1,is2) += propmat_clearsky(joker,iv,is1,is2).sum();
720 
721  // Add the absorption value to all the elements:
722  ext_mat += abs_total;
723 
724 }
725 
726 
727 /* Workspace method: Doxygen documentation will be auto-generated */
728 void abs_vecInit(Matrix& abs_vec,
729  const Vector& f_grid,
730  const Index& stokes_dim,
731  const Index& f_index,
732  const Verbosity& verbosity)
733 {
734  CREATE_OUT2;
735 
736  Index freq_dim;
737 
738  if( f_index < 0 )
739  freq_dim = f_grid.nelem();
740  else
741  freq_dim = 1;
742 
743  abs_vec.resize( freq_dim,
744  stokes_dim );
745  abs_vec = 0; // Initialize to zero!
746 
747  out2 << "Set dimensions of abs_vec as ["
748  << freq_dim << ","
749  << stokes_dim << "] and initialized to 0.\n";
750 }
751 
752 
753 /* Workspace method: Doxygen documentation will be auto-generated */
754 void abs_vecAddGas(Matrix& abs_vec,
755  const Tensor4& propmat_clearsky,
756  const Verbosity&)
757 {
758  // Number of frequencies:
759  const Index f_dim = abs_vec.nrows();
760  const Index stokes_dim = abs_vec.ncols();
761 
762  // This must be consistent with the second dimension of
763  // propmat_clearsky. Check this:
764  if ( f_dim != propmat_clearsky.npages() )
765  throw runtime_error("Frequency dimension of abs_vec and propmat_clearsky\n"
766  "are inconsistent in abs_vecAddGas.");
767  if ( stokes_dim != propmat_clearsky.ncols() )
768  throw runtime_error("Stokes dimension of abs_vec and propmat_clearsky\n"
769  "are inconsistent in abs_vecAddGas.");
770 
771  // Loop all frequencies. Of course this includes the special case
772  // that there is only one frequency.
773  for ( Index i=0; i<f_dim; ++i )
774  {
775  // Sum up the columns of propmat_clearsky and add to the first
776  // element of abs_vec.
777  for(Index is = 0; is < stokes_dim;is++)
778  abs_vec(i,is) += propmat_clearsky(joker,i,is,0).sum();
779  }
780 
781  // We don't have to do anything about higher elements of abs_vec,
782  // since scalar gas absorption only influences the first element.
783 }
784 
785 
786 /* Workspace method: Doxygen documentation will be auto-generated */
787 /*
788 void ext_matAddGasZeeman( Tensor3& ext_mat,
789  const Tensor3& ext_mat_zee,
790  const Verbosity&)
791 {
792  // Number of Stokes parameters:
793  const Index stokes_dim = ext_mat.ncols();
794 
795  // The second dimension of ext_mat must also match the number of
796  // Stokes parameters:
797  if ( stokes_dim != ext_mat.nrows() )
798  throw runtime_error("Row dimension of ext_mat inconsistent with "
799  "column dimension.");
800 
801  for ( Index i=0; i<stokes_dim; ++i )
802  {
803  for ( Index j=0; j<stokes_dim; ++j )
804  {
805  // Add the zeeman extinction to extinction matrix.
806  ext_mat(joker,i,j) += ext_mat_zee(joker, i, j);
807  }
808 
809  }
810 }
811 */
812 
813 
814 /* Workspace method: Doxygen documentation will be auto-generated */
815 /*
816 void abs_vecAddGasZeeman( Matrix& abs_vec,
817  const Matrix& abs_vec_zee,
818  const Verbosity&)
819 {
820  // Number of Stokes parameters:
821  const Index stokes_dim = abs_vec_zee.ncols();
822  // that there is only one frequency.
823  for ( Index j=0; j<stokes_dim; ++j )
824  {
825  abs_vec(joker,j) += abs_vec_zee(joker,j);
826  }
827 }
828 */
829 
830 
831 /* Workspace method: Doxygen documentation will be auto-generated */
832 void pha_matCalc(Tensor4& pha_mat,
833  const Tensor5& pha_mat_spt,
834  const Tensor4& pnd_field,
835  const Index& atmosphere_dim,
836  const Index& scat_p_index,
837  const Index& scat_lat_index,
838  const Index& scat_lon_index,
839  const Verbosity&)
840 {
841 
842  Index N_pt = pha_mat_spt.nshelves();
843  Index Nza = pha_mat_spt.nbooks();
844  Index Naa = pha_mat_spt.npages();
845  Index stokes_dim = pha_mat_spt.nrows();
846 
847  pha_mat.resize(Nza, Naa, stokes_dim, stokes_dim);
848 
849  // Initialisation
850  pha_mat = 0.0;
851 
852  if (atmosphere_dim == 1)
853  {
854  // this is a loop over the different particle types
855  for (Index pt_index = 0; pt_index < N_pt; ++ pt_index)
856  {
857  // these are loops over zenith angle and azimuth angle
858  for (Index za_index = 0; za_index < Nza; ++ za_index)
859  {
860  for (Index aa_index = 0; aa_index < Naa; ++ aa_index)
861  {
862 
863  // now the last two loops over the stokes dimension.
864  for (Index stokes_index_1 = 0; stokes_index_1 < stokes_dim;
865  ++ stokes_index_1)
866  {
867  for (Index stokes_index_2 = 0; stokes_index_2 < stokes_dim;
868  ++ stokes_index_2)
869  //summation of the product of pnd_field and
870  //pha_mat_spt.
871  pha_mat(za_index, aa_index,
872  stokes_index_1, stokes_index_2) +=
873 
874  (pha_mat_spt(pt_index, za_index, aa_index,
875  stokes_index_1, stokes_index_2) *
876  pnd_field(pt_index,scat_p_index, 0, 0));
877  }
878  }
879  }
880  }
881  }
882 
883  if (atmosphere_dim == 3)
884  {
885  // this is a loop over the different particle types
886  for (Index pt_index = 0; pt_index < N_pt; ++ pt_index)
887  {
888 
889  // these are loops over zenith angle and azimuth angle
890  for (Index za_index = 0; za_index < Nza; ++ za_index)
891  {
892  for (Index aa_index = 0; aa_index < Naa; ++ aa_index)
893  {
894 
895  // now the last two loops over the stokes dimension.
896  for (Index i = 0; i < stokes_dim; ++ i)
897  {
898  for (Index j = 0; j < stokes_dim; ++ j)
899  {
900  //summation of the product of pnd_field and
901  //pha_mat_spt.
902  pha_mat(za_index, aa_index, i,j ) +=
903  (pha_mat_spt(pt_index, za_index, aa_index, i, j) *
904  pnd_field(pt_index, scat_p_index,
905  scat_lat_index, scat_lon_index));
906 
907 
908  }
909  }
910  }
911  }
912  }
913  }
914 }
915 
916 
917 
918 /* Workspace method: Doxygen documentation will be auto-generated */
919 void scat_data_arrayCheck(//Input:
920  const ArrayOfSingleScatteringData& scat_data_array,
921  const Numeric& threshold,
922  const Verbosity& verbosity)
923 {
924  CREATE_OUT2;
925 
926 /* JM121024: we do not really need to write the scatt_data to file again. we
927  usually have just read them in a couple of commands before!?
928  if wanted/needed for debug cases, just uncomment the 2 lines below */
929 // xml_write_to_file("SingleScatteringData", scat_data_array, FILE_TYPE_ASCII,
930 // verbosity);
931 
932  const Index N_pt = scat_data_array.nelem();
933 
934  // Loop over the included particle_types
935  for (Index i_pt = 0; i_pt < N_pt; i_pt++)
936  {
937  out2 << " particle " << i_pt << "\n";
938 
939  switch (PART_TYPE){
940 
942  {
943  for (Index f = 0; f < F_DATAGRID.nelem(); f++)
944  {
945  out2 << "frequency " << F_DATAGRID[f] << "Hz\n";
946  for (Index t = 0; t < T_DATAGRID.nelem(); t++)
947  {
948  out2 << "Temperature " << T_DATAGRID[t] << "K\n";
949 
951  (PHA_MAT_DATA_RAW(f, t, joker, 0, 0, 0, 0), ZA_DATAGRID);
952 
953  Numeric Cext_data = EXT_MAT_DATA_RAW(f,t,0,0,0);
954 
955  Numeric Cabs = Cext_data - Csca;
956 
957  Numeric Cabs_data = ABS_VEC_DATA_RAW(f,t,0,0,0);
958 
959  Numeric Csca_data = Cext_data - Cabs_data;
960 
961 
962  out2 << " Coefficients in database: "
963  << "Cext: " << Cext_data << " Cabs: " << Cabs_data
964  << " Csca: " << Csca_data << "\n"
965  << " Calculated coefficients: "
966  << "Cabs calc: " << Cabs
967  << " Csca calc: " << Csca << "\n"
968  << " Deviations "
969  << "Cabs: " << 1e2*Cabs/Cabs_data-1e2
970  << "% Csca: " << 1e2*Csca/Csca_data-1e2
971  << "% Alb: " << (Csca-Csca_data)/Cext_data << "\n";
972 
973 
974 // if (abs(Csca/Csca_data-1.)*Csca_data/Cext_data > threshold)
975 // equivalent to the above (it's actually the (absolute) albedo
976 // deviation!)
977  if (abs(Csca-Csca_data)/Cext_data > threshold)
978  {
979  ostringstream os;
980  os << " Deviations in scat_data_array too large:\n"
981  << " scat dev [%] " << 1e2*Csca/Csca_data-1e2
982  << " at albedo of " << Csca_data/Cext_data << "\n"
983  << " Check entry for particle " << i_pt << " at "
984  << f << ".frequency and " << t << ".temperature!\n";
985  throw runtime_error( os.str() );
986  }
987  }
988  }
989  break;
990  }
991 
992  default:
993  {
994  CREATE_OUT0;
995  out0 << " WARNING:\n"
996  << " scat_data_array consistency check not implemented (yet?!) for\n"
997  << " particle type " << PART_TYPE << "!\n";
998  }
999  }
1000  }
1001 }
1002 
1003 
1004 
1005 /* Workspace method: Doxygen documentation will be auto-generated */
1007  ArrayOfTensor7& pha_mat_sptDOITOpt,
1008  ArrayOfSingleScatteringData& scat_data_array_mono,
1009  //Input:
1010  const Index& doit_za_grid_size,
1011  const Vector& scat_aa_grid,
1013  scat_data_array,
1014  const Vector& f_grid,
1015  const Index& f_index,
1016  const Index& atmosphere_dim,
1017  const Index& stokes_dim,
1018  const Verbosity& verbosity)
1019 {
1020 
1021  // Interpolate all the data in frequency
1022  scat_data_array_monoCalc(scat_data_array_mono, scat_data_array, f_grid, f_index, verbosity);
1023 
1024  // For 1D calculation the scat_aa dimension is not required:
1025  Index N_aa_sca;
1026  if(atmosphere_dim == 1)
1027  N_aa_sca = 1;
1028  else
1029  N_aa_sca = scat_aa_grid.nelem();
1030 
1031  Vector za_grid;
1032  nlinspace(za_grid, 0, 180, doit_za_grid_size);
1033 
1034  assert( scat_data_array.nelem() == scat_data_array_mono.nelem() );
1035 
1036  Index N_pt = scat_data_array.nelem();
1037 
1038  pha_mat_sptDOITOpt.resize(N_pt);
1039 
1040  for (Index i_pt = 0; i_pt < N_pt; i_pt++)
1041  {
1042  Index N_T = scat_data_array_mono[i_pt].T_grid.nelem();
1043  pha_mat_sptDOITOpt[i_pt].resize(N_T, doit_za_grid_size, N_aa_sca,
1044  doit_za_grid_size, scat_aa_grid.nelem(),
1045  stokes_dim, stokes_dim);
1046 
1047  // Initialize:
1048  pha_mat_sptDOITOpt[i_pt]= 0.;
1049 
1050  // Calculate all scattering angles for all combinations of incoming
1051  // and scattered directions and interpolation.
1052  for (Index t_idx = 0; t_idx < N_T; t_idx ++)
1053  {
1054  // These are the scattered directions as called in *scat_field_calc*
1055  for (Index za_sca_idx = 0; za_sca_idx < doit_za_grid_size; za_sca_idx ++)
1056  {
1057  for (Index aa_sca_idx = 0; aa_sca_idx < N_aa_sca; aa_sca_idx ++)
1058  {
1059  // Integration is performed over all incoming directions
1060  for (Index za_inc_idx = 0; za_inc_idx < doit_za_grid_size;
1061  za_inc_idx ++)
1062  {
1063  for (Index aa_inc_idx = 0; aa_inc_idx <
1064  scat_aa_grid.nelem();
1065  aa_inc_idx ++)
1066  {
1067  pha_matTransform(pha_mat_sptDOITOpt[i_pt]
1068  (t_idx, za_sca_idx,
1069  aa_sca_idx, za_inc_idx, aa_inc_idx,
1070  joker, joker),
1071  scat_data_array_mono[i_pt].
1072  pha_mat_data
1073  (0,t_idx,joker,joker,joker,
1074  joker,joker),
1075  scat_data_array_mono[i_pt].za_grid,
1076  scat_data_array_mono[i_pt].aa_grid,
1077  scat_data_array_mono[i_pt].particle_type,
1078  za_sca_idx,
1079  aa_sca_idx,
1080  za_inc_idx,
1081  aa_inc_idx,
1082  za_grid,
1083  scat_aa_grid,
1084  verbosity);
1085  }
1086  }
1087  }
1088  }
1089  }
1090  }
1091  }
1092 
1093 
1094 /* Workspace method: Doxygen documentation will be auto-generated */
1096  const ArrayOfSingleScatteringData& scat_data_array,
1097  const Vector& f_grid,
1098  const Index& f_index,
1099  const Verbosity&)
1100 {
1101  //Extrapolation factor:
1102  //const Numeric extpolfac = 0.5;
1103 
1104  // Check, whether single scattering data contains the right frequencies:
1105  // The check was changed to allow extrapolation at the boundaries of the
1106  // frequency grid.
1107  for (Index i = 0; i<scat_data_array.nelem(); i++)
1108  {
1109  // check with extrapolation
1110  chk_interpolation_grids("scat_data_array.f_grid to f_grid",
1111  scat_data_array[i].f_grid,
1112  f_grid[f_index]);
1113 
1114  // old check without extrapolation
1115  /*if (scat_data_array[i].f_grid[0] > f_grid[f_index] ||
1116  scat_data_array[i].f_grid[scat_data_array[i].f_grid.nelem()-1] < f_grid[f_index])
1117  {
1118  ostringstream os;
1119  os << "Frequency of the scattering calculation " << f_grid[f_index]
1120  << " GHz is not contained \nin the frequency grid of the " << i+1
1121  << "the single scattering data file \n(*ParticleTypeAdd*). "
1122  << "Range:" << scat_data_array[i].f_grid[0]/1e9 <<" - "
1123  << scat_data_array[i].f_grid[scat_data_array[i].f_grid.nelem()-1]/1e9
1124  <<" GHz \n";
1125  throw runtime_error( os.str() );
1126  }*/
1127  }
1128 
1129  const Index N_pt = scat_data_array.nelem();
1130 
1131  //Initialise scat_data_array_mono
1132  scat_data_array_mono.resize(N_pt);
1133 
1134  // Loop over the included particle_types
1135  for (Index i_pt = 0; i_pt < N_pt; i_pt++)
1136  {
1137  // Gridpositions:
1138  GridPos freq_gp;
1139  gridpos(freq_gp, F_DATAGRID, f_grid[f_index]);
1140 
1141  // Interpolationweights:
1142  Vector itw(2);
1143  interpweights(itw, freq_gp);
1144 
1145  //Stuff that doesn't need interpolating
1146  scat_data_array_mono[i_pt].particle_type=PART_TYPE;
1147  scat_data_array_mono[i_pt].f_grid.resize(1);
1148  scat_data_array_mono[i_pt].f_grid=f_grid[f_index];
1149  scat_data_array_mono[i_pt].T_grid=scat_data_array[i_pt].T_grid;
1150  scat_data_array_mono[i_pt].za_grid=ZA_DATAGRID;
1151  scat_data_array_mono[i_pt].aa_grid=AA_DATAGRID;
1152 
1153  //Phase matrix data
1154  scat_data_array_mono[i_pt].pha_mat_data.resize(1,
1155  PHA_MAT_DATA_RAW.nvitrines(),
1156  PHA_MAT_DATA_RAW.nshelves(),
1157  PHA_MAT_DATA_RAW.nbooks(),
1158  PHA_MAT_DATA_RAW.npages(),
1159  PHA_MAT_DATA_RAW.nrows(),
1160  PHA_MAT_DATA_RAW.ncols());
1161 
1162  for (Index t_index = 0; t_index < PHA_MAT_DATA_RAW.nvitrines(); t_index ++)
1163  {
1164  for (Index i_za_sca = 0; i_za_sca < PHA_MAT_DATA_RAW.nshelves();
1165  i_za_sca++)
1166  {
1167  for (Index i_aa_sca = 0; i_aa_sca < PHA_MAT_DATA_RAW.nbooks();
1168  i_aa_sca++)
1169  {
1170  for (Index i_za_inc = 0; i_za_inc <
1171  PHA_MAT_DATA_RAW.npages();
1172  i_za_inc++)
1173  {
1174  for (Index i_aa_inc = 0;
1175  i_aa_inc < PHA_MAT_DATA_RAW.nrows();
1176  i_aa_inc++)
1177  {
1178  for (Index i = 0; i < PHA_MAT_DATA_RAW.ncols(); i++)
1179  {
1180  scat_data_array_mono[i_pt].pha_mat_data(0, t_index,
1181  i_za_sca,
1182  i_aa_sca,
1183  i_za_inc,
1184  i_aa_inc, i) =
1185  interp(itw,
1186  PHA_MAT_DATA_RAW(joker, t_index,
1187  i_za_sca,
1188  i_aa_sca, i_za_inc,
1189  i_aa_inc, i),
1190  freq_gp);
1191  }
1192  }
1193  }
1194  }
1195  }
1196  //Extinction matrix data
1197  scat_data_array_mono[i_pt].ext_mat_data.resize(1, T_DATAGRID.nelem(),
1198  EXT_MAT_DATA_RAW.npages(),
1199  EXT_MAT_DATA_RAW.nrows(),
1200  EXT_MAT_DATA_RAW.ncols());
1201  for (Index i_za_sca = 0; i_za_sca < EXT_MAT_DATA_RAW.npages();
1202  i_za_sca++)
1203  {
1204  for(Index i_aa_sca = 0; i_aa_sca < EXT_MAT_DATA_RAW.nrows();
1205  i_aa_sca++)
1206  {
1207  //
1208  // Interpolation of extinction matrix:
1209  //
1210  for (Index i = 0; i < EXT_MAT_DATA_RAW.ncols(); i++)
1211  {
1212  scat_data_array_mono[i_pt].ext_mat_data(0, t_index,
1213  i_za_sca, i_aa_sca, i)
1214  = interp(itw, EXT_MAT_DATA_RAW(joker, t_index, i_za_sca,
1215  i_aa_sca, i),
1216  freq_gp);
1217  }
1218  }
1219  }
1220  //Absorption vector data
1221  scat_data_array_mono[i_pt].abs_vec_data.resize(1, T_DATAGRID.nelem(),
1222  ABS_VEC_DATA_RAW.npages(),
1223  ABS_VEC_DATA_RAW.nrows(),
1224  ABS_VEC_DATA_RAW.ncols());
1225  for (Index i_za_sca = 0; i_za_sca < ABS_VEC_DATA_RAW.npages() ;
1226  i_za_sca++)
1227  {
1228  for(Index i_aa_sca = 0; i_aa_sca < ABS_VEC_DATA_RAW.nrows();
1229  i_aa_sca++)
1230  {
1231  //
1232  // Interpolation of absorption vector:
1233  //
1234  for (Index i = 0; i < ABS_VEC_DATA_RAW.ncols(); i++)
1235  {
1236  scat_data_array_mono[i_pt].abs_vec_data(0, t_index, i_za_sca,
1237  i_aa_sca, i) =
1238  interp(itw, ABS_VEC_DATA_RAW(joker, t_index, i_za_sca,
1239  i_aa_sca, i),
1240  freq_gp);
1241  }
1242  }
1243  }
1244  }
1245  }
1246 }
1247 
1248 
1249 /* Workspace method: Doxygen documentation will be auto-generated */
1250 void opt_prop_sptFromMonoData(// Output and Input:
1251  Tensor3& ext_mat_spt,
1252  Matrix& abs_vec_spt,
1253  // Input:
1254  const ArrayOfSingleScatteringData& scat_data_array_mono,
1255  const Vector& scat_za_grid,
1256  const Vector& scat_aa_grid,
1257  const Index& scat_za_index, // propagation directions
1258  const Index& scat_aa_index,
1259  const Numeric& rtp_temperature,
1260  const Tensor4& pnd_field,
1261  const Index& scat_p_index,
1262  const Index& scat_lat_index,
1263  const Index& scat_lon_index,
1264  const Verbosity& verbosity)
1265 {
1266  const Index N_pt = scat_data_array_mono.nelem();
1267  const Index stokes_dim = ext_mat_spt.ncols();
1268  const Numeric za_sca = scat_za_grid[scat_za_index];
1269  const Numeric aa_sca = scat_aa_grid[scat_aa_index];
1270 
1271  if (stokes_dim > 4 || stokes_dim < 1){
1272  throw runtime_error("The dimension of the stokes vector \n"
1273  "must be 1,2,3 or 4");
1274  }
1275 
1276  assert( ext_mat_spt.npages() == N_pt );
1277  assert( abs_vec_spt.nrows() == N_pt );
1278 
1279  // Initialisation
1280  ext_mat_spt = 0.;
1281  abs_vec_spt = 0.;
1282 
1283  GridPos t_gp;
1284 
1285  Vector itw(2);
1286 
1287  // Loop over the included particle_types
1288  for (Index i_pt = 0; i_pt < N_pt; i_pt++)
1289  {
1290  // If the particle number density at a specific point in the atmosphere for
1291  // the i_pt particle type is zero, we don't need to do the transfromation!
1292  if (pnd_field(i_pt, scat_p_index, scat_lat_index, scat_lon_index) > PND_LIMIT)
1293  {
1294 
1295  // First we have to transform the data from the coordinate system
1296  // used in the database (depending on the kind of particle type
1297  // specified by *particle_type*) to the laboratory coordinate sytem.
1298 
1299  //
1300  // Do the transformation into the laboratory coordinate system.
1301  //
1302  // Extinction matrix:
1303  //
1304  Index ext_npages = scat_data_array_mono[i_pt].ext_mat_data.npages();
1305  Index ext_nrows = scat_data_array_mono[i_pt].ext_mat_data.nrows();
1306  Index ext_ncols = scat_data_array_mono[i_pt].ext_mat_data.ncols();
1307  Index abs_npages = scat_data_array_mono[i_pt].abs_vec_data.npages();
1308  Index abs_nrows = scat_data_array_mono[i_pt].abs_vec_data.nrows();
1309  Index abs_ncols = scat_data_array_mono[i_pt].abs_vec_data.ncols();
1310  Tensor3 ext_mat_data1temp(ext_npages,ext_nrows,ext_ncols,0.0);
1311  Tensor3 abs_vec_data1temp(abs_npages,abs_nrows,abs_ncols,0.0);
1312 
1313  //Check that scattering data temperature range covers required temperature
1314  ConstVectorView t_grid = scat_data_array_mono[i_pt].T_grid;
1315 
1316  if (t_grid.nelem() > 1)
1317  {
1318  // if ((rtp_temperature<t_grid[0])||(rtp_temperature>t_grid[t_grid.nelem()-1]))
1319  // {
1320  // throw runtime_error("rtp_temperature outside scattering data temperature range");
1321  // }
1322 
1323  //interpolate over temperature
1324  gridpos(t_gp, scat_data_array_mono[i_pt].T_grid, rtp_temperature);
1325  interpweights(itw, t_gp);
1326  for (Index i_p = 0; i_p < ext_npages ; i_p++)
1327  {
1328  for (Index i_r = 0; i_r < ext_nrows ; i_r++)
1329  {
1330  for (Index i_c = 0; i_c < ext_ncols ; i_c++)
1331  {
1332  ext_mat_data1temp(i_p,i_r,i_c)=interp(itw,
1333  scat_data_array_mono[i_pt].ext_mat_data(0,joker,i_p,i_r,i_c),t_gp);
1334  }
1335  }
1336  }
1337  }
1338  else
1339  {
1340  ext_mat_data1temp =
1341  scat_data_array_mono[i_pt].ext_mat_data(0,0,joker,joker,joker);
1342  }
1343 
1344  ext_matTransform(ext_mat_spt(i_pt, joker, joker),
1345  ext_mat_data1temp,
1346  scat_data_array_mono[i_pt].za_grid,
1347  scat_data_array_mono[i_pt].aa_grid,
1348  scat_data_array_mono[i_pt].particle_type,
1349  za_sca, aa_sca,
1350  verbosity);
1351  //
1352  // Absorption vector:
1353  //
1354 
1355  if (t_grid.nelem() > 1)
1356  {
1357  //interpolate over temperature
1358  for (Index i_p = 0; i_p < abs_npages ; i_p++)
1359  {
1360  for (Index i_r = 0; i_r < abs_nrows ; i_r++)
1361  {
1362  for (Index i_c = 0; i_c < abs_ncols ; i_c++)
1363  {
1364  abs_vec_data1temp(i_p,i_r,i_c)=interp(itw,
1365  scat_data_array_mono[i_pt].abs_vec_data(0,joker,i_p,i_r,i_c),t_gp);
1366  }
1367  }
1368  }
1369  }
1370  else
1371  {
1372  abs_vec_data1temp =
1373  scat_data_array_mono[i_pt].abs_vec_data(0,0,joker,joker,joker);
1374  }
1375 
1376  abs_vecTransform(abs_vec_spt(i_pt, joker),
1377  abs_vec_data1temp,
1378  scat_data_array_mono[i_pt].za_grid,
1379  scat_data_array_mono[i_pt].aa_grid,
1380  scat_data_array_mono[i_pt].particle_type,
1381  za_sca, aa_sca,
1382  verbosity);
1383  }
1384 
1385  }
1386 }
1387 
1388 
1389 /* Workspace method: Doxygen documentation will be auto-generated */
1391  Tensor5& pha_mat_spt,
1392  // Input:
1393  const ArrayOfSingleScatteringData& scat_data_array_mono,
1394  const Index& doit_za_grid_size,
1395  const Vector& scat_aa_grid,
1396  const Index& scat_za_index, // propagation directions
1397  const Index& scat_aa_index,
1398  const Numeric& rtp_temperature,
1399  const Tensor4& pnd_field,
1400  const Index& scat_p_index,
1401  const Index& scat_lat_index,
1402  const Index& scat_lon_index,
1403  const Verbosity& verbosity)
1404 {
1405  CREATE_OUT3;
1406 
1407  out3 << "Calculate *pha_mat_spt* from scat_data_array_mono. \n";
1408 
1409  Vector za_grid;
1410  nlinspace(za_grid, 0, 180, doit_za_grid_size);
1411 
1412  const Index N_pt = scat_data_array_mono.nelem();
1413  const Index stokes_dim = pha_mat_spt.ncols();
1414 
1415 
1416 
1417  if (stokes_dim > 4 || stokes_dim < 1){
1418  throw runtime_error("The dimension of the stokes vector \n"
1419  "must be 1,2,3 or 4");
1420  }
1421 
1422  assert( pha_mat_spt.nshelves() == N_pt );
1423 
1424  GridPos T_gp;
1425  Vector itw(2);
1426 
1427  // Initialisation
1428  pha_mat_spt = 0.;
1429 
1430  for (Index i_pt = 0; i_pt < N_pt; i_pt ++)
1431  {
1432  // If the particle number density at a specific point in the atmosphere
1433  // for the i_pt particle type is zero, we don't need to do the
1434  // transfromation!
1435  if (pnd_field(i_pt, scat_p_index, scat_lat_index, scat_lon_index) >
1436  PND_LIMIT)
1437  {
1438  // Temporary phase matrix wich icludes the all temperatures.
1439  Tensor3 pha_mat_spt_tmp(scat_data_array_mono[i_pt].T_grid.nelem(),
1440  pha_mat_spt.nrows(), pha_mat_spt.ncols());
1441 
1442  pha_mat_spt_tmp = 0.;
1443 
1444  if( scat_data_array_mono[i_pt].T_grid.nelem() > 1)
1445  {
1446  ostringstream os;
1447  os << "The temperature grid of the scattering data does not cover the \n"
1448  "atmospheric temperature at cloud location. The data should \n"
1449  "include the value T="<< rtp_temperature << " K. \n";
1450  chk_interpolation_grids(os.str(), scat_data_array_mono[i_pt].T_grid, rtp_temperature);
1451 
1452  // Gridpositions:
1453  gridpos(T_gp, scat_data_array_mono[i_pt].T_grid, rtp_temperature);
1454  // Interpolationweights:
1455  interpweights(itw, T_gp);
1456  }
1457 
1458  // Do the transformation into the laboratory coordinate system.
1459  for (Index za_inc_idx = 0; za_inc_idx < doit_za_grid_size;
1460  za_inc_idx ++)
1461  {
1462  for (Index aa_inc_idx = 0; aa_inc_idx < scat_aa_grid.nelem();
1463  aa_inc_idx ++)
1464  {
1465  for (Index t_idx = 0; t_idx <
1466  scat_data_array_mono[i_pt].T_grid.nelem();
1467  t_idx ++)
1468  {
1469  pha_matTransform( pha_mat_spt_tmp(t_idx, joker, joker),
1470  scat_data_array_mono[i_pt].
1471  pha_mat_data
1472  (0,0,joker,joker,joker,
1473  joker,joker),
1474  scat_data_array_mono[i_pt].za_grid,
1475  scat_data_array_mono[i_pt].aa_grid,
1476  scat_data_array_mono[i_pt].particle_type,
1477  scat_za_index, scat_aa_index,
1478  za_inc_idx,
1479  aa_inc_idx, za_grid, scat_aa_grid,
1480  verbosity );
1481  }
1482  // Temperature interpolation
1483  if( scat_data_array_mono[i_pt].T_grid.nelem() > 1)
1484  {
1485  for (Index i = 0; i< stokes_dim; i++)
1486  {
1487  for (Index j = 0; j< stokes_dim; j++)
1488  {
1489  pha_mat_spt(i_pt, za_inc_idx, aa_inc_idx, i, j)=
1490  interp(itw, pha_mat_spt_tmp(joker, i, j), T_gp);
1491  }
1492  }
1493  }
1494  else // no temperatue interpolation required
1495  {
1496  pha_mat_spt(i_pt, za_inc_idx, aa_inc_idx, joker, joker) =
1497  pha_mat_spt_tmp(0, joker, joker);
1498  }
1499  }
1500  }
1501  }
1502  }
1503 }
1504 
1505 
1506 /* Workspace method: Doxygen documentation will be auto-generated */
1508  Tensor4& pnd_field,
1509  ArrayOfSingleScatteringData& scat_data_array,
1510  //WS Input:
1511  const Index& atmosphere_dim,
1512  const Index& cloudbox_on,
1513  const ArrayOfIndex& cloudbox_limits,
1514  const Tensor3& t_field,
1515  const Tensor3& z_field,
1516  const Matrix& z_surface,
1517  const Index& cloudbox_checked,
1518  const Verbosity& /*verbosity*/)
1519 {
1520  if (!cloudbox_checked)
1521  throw std::runtime_error("You must call *cloudbox_checkedCalc* before this method.");
1522 
1523  if (atmosphere_dim != 1)
1524  throw std::runtime_error("Merging particles only works with a 1D atmoshere");
1525 
1526  // Cloudbox on/off?
1527  if ( !cloudbox_on )
1528  {
1529  /* Must initialise pnd_field anyway; but empty */
1530  pnd_field.resize(0, 0, 0, 0);
1531  return;
1532  }
1533 
1534  // ------- setup new pnd_field and scat_data_array -------------------
1535  ArrayOfIndex limits(2);
1536  //pressure
1537  limits[0] = cloudbox_limits[0];
1538  limits[1] = cloudbox_limits[1] + 1;
1539 
1540  Tensor4 pnd_field_merged(limits[1] - limits[0],
1541  limits[1] - limits[0],
1542  1,
1543  1,
1544  0.);
1545 
1546  ArrayOfSingleScatteringData scat_data_array_merged;
1547 
1548  scat_data_array_merged.resize(pnd_field_merged.nbooks());
1549  for (Index sp = 0; sp < scat_data_array_merged.nelem(); sp++)
1550  {
1551  SingleScatteringData &this_part = scat_data_array_merged[sp];
1552  this_part.particle_type = scat_data_array[0].particle_type;
1553  this_part.description = "Merged particles";
1554  this_part.f_grid = scat_data_array[0].f_grid;
1555  this_part.za_grid = scat_data_array[0].za_grid;
1556  this_part.aa_grid = scat_data_array[0].aa_grid;
1557  this_part.description = "Merged particles";
1558  this_part.pha_mat_data.resize(scat_data_array[0].pha_mat_data.nlibraries(),
1559  1,
1560  scat_data_array[0].pha_mat_data.nshelves(),
1561  scat_data_array[0].pha_mat_data.nbooks(),
1562  scat_data_array[0].pha_mat_data.npages(),
1563  scat_data_array[0].pha_mat_data.nrows(),
1564  scat_data_array[0].pha_mat_data.ncols());
1565  this_part.ext_mat_data.resize(scat_data_array[0].ext_mat_data.nshelves(),
1566  1,
1567  scat_data_array[0].ext_mat_data.npages(),
1568  scat_data_array[0].ext_mat_data.nrows(),
1569  scat_data_array[0].ext_mat_data.ncols());
1570  this_part.abs_vec_data.resize(scat_data_array[0].abs_vec_data.nshelves(),
1571  1,
1572  scat_data_array[0].abs_vec_data.npages(),
1573  scat_data_array[0].abs_vec_data.nrows(),
1574  scat_data_array[0].abs_vec_data.ncols());
1575  this_part.pha_mat_data = 0.;
1576  this_part.ext_mat_data = 0.;
1577  this_part.abs_vec_data = 0.;
1578  this_part.T_grid.resize(1);
1579  this_part.T_grid[0] = t_field(sp, 0, 0);
1580  }
1581 
1582  // Check that all particles have same type and data dimensions
1583  SingleScatteringData &first_part = scat_data_array[0];
1584  for (Index i_pt = 1; i_pt < pnd_field.nbooks(); i_pt++)
1585  {
1586  SingleScatteringData &orig_part = scat_data_array[i_pt];
1587 
1588  if (orig_part.particle_type != first_part.particle_type)
1589  throw std::runtime_error("All particles must have the same type");
1590 
1591  if (orig_part.f_grid.nelem() != first_part.f_grid.nelem())
1592  throw std::runtime_error("All particles must have the same f_grid");
1593 
1594  if (!is_size(orig_part.pha_mat_data(joker, 0, joker, joker, joker, joker, joker),
1595  first_part.pha_mat_data.nlibraries(),
1596  first_part.pha_mat_data.nshelves(),
1597  first_part.pha_mat_data.nbooks(),
1598  first_part.pha_mat_data.npages(),
1599  first_part.pha_mat_data.nrows(),
1600  first_part.pha_mat_data.ncols()
1601  ))
1602  throw std::runtime_error("All particles must have the same pha_mat_data size"
1603  " (except for temperature).");
1604 
1605  if (!is_size(orig_part.ext_mat_data(joker, 0, joker, joker, joker),
1606  first_part.ext_mat_data.nshelves(),
1607  first_part.ext_mat_data.npages(),
1608  first_part.ext_mat_data.nrows(),
1609  first_part.ext_mat_data.ncols()
1610  ))
1611  throw std::runtime_error("All particles must have the same ext_mat_data size"
1612  " (except for temperature).");
1613 
1614  if (!is_size(orig_part.abs_vec_data(joker, 0, joker, joker, joker),
1615  first_part.abs_vec_data.nshelves(),
1616  first_part.abs_vec_data.npages(),
1617  first_part.abs_vec_data.nrows(),
1618  first_part.abs_vec_data.ncols()
1619  ))
1620  throw std::runtime_error("All particles must have the same abs_vec_data size"
1621  " (except for temperature).");
1622  }
1623 
1624  //-------- Start pnd_field_merged and scat_data_array_merged calculations--------------------
1625 
1626  GridPos T_gp;
1627  Vector itw(2);
1628 
1629  Index nlevels = pnd_field_merged.nbooks();
1630  // loop over pressure levels in cloudbox
1631  for (Index i_lv = 0; i_lv < nlevels-1; i_lv++)
1632  {
1633  pnd_field_merged(i_lv,i_lv,0,0) = 1.;
1634 
1635  SingleScatteringData &this_part = scat_data_array_merged[i_lv];
1636  for (Index i_pt = 0; i_pt < pnd_field.nbooks(); i_pt++)
1637  {
1638  SingleScatteringData &orig_part = scat_data_array[i_pt];
1639 
1640  // If the particle number density at a specific point in the atmosphere
1641  // for the i_pt particle type is zero, we don't need to do the
1642  // transformation!
1643  if (pnd_field(i_pt, i_lv, 0, 0) > PND_LIMIT) //TRS
1644  {
1645  Numeric temperature = this_part.T_grid[0];
1646  if( scat_data_array[i_pt].T_grid.nelem() > 1)
1647  {
1648  ostringstream os;
1649  os << "The temperature grid of the scattering data does not cover the \n"
1650  "atmospheric temperature at cloud location. The data should \n"
1651  "include the value T="<< this_part.T_grid[0] << " K. \n";
1652  chk_interpolation_grids(os.str(), orig_part.T_grid, temperature);
1653 
1654  // Gridpositions:
1655  gridpos(T_gp, orig_part.T_grid, temperature);
1656  // Interpolationweights:
1657  interpweights(itw, T_gp);
1658  }
1659 
1660  // Loop over frequencies
1661  for (Index i_f = 0; i_f < orig_part.pha_mat_data.nlibraries(); i_f++)
1662  {
1663  // Weighted sum of ext_mat_data and abs_vec_data
1664  if( scat_data_array[i_pt].T_grid.nelem() == 1)
1665  {
1666  this_part.ext_mat_data(i_f, 0, 0, 0, joker) +=
1667  pnd_field(i_pt, i_lv, 0, 0)
1668  * orig_part.ext_mat_data(i_f, 0, 0, 0, joker);
1669 
1670  this_part.abs_vec_data(i_f, 0, 0, 0, joker) +=
1671  pnd_field(i_pt, i_lv, 0, 0)
1672  * orig_part.abs_vec_data(i_f, 0, 0, 0, joker);
1673  }
1674  else
1675  {
1676  for (Index i = 0; i < orig_part.ext_mat_data.ncols(); i++)
1677  {
1678  // Temperature interpolation
1679  this_part.ext_mat_data(i_f, 0, 0, 0, i) +=
1680  pnd_field(i_pt, i_lv, 0, 0)
1681  * interp(itw,
1682  orig_part.ext_mat_data(i_f, joker, 0, 0, i),
1683  T_gp);
1684  }
1685  for (Index i = 0; i < orig_part.abs_vec_data.ncols(); i++)
1686  {
1687  // Temperature interpolation
1688  this_part.abs_vec_data(i_f, 0, 0, 0, i) +=
1689  pnd_field(i_pt, i_lv, 0, 0)
1690  * interp(itw,
1691  orig_part.abs_vec_data(i_f, joker, 0, 0, i),
1692  T_gp);
1693  }
1694  }
1695 
1696  // Loop over zenith angles
1697  for (Index i_za = 0; i_za < orig_part.pha_mat_data.nshelves(); i_za++)
1698  {
1699  // Weighted sum of pha_mat_data
1700  if( scat_data_array[i_pt].T_grid.nelem() == 1)
1701  {
1702  const Numeric pnd = pnd_field(i_pt, i_lv, 0, 0);
1703  for (Index i_s = 0; i_s < orig_part.pha_mat_data.ncols(); i_s++)
1704  {
1705  this_part.pha_mat_data(i_f, 0, i_za, 0, 0, 0, i_s) =
1706  pnd * orig_part.pha_mat_data(i_f, 0, i_za, 0, 0, 0, i_s);
1707  }
1708  }
1709  else
1710  {
1711  // Temperature interpolation
1712  for (Index i = 0; i < orig_part.pha_mat_data.ncols(); i++)
1713  {
1714  this_part.pha_mat_data(i_f, 0, i_za, 0, 0, 0, i) +=
1715  pnd_field(i_pt, i_lv, 0, 0)
1716  * interp(itw,
1717  orig_part.pha_mat_data(i_f, joker, i_za, 0, 0, 0, i),
1718  T_gp);
1719  }
1720  }
1721  }
1722  }
1723  }
1724  }
1725  }
1726 
1727  // Set new pnd_field at lowest altitude to 0 if the cloudbox doesn't touch the ground
1728  // The consistency for the original pnd_field has already been ensured by
1729  // cloudbox_checkedCalc
1730  if (z_field(cloudbox_limits[0], 0, 0) > z_surface(0, 0))
1731  pnd_field_merged(0, 0, 0, 0) = 0.;
1732 
1733  pnd_field = pnd_field_merged;
1734  scat_data_array = scat_data_array_merged;
1735 }
1736 
1737 
1738 /* Workspace method: Doxygen documentation will be auto-generated */
1740  //WS Output:
1741  Vector& meta_param,
1742  //WS Input:
1743  const ArrayOfScatteringMetaData& scat_meta_array,
1744  const ArrayOfIndex& scat_data_per_part_species,
1745  const String& meta_name,
1746  const Index& part_species_index,
1747  const Verbosity& ) //verbosity
1748 {
1749  if ( part_species_index<0 )
1750  {
1751  ostringstream os;
1752  os << "part_species_index can't be <0!";
1753  throw runtime_error( os.str() );
1754  }
1755 
1756  const Index nps = scat_data_per_part_species.nelem();
1757 
1758  // check that scat_data_per_part_species actually has at least
1759  // part_species_index elements
1760  if ( !(nps>part_species_index) )
1761  {
1762  ostringstream os;
1763  os << "Can not extract data for scattering species #"
1764  << part_species_index << "\n"
1765  << "because scat_data_per_part_species has only " << nps << " elements.";
1766  throw runtime_error( os.str() );
1767  }
1768 
1769  const Index npe = scat_data_per_part_species[part_species_index];
1770  meta_param.resize(npe);
1771 
1772  // calculate offset index scattering species part_species_index within
1773  // scat_meta_array
1774  Index indoff=0;
1775  for ( Index i=0; i<part_species_index; i++ )
1776  indoff+=scat_data_per_part_species[i];
1777 
1778  // eventually, copy the meta data into the output Vector
1779  for ( Index i=0; i<npe; i++ )
1780  {
1781  if ( meta_name=="volume" )
1782  meta_param[i] = scat_meta_array[indoff+i].volume;
1783  else if ( meta_name=="diameter_max" )
1784  meta_param[i] = scat_meta_array[indoff+i].diameter_max;
1785  else if ( meta_name=="density" )
1786  meta_param[i] = scat_meta_array[indoff+i].density;
1787  else if ( meta_name=="area_projected" )
1788  meta_param[i] = scat_meta_array[indoff+i].area_projected;
1789  else if ( meta_name=="aspect_ratio" )
1790  meta_param[i] = scat_meta_array[indoff+i].aspect_ratio;
1791  else
1792  {
1793  ostringstream os;
1794  os << "Meta parameter \"" << meta_name << "\"is unknown.";
1795  throw runtime_error( os.str() );
1796  }
1797  }
1798 }
Index npages() const
Returns the number of pages.
Definition: matpackV.cc:47
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:35
void opt_prop_sptFromMonoData(Tensor3 &ext_mat_spt, Matrix &abs_vec_spt, const ArrayOfSingleScatteringData &scat_data_array_mono, const Vector &scat_za_grid, const Vector &scat_aa_grid, const Index &scat_za_index, const Index &scat_aa_index, const Numeric &rtp_temperature, const Tensor4 &pnd_field, const Index &scat_p_index, const Index &scat_lat_index, const Index &scat_lon_index, const Verbosity &verbosity)
WORKSPACE METHOD: opt_prop_sptFromMonoData.
Index nrows() const
Returns the number of rows.
Definition: matpackV.cc:53
void ScatteringMergeParticles1D(Tensor4 &pnd_field, ArrayOfSingleScatteringData &scat_data_array, const Index &atmosphere_dim, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Tensor3 &t_field, const Tensor3 &z_field, const Matrix &z_surface, const Index &cloudbox_checked, const Verbosity &)
WORKSPACE METHOD: ScatteringMergeParticles1D.
Index nelem() const
Number of elements.
Definition: array.h:176
void resize(Index s, Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackV.cc:2430
Declarations having to do with the four output streams.
void ext_matAddPart(Tensor3 &ext_mat, const Tensor3 &ext_mat_spt, const Tensor4 &pnd_field, const Index &atmosphere_dim, const Index &scat_p_index, const Index &scat_lat_index, const Index &scat_lon_index, const Verbosity &)
WORKSPACE METHOD: ext_matAddPart.
The Vector class.
Definition: matpackI.h:556
void ext_matTransform(MatrixView ext_mat_lab, ConstTensor3View ext_mat_data, ConstVectorView za_datagrid, ConstVectorView aa_datagrid, const ParticleType &particle_type, const Numeric &za_sca, const Numeric &aa_sca, const Verbosity &verbosity)
Transformation of extinction matrix.
void abs_vecTransform(VectorView abs_vec_lab, ConstTensor3View abs_vec_data, ConstVectorView za_datagrid, ConstVectorView aa_datagrid, const ParticleType &particle_type, const Numeric &za_sca, const Numeric &aa_sca, const Verbosity &verbosity)
Transformation of absorption vector.
Numeric interp(ConstVectorView itw, ConstVectorView a, const GridPos &tc)
Red 1D Interpolate.
#define PND_LIMIT
The Tensor4 class.
Definition: matpackIV.h:383
The range class.
Definition: matpackI.h:148
Index nrows() const
Returns the number of rows.
Definition: matpackVII.cc:62
Index npages() const
Returns the number of pages.
Definition: matpackIV.cc:69
#define T_DATAGRID
#define ZA_DATAGRID
This file contains basic functions to handle XML data files.
Structure which describes the single scattering properties of a particle or a particle distribution...
Definition: optproperties.h:84
Header file for interpolation.cc.
Index nbooks() const
Returns the number of books.
Definition: matpackV.cc:41
const Numeric RAD2DEG
Index nrows() const
Returns the number of rows.
Definition: matpackIII.h:146
Index nelem() const
Returns the number of elements.
Definition: matpackI.cc:180
Contains sorting routines.
void chk_interpolation_grids(const String &which_interpolation, ConstVectorView old_grid, ConstVectorView new_grid, const Index order, const Numeric &extpolfac, const bool islog)
Check interpolation grids.
void ExtractFromMetaSinglePartSpecies(Vector &meta_param, const ArrayOfScatteringMetaData &scat_meta_array, const ArrayOfIndex &scat_data_per_part_species, const String &meta_name, const Index &part_species_index, const Verbosity &)
WORKSPACE METHOD: ExtractFromMetaSinglePartSpecies.
Structure to store a grid position.
Definition: interpolation.h:74
This file contains the definition of Array.
The implementation for String, the ARTS string class.
Definition: mystring.h:63
Index ncols() const
Returns the number of columns.
Definition: matpackI.cc:838
#define PHA_MAT_DATA_RAW
The Tensor3 class.
Definition: matpackIII.h:348
The global header file for ARTS.
Index nshelves() const
Returns the number of shelves.
Definition: matpackV.cc:35
void opt_prop_sptFromData(Tensor3 &ext_mat_spt, Matrix &abs_vec_spt, const ArrayOfSingleScatteringData &scat_data_array, const Vector &scat_za_grid, const Vector &scat_aa_grid, const Index &scat_za_index, const Index &scat_aa_index, const Index &f_index, const Vector &f_grid, const Numeric &rtp_temperature, const Tensor4 &pnd_field, const Index &scat_p_index, const Index &scat_lat_index, const Index &scat_lon_index, const Verbosity &verbosity)
WORKSPACE METHOD: opt_prop_sptFromData.
const Numeric DEG2RAD
void scat_data_array_monoCalc(ArrayOfSingleScatteringData &scat_data_array_mono, const ArrayOfSingleScatteringData &scat_data_array, const Vector &f_grid, const Index &f_index, const Verbosity &)
WORKSPACE METHOD: scat_data_array_monoCalc.
#define EXT_MAT_DATA_RAW
Index ncols() const
Returns the number of columns.
Definition: matpackIII.h:149
#define AA_DATAGRID
void pha_matTransform(MatrixView pha_mat_lab, ConstTensor5View pha_mat_data, ConstVectorView za_datagrid, ConstVectorView aa_datagrid, const ParticleType &particle_type, const Index &za_sca_idx, const Index &aa_sca_idx, const Index &za_inc_idx, const Index &aa_inc_idx, ConstVectorView scat_za_grid, ConstVectorView scat_aa_grid, const Verbosity &verbosity)
Transformation of phase matrix.
void pha_mat_sptFromMonoData(Tensor5 &pha_mat_spt, const ArrayOfSingleScatteringData &scat_data_array_mono, const Index &doit_za_grid_size, const Vector &scat_aa_grid, const Index &scat_za_index, const Index &scat_aa_index, const Numeric &rtp_temperature, const Tensor4 &pnd_field, const Index &scat_p_index, const Index &scat_lat_index, const Index &scat_lon_index, const Verbosity &verbosity)
WORKSPACE METHOD: pha_mat_sptFromMonoData.
The declarations of all the exception classes.
void nlinspace(Vector &x, const Numeric start, const Numeric stop, const Index n)
nlinspace
Definition: math_funcs.cc:261
#define ABS_VEC_DATA_RAW
void gridpos(ArrayOfGridPos &gp, ConstVectorView old_grid, ConstVectorView new_grid, const Numeric &extpolfac)
Set up a grid position Array.
#define PART_TYPE
#define abs(x)
Definition: continua.cc:20458
const Joker joker
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:29
The Matrix class.
Definition: matpackI.h:788
Index nlibraries() const
Returns the number of libraries.
Definition: matpackVII.cc:32
void abs_vecInit(Matrix &abs_vec, const Vector &f_grid, const Index &stokes_dim, const Index &f_index, const Verbosity &verbosity)
WORKSPACE METHOD: abs_vecInit.
void DoitScatteringDataPrepare(ArrayOfTensor7 &pha_mat_sptDOITOpt, ArrayOfSingleScatteringData &scat_data_array_mono, const Index &doit_za_grid_size, const Vector &scat_aa_grid, const ArrayOfSingleScatteringData &scat_data_array, const Vector &f_grid, const Index &f_index, const Index &atmosphere_dim, const Index &stokes_dim, const Verbosity &verbosity)
WORKSPACE METHOD: DoitScatteringDataPrepare.
void resize(Index p, Index r, Index c)
Resize function.
Definition: matpackIII.cc:862
Header file for logic.cc.
Index npages() const
Returns the number of pages.
Definition: matpackIII.h:143
This can be used to make arrays out of anything.
Definition: array.h:40
void resize(Index n)
Assignment operator from VectorView.
Definition: matpackI.cc:798
void pha_matCalc(Tensor4 &pha_mat, const Tensor5 &pha_mat_spt, const Tensor4 &pnd_field, const Index &atmosphere_dim, const Index &scat_p_index, const Index &scat_lat_index, const Index &scat_lon_index, const Verbosity &)
WORKSPACE METHOD: pha_matCalc.
A constant view of a Vector.
Definition: matpackI.h:292
Index nshelves() const
Returns the number of shelves.
Definition: matpackVII.cc:44
Index npages() const
Returns the number of pages.
Definition: matpackVII.cc:56
Index nbooks() const
Returns the number of books.
Definition: matpackIV.cc:63
void ext_matAddGas(Tensor3 &ext_mat, const Tensor4 &propmat_clearsky, const Verbosity &)
WORKSPACE METHOD: ext_matAddGas.
void abs_vecAddPart(Matrix &abs_vec, const Matrix &abs_vec_spt, const Tensor4 &pnd_field, const Index &atmosphere_dim, const Index &scat_p_index, const Index &scat_lat_index, const Index &scat_lon_index, const Verbosity &)
WORKSPACE METHOD: abs_vecAddPart.
Index ncols() const
Returns the number of columns.
Definition: matpackV.cc:59
#define CREATE_OUT0
Definition: messages.h:211
#define CREATE_OUT3
Definition: messages.h:214
bool is_size(ConstVectorView x, const Index &n)
Verifies that the size of x is l.
Definition: logic.cc:91
Index ncols() const
Returns the number of columns.
Definition: matpackVII.cc:68
ParticleType particle_type
Definition: optproperties.h:85
const Numeric PI
Index ncols() const
Returns the number of columns.
Definition: matpackIV.cc:81
void pha_mat_sptFromData(Tensor5 &pha_mat_spt, const ArrayOfSingleScatteringData &scat_data_array, const Vector &scat_za_grid, const Vector &scat_aa_grid, const Index &scat_za_index, const Index &scat_aa_index, const Index &f_index, const Vector &f_grid, const Numeric &rtp_temperature, const Tensor4 &pnd_field, const Index &scat_p_index, const Index &scat_lat_index, const Index &scat_lon_index, const Verbosity &verbosity)
WORKSPACE METHOD: pha_mat_sptFromData.
void scat_data_arrayCheck(const ArrayOfSingleScatteringData &scat_data_array, const Numeric &threshold, const Verbosity &verbosity)
WORKSPACE METHOD: scat_data_arrayCheck.
void interpweights(VectorView itw, const GridPos &tc)
Red 1D interpolation weights.
void resize(Index l, Index v, Index s, Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackVII.cc:5379
#define CREATE_OUT2
Definition: messages.h:213
void abs_vecAddGas(Matrix &abs_vec, const Tensor4 &propmat_clearsky, const Verbosity &)
WORKSPACE METHOD: abs_vecAddGas.
Scattering database structure and functions.
void pha_mat_sptFromDataDOITOpt(Tensor5 &pha_mat_spt, const ArrayOfTensor7 &pha_mat_sptDOITOpt, const ArrayOfSingleScatteringData &scat_data_array_mono, const Index &doit_za_grid_size, const Vector &scat_aa_grid, const Index &scat_za_index, const Index &scat_aa_index, const Numeric &rtp_temperature, const Tensor4 &pnd_field, const Index &scat_p_index, const Index &scat_lat_index, const Index &scat_lon_index, const Verbosity &)
WORKSPACE METHOD: pha_mat_sptFromDataDOITOpt.
Index nrows() const
Returns the number of rows.
Definition: matpackI.cc:832
Index nbooks() const
Returns the number of books.
Definition: matpackVII.cc:50
#define F_DATAGRID
void resize(Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackIV.cc:1403
Numeric AngIntegrate_trapezoid(ConstMatrixView Integrand, ConstVectorView za_grid, ConstVectorView aa_grid)
AngIntegrate_trapezoid.
Definition: math_funcs.cc:327
The Tensor5 class.
Definition: matpackV.h:451
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1580
void ext_matInit(Tensor3 &ext_mat, const Vector &f_grid, const Index &stokes_dim, const Index &f_index, const Verbosity &verbosity)
WORKSPACE METHOD: ext_matInit.