ARTS  2.3.1285(git:92a29ea9-dirty)
m_fos.cc
Go to the documentation of this file.
1 /* Copyright (C) 2013
2  Patrick Eriksson <patrick.eriksson@chalmers.se>
3 
4  This program is free software; you can redistribute it and/or modify it
5  under the terms of the GNU General Public License as published by the
6  Free Software Foundation; either version 2, or (at your option) any
7  later version.
8 
9  This program is distributed in the hope that it will be useful,
10  but WITHOUT ANY WARRANTY; without even the implied warranty of
11  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12  GNU General Public License for more details.
13 
14  You should have received a copy of the GNU General Public License
15  along with this program; if not, write to the Free Software
16  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
17  USA. */
18 
19 /*===========================================================================
20  === File description
21  ===========================================================================*/
22 
34 /*===========================================================================
35  === External declarations
36  ===========================================================================*/
37 
38 #include <cmath>
39 #include <stdexcept>
40 #include "arts.h"
41 #include "arts_omp.h"
42 #include "auto_md.h"
43 #include "doit.h"
44 #include "m_xml.h"
45 #include "math_funcs.h"
46 #include "montecarlo.h"
47 #include "rte.h"
48 
49 extern const Numeric DEG2RAD;
50 extern const Numeric RAD2DEG;
51 extern const Numeric PI;
52 
53 // FOS implemented as an internal function, to allow an recursive algorithm
54 /*
55 void fos(
56  Workspace& ws,
57  Matrix& iy,
58  ArrayOfTensor4& iy_aux,
59  Ppath& ppath,
60  ArrayOfTensor3& diy_dx,
61  const Index& stokes_dim,
62  const Vector& f_grid,
63  const Index& atmosphere_dim,
64  const Vector& p_grid,
65  const Tensor3& z_field,
66  const Tensor3& t_field,
67  const Tensor4& vmr_field,
68  const ArrayOfArrayOfSpeciesTag& abs_species,
69  const Tensor3& wind_u_field,
70  const Tensor3& wind_v_field,
71  const Tensor3& wind_w_field,
72  const Tensor3& mag_u_field,
73  const Tensor3& mag_v_field,
74  const Tensor3& mag_w_field,
75  const Index& cloudbox_on,
76  const ArrayOfIndex& cloudbox_limits,
77  const Tensor4& pnd_field,
78  const ArrayOfArrayOfSingleScatteringData& scat_data,
79  const Matrix& particle_masses,
80  const String& iy_unit,
81  const ArrayOfString& iy_aux_vars,
82  const Index& jacobian_do,
83  const Agenda& ppath_agenda,
84  const Agenda& propmat_clearsky_agenda,
85  const Agenda& iy_main_agenda,
86  const Agenda& iy_space_agenda,
87  const Agenda& iy_surface_agenda,
88  const Index& iy_agenda_call1,
89  const Tensor3& iy_transmission,
90  const Vector& rte_pos,
91  const Vector& rte_los,
92  const Vector& rte_pos2,
93  const Numeric& rte_alonglos_v,
94  const Numeric& ppath_lmax,
95  const Numeric& ppath_lraytrace,
96  const Matrix& fos_scatint_angles,
97  const Vector& fos_iyin_za_angles,
98  const Index& fos_za_interporder,
99  const Index& fos_n,
100  const Index& fos_i,
101  const Verbosity& verbosity )
102 {
103  // A temporary restriction
104  if( atmosphere_dim > 1 )
105  throw runtime_error( "FOS is so far only handling 1D atmospheres." );
106 
107  assert( fos_i >= 0 && fos_i <= fos_n );
108 
109 
110  // Determine propagation path
111  //
112  ppath_agendaExecute( ws, ppath, ppath_lmax, ppath_lraytrace,
113  rte_pos, rte_los, rte_pos2,
114  0, 0, t_field, z_field, vmr_field, f_grid,
115  ppath_agenda );
116 
117  // Some basic sizes
118  //
119  const Index nf = f_grid.nelem();
120  const Index ns = stokes_dim;
121  const Index np = ppath.np;
122 
123  // The below copied from iyEmission. Activate for-loop when jacobians
124  // introduced.
125 
126  // Set up variable with index of species where we want abs_per_species.
127  // This variable can below be extended in iy_aux part.
128  //
129  ArrayOfIndex iaps(0);
130  //
131  //for( Index i=0; i<jac_species_i.nelem(); i++ )
132  // {
133  // if( jac_species_i[i] >= 0 )
134  // { iaps.push_back( jac_species_i[i] ); }
135  // }
136 
137  //=== iy_aux part ===========================================================
138  Index auxPressure = -1,
139  auxTemperature = -1,
140  auxAbsSum = -1,
141  auxBackground = -1,
142  auxIy = -1,
143  auxOptDepth = -1;
144  ArrayOfIndex auxAbsSpecies(0), auxAbsIsp(0);
145  ArrayOfIndex auxVmrSpecies(0), auxVmrIsp(0);
146  ArrayOfIndex auxPartCont(0), auxPartContI(0);
147  ArrayOfIndex auxPartField(0), auxPartFieldI(0);
148  //
149  if( !iy_agenda_call1 )
150  { iy_aux.resize( 0 ); }
151  else
152  {
153  const Index naux = iy_aux_vars.nelem();
154  iy_aux.resize( naux );
155  //
156  for( Index i=0; i<naux; i++ )
157  {
158  if( iy_aux_vars[i] == "Pressure" )
159  { auxPressure = i; iy_aux[i].resize( 1, 1, 1, np ); }
160  else if( iy_aux_vars[i] == "Temperature" )
161  { auxTemperature = i; iy_aux[i].resize( 1, 1, 1, np ); }
162  else if( iy_aux_vars[i].substr(0,13) == "VMR, species " )
163  {
164  Index ispecies;
165  istringstream is(iy_aux_vars[i].substr(13,2));
166  is >> ispecies;
167  if( ispecies < 0 || ispecies>=abs_species.nelem() )
168  {
169  ostringstream os;
170  os << "You have selected VMR of species with index "
171  << ispecies << ".\nThis species does not exist!";
172  throw runtime_error( os.str() );
173  }
174  auxVmrSpecies.push_back(i);
175  auxVmrIsp.push_back(ispecies);
176  iy_aux[i].resize( 1, 1, 1, np );
177  }
178  else if( iy_aux_vars[i] == "Absorption, summed" )
179  { auxAbsSum = i; iy_aux[i].resize( nf, ns, ns, np ); }
180  else if( iy_aux_vars[i].substr(0,20) == "Absorption, species " )
181  {
182  Index ispecies;
183  istringstream is(iy_aux_vars[i].substr(20,2));
184  is >> ispecies;
185  if( ispecies < 0 || ispecies>=abs_species.nelem() )
186  {
187  ostringstream os;
188  os << "You have selected absorption species with index "
189  << ispecies << ".\nThis species does not exist!";
190  throw runtime_error( os.str() );
191  }
192  auxAbsSpecies.push_back(i);
193  const Index ihit = find_first( iaps, ispecies );
194  if( ihit >= 0 )
195  { auxAbsIsp.push_back( ihit ); }
196  else
197  {
198  iaps.push_back(ispecies);
199  auxAbsIsp.push_back( iaps.nelem()-1 );
200  }
201  iy_aux[i].resize( nf, ns, ns, np );
202  }
203  else if( iy_aux_vars[i] == "Radiative background" )
204  { auxBackground = i; iy_aux[i].resize( nf, 1, 1, 1 ); }
205  else if( iy_aux_vars[i] == "iy" && auxIy < 0 )
206  { auxIy = i; iy_aux[i].resize( nf, ns, 1, np ); }
207  else if( iy_aux_vars[i] == "Optical depth" )
208  { auxOptDepth = i; iy_aux[i].resize( nf, 1, 1, 1 ); }
209  else if( iy_aux_vars[i].substr(0,14) == "Mass content, " )
210  {
211  Index icont;
212  istringstream is(iy_aux_vars[i].substr(14,2));
213  is >> icont;
214  if( icont < 0 || icont>=particle_masses.ncols() )
215  {
216  ostringstream os;
217  os << "You have selected particle mass content category with "
218  << "index " << icont << ".\nThis category is not defined!";
219  throw runtime_error( os.str() );
220  }
221  auxPartCont.push_back(i);
222  auxPartContI.push_back(icont);
223  iy_aux[i].resize( 1, 1, 1, np );
224  }
225  else if( iy_aux_vars[i].substr(0,10) == "PND, type " )
226  {
227  Index ip;
228  istringstream is(iy_aux_vars[i].substr(10,2));
229  is >> ip;
230  if( ip < 0 || ip>=pnd_field.nbooks() )
231  {
232  ostringstream os;
233  os << "You have selected particle number density field with "
234  << "index " << ip << ".\nThis field is not defined!";
235  throw runtime_error( os.str() );
236  }
237  auxPartField.push_back(i);
238  auxPartFieldI.push_back(ip);
239  iy_aux[i].resize( 1, 1, 1, np );
240  }
241  else
242  {
243  ostringstream os;
244  os << "In *iy_aux_vars* you have included: \"" << iy_aux_vars[i]
245  << "\"\nThis choice is not recognised.";
246  throw runtime_error( os.str() );
247  }
248  }
249  }
250  //===========================================================================
251 
252  // Get atmospheric and RT quantities for each ppath point/step
253  //
254  Vector ppath_p, ppath_t;
255  Matrix ppath_vmr, ppath_pnd, ppath_wind, ppath_mag, ppath_f, ppath_nlte;
256  Matrix ppath_blackrad;
257  ArrayOfArrayOfPropagationMatrix abs_per_species, dummy_dppath_ext_dx;
258  ArrayOfPropagationMatrix ppath_ext, pnd_ext_mat;
259  Tensor4 trans_partial, trans_cumulat;
260  ArrayOfArrayOfStokesVector dummy_dppath_nlte_dx;
261  Tensor3 pnd_abs_vec;
262  ArrayOfStokesVector ppath_nlte_source;
263  Vector scalar_tau;
264  ArrayOfIndex clear2cloudy, lte;
265  ArrayOfMatrix dummy_ppath_dpnd_dx;
266  ArrayOfTensor4 dummy_dpnd_field_dx;
267  const Tensor4 nlte_field_empty(0,0,0,0);
268  //
269  Array<ArrayOfArrayOfSingleScatteringData> scat_data_single;
270  ArrayOfArrayOfIndex extmat_case;
271  //
272  if( np > 1 )
273  {
274  get_ppath_atmvars( ppath_p, ppath_t, ppath_nlte, ppath_vmr,
275  ppath_wind, ppath_mag,
276  ppath, atmosphere_dim, p_grid, t_field,
277  nlte_field_empty, vmr_field,
278  wind_u_field, wind_v_field, wind_w_field,
279  mag_u_field, mag_v_field, mag_w_field );
280 
281  get_ppath_f( ppath_f, ppath, f_grid, atmosphere_dim,
282  rte_alonglos_v, ppath_wind );
283 
284  get_ppath_pmat( ws, ppath_ext, ppath_nlte_source, lte, abs_per_species,
285  dummy_dppath_ext_dx, dummy_dppath_nlte_dx,
286  propmat_clearsky_agenda, ArrayOfRetrievalQuantity(0), ppath,
287  ppath_p, ppath_t, ppath_nlte, ppath_vmr, ppath_f,
288  ppath_mag, f_grid, stokes_dim, iaps );
289 
290  for( Index i=0; i<lte.nelem(); i++ )
291  {
292  if( lte[i] == 0 )
293  throw runtime_error( "FOS can so far only handle LTE conditions." );
294  }
295 
296  get_ppath_blackrad( ppath_blackrad, ppath, ppath_t, ppath_f );
297 
298  if( !cloudbox_on )
299  {
300  get_ppath_trans( trans_partial, extmat_case, trans_cumulat,
301  scalar_tau, ppath, ppath_ext, f_grid, stokes_dim );
302  }
303  else
304  {
305  get_ppath_cloudvars( clear2cloudy, ppath_pnd, dummy_ppath_dpnd_dx,
306  ppath, atmosphere_dim, cloudbox_limits,
307  pnd_field, dummy_dpnd_field_dx );
308  get_ppath_partopt( pnd_abs_vec, pnd_ext_mat, scat_data_single,
309  clear2cloudy, ppath_pnd, ppath, ppath_t,
310  stokes_dim, ppath_f, atmosphere_dim, scat_data,
311  verbosity );
312  get_ppath_trans2( trans_partial, extmat_case, trans_cumulat,
313  scalar_tau, ppath, ppath_ext, f_grid, stokes_dim,
314  clear2cloudy, pnd_ext_mat );
315  }
316  }
317  else // For cases totally outside the atmosphere,
318  { // set zero optical thickness and unit transmission
319  scalar_tau.resize( nf );
320  scalar_tau = 0;
321  trans_cumulat.resize( nf, ns, ns, np );
322  for( Index iv=0; iv<nf; iv++ )
323  { id_mat( trans_cumulat(iv,joker,joker,np-1) ); }
324  }
325 
326  // iy_transmission
327  //
328  Tensor3 iy_trans_new;
329  //
330  if( iy_agenda_call1 )
331  { iy_trans_new = trans_cumulat(joker,joker,joker,np-1); }
332  else
333  {
334  iy_transmission_mult( iy_trans_new, iy_transmission,
335  trans_cumulat(joker,joker,joker,np-1) ); }
336 
337  // Radiative background
338  //
339  {
340  Agenda iy_cbox_agenda;
341  const Index iy_id = 0;
342  get_iy_of_background( ws, iy, diy_dx,
343  iy_trans_new, iy_id, jacobian_do, ppath, rte_pos2,
344  atmosphere_dim, t_field, z_field, vmr_field,
345  cloudbox_on, stokes_dim, f_grid, iy_unit,
346  iy_main_agenda, iy_space_agenda, iy_surface_agenda,
347  iy_cbox_agenda, verbosity );
348  }
349 
350  //=== iy_aux part ===========================================================
351  // Fill parts of iy_aux that are defined even for np=1.
352  // Radiative background
353  if( auxBackground >= 0 )
354  { iy_aux[auxBackground](0,0,0,0) = (Numeric)min( (Index)2,
355  ppath_what_background(ppath)-1); }
356  // Radiance
357  if( auxIy >= 0 )
358  { iy_aux[auxIy](joker,joker,0,np-1) = iy; }
359  // Transmission, total
360  if( auxOptDepth >= 0 )
361  { iy_aux[auxOptDepth](joker,0,0,0) = scalar_tau; }
362  //===========================================================================
363 
364 
365  // Do RT calculations
366  //
367  if( np > 1 )
368  {
369  //=== iy_aux part =======================================================
370  // iy_aux for point np-1:
371  // Pressure
372  if( auxPressure >= 0 )
373  { iy_aux[auxPressure](0,0,0,np-1) = ppath_p[np-1]; }
374  // Temperature
375  if( auxTemperature >= 0 )
376  { iy_aux[auxTemperature](0,0,0,np-1) = ppath_t[np-1]; }
377  // VMR
378  for( Index j=0; j<auxVmrSpecies.nelem(); j++ )
379  { iy_aux[auxVmrSpecies[j]](0,0,0,np-1) = ppath_vmr(auxVmrIsp[j],np-1); }
380  // Absorption
381  if( auxAbsSum >= 0 )
382  { for( Index iv=0; iv<nf; iv++ ) {
383  for( Index is1=0; is1<ns; is1++ ){
384  for( Index is2=0; is2<ns; is2++ ){
385  iy_aux[auxAbsSum](iv,is1,is2,np-1) =
386  ppath_ext[np-1](iv,is1,is2); } } } }
387  for( Index j=0; j<auxAbsSpecies.nelem(); j++ )
388  { for( Index iv=0; iv<nf; iv++ ) {
389  for( Index is1=0; is1<stokes_dim; is1++ ){
390  for( Index is2=0; is2<stokes_dim; is2++ ){
391  iy_aux[auxAbsSpecies[j]](iv,is1,is2,np-1) =
392  abs_per_species[np-1][auxAbsIsp[j]](iv,is1,is2); } } } }
393  // Particle properties
394  if( cloudbox_on )
395  {
396  // Particle mass content
397  for( Index j=0; j<auxPartCont.nelem(); j++ )
398  { iy_aux[auxPartCont[j]](0,0,0,np-1) = ppath_pnd(joker,np-1) *
399  particle_masses(joker,auxPartContI[j]); }
400  // Particle number density
401  for( Index j=0; j<auxPartField.nelem(); j++ )
402  { iy_aux[auxPartField[j]](0,0,0,np-1) =
403  ppath_pnd(auxPartFieldI[j],np-1); }
404  }
405  // Radiance for this point is handled above
406  //=======================================================================
407 
408 
409  // Scattering source term at ip (0) and ip+1 (1)
410  // (If any particles at ip=np-1, this is ignored. Could happen if
411  // particles at surface level.)
412  Matrix ssource0(nf,ns,0), ssource1(nf,ns);
413 
414  // Dummy variables for non-LTE
415  const bool nonlte = false;
416  const Matrix sourcebar_dummy(0,0);
417  const Tensor3 extbar_dummy(0,0,0);
418 
419  // Loop ppath steps
420  for( Index ip=np-2; ip>=0; ip-- )
421  {
422  // Path step average of emission source function: Bbar
423  Vector bbar(nf);
424  for( Index iv=0; iv<nf; iv++ )
425  { bbar[iv] = 0.5 * ( ppath_blackrad(iv,ip) +
426  ppath_blackrad(iv,ip+1) ); }
427 
428  // Check if any particles to consider
429  bool any_particles = clear2cloudy[ip] >= 0 ||
430  clear2cloudy[ip+1] >= 0;
431 
432  // -----------------------------------------------------------------
433  // i = N (only absorption/emission)
434  // -----------------------------------------------------------------
435  if( fos_i == fos_n )
436  {
437  // No particle absorption to consider
438  if( !any_particles )
439  {
440  emission_rtstep( iy, stokes_dim, bbar, extmat_case[ip],
441  trans_partial(joker,joker,joker,ip),
442  nonlte, extbar_dummy, sourcebar_dummy );
443  }
444 
445  else // We want to include particle absorption, but not
446  { // extinction. trans_partial is then not valid.
447  Tensor3 t(nf,ns,ns);
448  ArrayOfIndex extmat_cas2(nf);
449  //
450  for( Index iv=0; iv<nf; iv++ )
451  {
452  // Particle absorption
453  //
454  Matrix pabs_mat(ns,ns,0);
455  //
456  if( clear2cloudy[ip] >= 0 )
457  { ext_matFromabs_vec( pabs_mat, pnd_abs_vec(iv,joker,
458  clear2cloudy[ip]), stokes_dim ); }
459  if( clear2cloudy[ip+1] >= 0 )
460  { ext_matFromabs_vec( pabs_mat, pnd_abs_vec(iv,joker,
461  clear2cloudy[ip+1]), stokes_dim ); }
462 
463  // Transmission of step
464  Matrix ext_mat(stokes_dim,stokes_dim);
465  for( Index is1=0; is1<stokes_dim; is1++ ) {
466  for( Index is2=0; is2<stokes_dim; is2++ ) {
467  ext_mat(is1,is2) = 0.5 * ( pabs_mat(is1,is2) +
468  ppath_ext[ip](iv,is1,is2) +
469  ppath_ext[ip+1](iv,is1,is2) );
470  } }
471  //
472  extmat_cas2[iv] = 0;
473  ext2trans( t(iv,joker,joker), extmat_cas2[iv],
474  ext_mat, ppath.lstep[ip] );
475  }
476 
477  // Perform RT
478  emission_rtstep( iy, stokes_dim, bbar, extmat_cas2, t,
479  nonlte, extbar_dummy, sourcebar_dummy );
480  }
481  }
482 
483 
484  // -----------------------------------------------------------------
485  // i < N
486  // -----------------------------------------------------------------
487  else
488  {
489  // Shift scattering source term (new 1 is old 0)
490  ssource1 = ssource0;
491 
492  // Clear-sky path step
493  if( !any_particles )
494  {
495  // Perform RT
496  emission_rtstep( iy, stokes_dim, bbar, extmat_case[ip],
497  trans_partial(joker,joker,joker,ip),
498  nonlte, extbar_dummy, sourcebar_dummy );
499 
500  // Scattering source term at ip is zero:
501  ssource0 = 0;
502  }
503 
504  // Include scattering
505  else
506  {
507  // Determine scattering source term at ip
508  if( clear2cloudy[ip] < 0 )
509  { ssource0 = 0; }
510  else
511  {
512  // Present position
513  // (Note that the Ppath positions (ppath.pos) for 1D have
514  // one column more than expected by most functions. Only
515  // the first atmosphere_dim values shall be copied.)
516  Vector pos = ppath.pos(ip,Range(0,atmosphere_dim));
517 
518  // Determine incoming radiation
519  //
520  const Index nin = fos_scatint_angles.nrows();
521  Tensor3 Y(nin,nf,ns);
522  {
523  // Do RT
524  const Index nza = fos_iyin_za_angles.nelem();
525  Tensor3 Y1(nza,nf,ns);
526  //
527  for( Index i=0; i<nza; i++ )
528  {
529  // LOS
530  Vector los( 1, fos_iyin_za_angles[i] );
531 
532  // Call recursively, with fos_i increased
533  //
534  Matrix iyl;
535  ArrayOfTensor4 iy_auxl;
536  Ppath ppathl;
537  ArrayOfTensor3 diy_dxl;
538  //
539  fos( ws, iyl, iy_auxl, ppathl, diy_dxl,
540  stokes_dim, f_grid, atmosphere_dim,
541  p_grid, z_field, t_field, vmr_field,
542  abs_species, wind_u_field, wind_v_field,
543  wind_w_field, mag_u_field, mag_v_field,
544  mag_w_field, cloudbox_on, cloudbox_limits,
545  pnd_field,
546  scat_data,
547  particle_masses, iy_unit, iy_aux_vars,
548  jacobian_do, ppath_agenda,
549  propmat_clearsky_agenda, iy_main_agenda,
550  iy_space_agenda, iy_surface_agenda, 0,
551  iy_trans_new, pos, los, rte_pos2,
552  rte_alonglos_v, ppath_lmax, ppath_lraytrace,
553  fos_scatint_angles, fos_iyin_za_angles,
554  fos_za_interporder, fos_n, fos_i+1,
555  verbosity );
556 
557  Y1(i,joker,joker) = iyl;
558  }
559 
560  // Zenith angle interpolation of Y
561  ArrayOfGridPosPoly gp( nin );
562  gridpos_poly( gp, fos_iyin_za_angles,
563  fos_scatint_angles(joker,0),
564  fos_za_interporder );
565  Matrix itw( nin, fos_za_interporder+1 );
566  interpweights( itw, gp );
567  //
568  for( Index iv=0; iv<nf; iv++ )
569  {
570  for( Index is1=0; is1<stokes_dim; is1++ )
571  {
572  interp( Y(joker,iv,is1), itw,
573  Y1(joker,iv,is1), gp );
574  }
575  }
576  }
577 
578  // Direction of outgoing scattered radiation (which is
579  // reversed to LOS). Note that this outlos is only used
580  // for extracting scattering properties.
581  Vector outlos;
582  mirror_los( outlos, ppath.los(ip,joker), atmosphere_dim );
583 
584  // Determine phase matrix
585  Tensor4 P( nin, nf, stokes_dim, stokes_dim );
586  Matrix P1( stokes_dim, stokes_dim );
587  //
588  for( Index ii=0; ii<nin; ii++ )
589  {
590  for( Index iv=0; iv<nf; iv++ )
591  {
592  pha_mat_singleCalc( P1, outlos[0], outlos[1],
593  fos_scatint_angles(ii,0),
594  fos_scatint_angles(ii,1),
595  scat_data_single[iv], stokes_dim,
596  ppath_pnd(joker,ip),
597  ppath_t[ip], verbosity );
598  P(ii,iv,joker,joker) = P1;
599  }
600  }
601 
602 
603  // Scattering source term
604  ssource0 = 0.0;
605  for( Index iv=0; iv<nf; iv++ )
606  {
607  Vector sp(stokes_dim);
608  for( Index ii=0; ii<nin; ii++ )
609  {
610  mult( sp, P(ii,iv,joker,joker),
611  Y(ii,iv,joker) );
612  ssource0(iv,joker) += sp;
613  }
614  }
615  ssource0 *= 4*PI/(Numeric)nin;
616  }
617 
618  // RT of ppath step
619  for( Index iv=0; iv<nf; iv++ )
620  {
621  // Calculate average of absorption, extinction etc.
622  Matrix ext_mat( stokes_dim, stokes_dim );
623  Vector abs_vec( stokes_dim );
624  Vector sbar( stokes_dim, 0 );
625  //
626  // Contribution from abs_species
627  for( Index is1=0; is1<stokes_dim; is1++ )
628  {
629  abs_vec[is1] = 0.5 * (
630  ppath_ext[ip](iv,is1,0) +
631  ppath_ext[ip+1](iv,is1,0) );
632  for( Index is2=0; is2<stokes_dim; is2++ )
633  {
634  ext_mat(is1,is2) = 0.5 * (
635  ppath_ext[ip](iv,is1,is2) +
636  ppath_ext[ip+1](iv,is1,is2) );
637  }
638  }
639  // Particle contribution
640  if( clear2cloudy[ip] >= 0 )
641  {
642  for( Index is1=0; is1<stokes_dim; is1++ )
643  {
644  sbar[is1] += 0.5 * ssource0(iv,is1);
645  abs_vec[is1] += 0.5 * (
646  pnd_abs_vec(iv,is1,clear2cloudy[ip]) );
647  for( Index is2=0; is2<stokes_dim; is2++ )
648  {
649  ext_mat(is1,is2) += 0.5 * (
650  pnd_ext_mat[clear2cloudy[ip]](iv,is1,is2) );
651  }
652  }
653  }
654  if( clear2cloudy[ip+1] >= 0 )
655  {
656  for( Index is1=0; is1<stokes_dim; is1++ )
657  {
658  sbar[is1] += 0.5 * ssource1(iv,is1);
659  abs_vec[is1] += 0.5 * (
660  pnd_abs_vec(iv,is1,clear2cloudy[ip+1]) );
661  for( Index is2=0; is2<stokes_dim; is2++ )
662  {
663  ext_mat(is1,is2) += 0.5 * (
664  pnd_ext_mat[clear2cloudy[ip+1]](iv,is1,is2));
665  }
666  }
667  }
668 
669  // Perform RT
670  //
671  Matrix trans_mat = trans_partial(iv,joker,joker,ip);
672  rte_step_doit( iy(iv,joker), trans_mat, ext_mat, abs_vec,
673  sbar, ppath.lstep[ip], bbar[iv], true );
674  }
675  }
676  }
677 
678  //=== iy_aux part ===================================================
679  // Pressure
680  if( auxPressure >= 0 )
681  { iy_aux[auxPressure](0,0,0,ip) = ppath_p[ip]; }
682  // Temperature
683  if( auxTemperature >= 0 )
684  { iy_aux[auxTemperature](0,0,0,ip) = ppath_t[ip]; }
685  // VMR
686  for( Index j=0; j<auxVmrSpecies.nelem(); j++ )
687  { iy_aux[auxVmrSpecies[j]](0,0,0,ip) = ppath_vmr(auxVmrIsp[j],ip);}
688  // Absorption
689  if( auxAbsSum >= 0 )
690  { for( Index iv=0; iv<nf; iv++ ) {
691  for( Index is1=0; is1<ns; is1++ ){
692  for( Index is2=0; is2<ns; is2++ ){
693  iy_aux[auxAbsSum](iv,is1,is2,ip) =
694  ppath_ext[ip](iv,is1,is2); } } } }
695  for( Index j=0; j<auxAbsSpecies.nelem(); j++ )
696  { for( Index iv=0; iv<nf; iv++ ) {
697  for( Index is1=0; is1<stokes_dim; is1++ ){
698  for( Index is2=0; is2<stokes_dim; is2++ ){
699  iy_aux[auxAbsSpecies[j]](iv,is1,is2,ip) =
700  abs_per_species[ip][auxAbsIsp[j]](iv,is1,is2); } } } }
701  // Particle properties
702  if( cloudbox_on )
703  {
704  // Particle mass content
705  for( Index j=0; j<auxPartCont.nelem(); j++ )
706  { iy_aux[auxPartCont[j]](0,0,0,ip) = ppath_pnd(joker,ip) *
707  particle_masses(joker,auxPartContI[j]); }
708  // Particle number density
709  for( Index j=0; j<auxPartField.nelem(); j++ )
710  { iy_aux[auxPartField[j]](0,0,0,ip) =
711  ppath_pnd(auxPartFieldI[j],ip); }
712  }
713  // Radiance
714  if( auxIy >= 0 )
715  { iy_aux[auxIy](joker,joker,0,ip) = iy; }
716  //===================================================================
717  }
718  } // if np>1
719 
720 
721  // Unit conversions
722  if( iy_agenda_call1 )
723  {
724  // Determine refractive index to use for the n2 radiance law
725  Numeric n = 1.0; // First guess is that sensor is in space
726  //
727  if( ppath.end_lstep == 0 ) // If true, sensor inside the atmosphere
728  { n = ppath.nreal[np-1]; }
729 
730  // Polarisation index variable
731  ArrayOfIndex i_pol(stokes_dim);
732  for( Index is=0; is<stokes_dim; is++ )
733  { i_pol[is] = is + 1; }
734 
735  // iy
736  apply_iy_unit( iy, iy_unit, f_grid, n, i_pol );
737 
738  // iy_aux
739  for( Index q=0; q<iy_aux.nelem(); q++ )
740  {
741  if( iy_aux_vars[q] == "iy")
742  {
743  for( Index ip=0; ip<ppath.np; ip++ )
744  { apply_iy_unit( iy_aux[q](joker,joker,0,ip), iy_unit, f_grid,
745  ppath.nreal[ip], i_pol ); }
746  }
747  }
748  }
749 }
750 */
751 
752 /*===========================================================================
753  === The functions (in alphabetical order)
754  ===========================================================================*/
755 
756 /* Workspace method: Doxygen documentation will be auto-generated */
757 /*
758 void iyFOS(
759  Workspace& ws,
760  Matrix& iy,
761  ArrayOfTensor4& iy_aux,
762  Ppath& ppath,
763  ArrayOfTensor3& diy_dx,
764  const Index& stokes_dim,
765  const Vector& f_grid,
766  const Index& atmosphere_dim,
767  const Vector& p_grid,
768  const Tensor3& z_field,
769  const Tensor3& t_field,
770  const Tensor4& vmr_field,
771  const ArrayOfArrayOfSpeciesTag& abs_species,
772  const Tensor3& wind_u_field,
773  const Tensor3& wind_v_field,
774  const Tensor3& wind_w_field,
775  const Tensor3& mag_u_field,
776  const Tensor3& mag_v_field,
777  const Tensor3& mag_w_field,
778  const Index& cloudbox_on,
779  const ArrayOfIndex& cloudbox_limits,
780  const Tensor4& pnd_field,
781  const ArrayOfArrayOfSingleScatteringData& scat_data,
782  const Matrix& particle_masses,
783  const String& iy_unit,
784  const ArrayOfString& iy_aux_vars,
785  const Index& jacobian_do,
786  const Agenda& ppath_agenda,
787  const Agenda& propmat_clearsky_agenda,
788  const Agenda& iy_main_agenda,
789  const Agenda& iy_space_agenda,
790  const Agenda& iy_surface_agenda,
791  const Index& iy_agenda_call1,
792  const Tensor3& iy_transmission,
793  const Vector& rte_pos,
794  const Vector& rte_los,
795  const Vector& rte_pos2,
796  const Numeric& rte_alonglos_v,
797  const Numeric& ppath_lmax,
798  const Numeric& ppath_lraytrace,
799  const Matrix& fos_scatint_angles,
800  const Vector& fos_iyin_za_angles,
801  const Index& fos_za_interporder,
802  const Index& fos_n,
803  const Verbosity& verbosity )
804 {
805  // Input checks
806  if( jacobian_do )
807  throw runtime_error(
808  "This method does not yet provide any jacobians (jacobian_do must be 0)" );
809  if( fos_scatint_angles.ncols() != 2 )
810  throw runtime_error( "The WSV *fos_scatint_angles* must have two columns." );
811  if( min(fos_scatint_angles(joker,0))<0 ||
812  max(fos_scatint_angles(joker,0))>180 )
813  throw runtime_error(
814  "The zenith angles in *fos_scatint_angles* shall be inside [0,180]." );
815  if( min(fos_scatint_angles(joker,1))<-180 ||
816  max(fos_scatint_angles(joker,1))>180 )
817  throw runtime_error(
818  "The azimuth angles in *fos_scatint_angles* shall be inside [-180,180]." );
819  if( min(fos_iyin_za_angles)<0 || max(fos_iyin_za_angles)>180 )
820  throw runtime_error(
821  "The zenith angles in *fos_iyin_za_angles* shall be inside [0,180]." );
822  if( fos_iyin_za_angles[0] != 0 )
823  throw runtime_error( "The first value in *fos_iyin_za_angles* must be 0." );
824  if( last(fos_iyin_za_angles) != 180 )
825  throw runtime_error( "The last value in *fos_iyin_za_angles* must be 180." );
826  if( fos_za_interporder < 1 )
827  throw runtime_error( "The argument *fos_za_interporder* must be >= 1." );
828  if( fos_iyin_za_angles.nelem() <= fos_za_interporder )
829  throw runtime_error( "The length of *fos_iyin_za_angles* must at least "
830  "be *fos_za_interporder* + 1." );
831  if( fos_n < 0 )
832  throw runtime_error( "The argument *fos_n* must be >= 0." );
833 
834  // Switch to order 0 if not primary call
835  // (This happens after surface reflection. If fos_n used (and >=1), new
836  // surface relections are created ..., and recursion never ends.)
837  Index n = fos_n;
838  if( !iy_agenda_call1 )
839  { n = 0; }
840 
841  fos( ws, iy, iy_aux, ppath, diy_dx, stokes_dim, f_grid, atmosphere_dim,
842  p_grid, z_field, t_field, vmr_field, abs_species, wind_u_field,
843  wind_v_field, wind_w_field, mag_u_field, mag_v_field, mag_w_field,
844  cloudbox_on, cloudbox_limits, pnd_field,
845  scat_data, particle_masses, iy_unit, iy_aux_vars, jacobian_do,
846  ppath_agenda, propmat_clearsky_agenda,
847  iy_main_agenda, iy_space_agenda, iy_surface_agenda, iy_agenda_call1,
848  iy_transmission, rte_pos, rte_los, rte_pos2, rte_alonglos_v,
849  ppath_lmax, ppath_lraytrace, fos_scatint_angles, fos_iyin_za_angles,
850  fos_za_interporder, n, 0, verbosity );
851 }
852 */
853 
854 /* Workspace method: Doxygen documentation will be auto-generated */
855 void iyHybrid(Workspace& ws,
856  Matrix& iy,
857  ArrayOfMatrix& iy_aux,
858  ArrayOfTensor3& diy_dx,
859  Vector& ppvar_p,
860  Vector& ppvar_t,
861  EnergyLevelMap& ppvar_nlte,
862  Matrix& ppvar_vmr,
863  Matrix& ppvar_wind,
864  Matrix& ppvar_mag,
865  Matrix& ppvar_pnd,
866  Matrix& ppvar_f,
867  Tensor3& ppvar_iy,
868  Tensor4& ppvar_trans_cumulat,
869  const Index& iy_id,
870  const Index& stokes_dim,
871  const Vector& f_grid,
872  const Index& atmosphere_dim,
873  const Vector& p_grid,
874  const Tensor3& t_field,
875  const EnergyLevelMap& nlte_field,
876  const Tensor4& vmr_field,
877  const ArrayOfArrayOfSpeciesTag& abs_species,
878  const Tensor3& wind_u_field,
879  const Tensor3& wind_v_field,
880  const Tensor3& wind_w_field,
881  const Tensor3& mag_u_field,
882  const Tensor3& mag_v_field,
883  const Tensor3& mag_w_field,
884  const Index& cloudbox_on,
885  const ArrayOfIndex& cloudbox_limits,
886  const Tensor4& pnd_field,
887  const ArrayOfTensor4& dpnd_field_dx,
888  const ArrayOfString& scat_species,
889  const ArrayOfArrayOfSingleScatteringData& scat_data,
890  const String& iy_unit,
891  const ArrayOfString& iy_aux_vars,
892  const Index& jacobian_do,
893  const ArrayOfRetrievalQuantity& jacobian_quantities,
894  const Agenda& propmat_clearsky_agenda,
895  const Agenda& water_p_eq_agenda,
896  const Agenda& iy_main_agenda,
897  const Agenda& iy_space_agenda,
898  const Agenda& iy_surface_agenda,
899  const Agenda& iy_cloudbox_agenda,
900  const Index& iy_agenda_call1,
901  const Tensor3& iy_transmission,
902  const Ppath& ppath,
903  const Vector& rte_pos2,
904  const Numeric& rte_alonglos_v,
905  const Tensor3& surface_props_data,
906  const Tensor7& cloudbox_field,
907  const Vector& za_grid,
908  const Index& Naa,
909  const Index& t_interp_order,
910  const Verbosity& verbosity) {
911  // If cloudbox off, switch to use clearsky method
912  if (!cloudbox_on) {
913  Tensor4 dummy;
914  iyEmissionStandard(ws,
915  iy,
916  iy_aux,
917  diy_dx,
918  ppvar_p,
919  ppvar_t,
920  ppvar_nlte,
921  ppvar_vmr,
922  ppvar_wind,
923  ppvar_mag,
924  ppvar_f,
925  ppvar_iy,
926  ppvar_trans_cumulat,
927  dummy,
928  iy_id,
929  stokes_dim,
930  f_grid,
931  atmosphere_dim,
932  p_grid,
933  t_field,
934  nlte_field,
935  vmr_field,
936  abs_species,
937  wind_u_field,
938  wind_v_field,
939  wind_w_field,
940  mag_u_field,
941  mag_v_field,
942  mag_w_field,
943  cloudbox_on,
944  iy_unit,
945  iy_aux_vars,
946  jacobian_do,
947  jacobian_quantities,
948  ppath,
949  rte_pos2,
950  propmat_clearsky_agenda,
951  water_p_eq_agenda,
952  iy_main_agenda,
953  iy_space_agenda,
954  iy_surface_agenda,
955  iy_cloudbox_agenda,
956  iy_agenda_call1,
957  iy_transmission,
958  rte_alonglos_v,
959  surface_props_data,
960  verbosity);
961  return;
962  }
963 
964  // Some basic sizes
965  const Index nf = f_grid.nelem();
966  const Index ns = stokes_dim;
967  const Index np = ppath.np;
968 
969  // Radiative background index
970  const Index rbi = ppath_what_background(ppath);
971 
972  // Throw error if unsupported features are requested
973  if (atmosphere_dim != 1)
974  throw runtime_error(
975  "With cloudbox on, this method handles only 1D calculations.");
976  if (Naa < 3) throw runtime_error("Naa must be > 2.");
977  if (jacobian_do) {
978  if (dpnd_field_dx.nelem() != jacobian_quantities.nelem())
979  throw runtime_error(
980  "*dpnd_field_dx* not properly initialized:\n"
981  "Number of elements in dpnd_field_dx must be equal number of jacobian"
982  " quantities.\n(Note: jacobians have to be defined BEFORE *pnd_field*"
983  " is calculated/set.");
984  }
985  if (rbi < 1 || rbi > 9)
986  throw runtime_error(
987  "ppath.background is invalid. Check your "
988  "calculation of *ppath*?");
989  if (rbi == 3 || rbi == 4)
990  throw runtime_error(
991  "The propagation path ends inside or at boundary of "
992  "the cloudbox.\nFor this method, *ppath* must be "
993  "calculated in this way:\n ppathCalc( cloudbox_on = 0 ).");
994  // iy_aux_vars checked below
995  // Checks of i_field
996  if (cloudbox_field.ncols() != stokes_dim)
997  throw runtime_error(
998  "Obtained *cloudbox_field* number of Stokes elements inconsistent with "
999  "*stokes_dim*.");
1000  if (cloudbox_field.nrows() != 1)
1001  throw runtime_error(
1002  "Obtained *cloudbox_field* has wrong number of azimuth angles.");
1003  if (cloudbox_field.npages() != za_grid.nelem())
1004  throw runtime_error(
1005  "Obtained *cloudbox_field* number of zenith angles inconsistent with "
1006  "*za_grid*.");
1007  if (cloudbox_field.nbooks() != 1)
1008  throw runtime_error(
1009  "Obtained *cloudbox_field* has wrong number of longitude points.");
1010  if (cloudbox_field.nshelves() != 1)
1011  throw runtime_error(
1012  "Obtained *cloudbox_field* has wrong number of latitude points.");
1013  if (cloudbox_field.nvitrines() != cloudbox_limits[1] - cloudbox_limits[0] + 1)
1014  throw runtime_error(
1015  "Obtained *cloudbox_field* number of pressure points inconsistent with "
1016  "*cloudbox_limits*.");
1017  if (cloudbox_field.nlibraries() != nf)
1018  throw runtime_error(
1019  "Obtained *cloudbox_field* number of frequency points inconsistent with "
1020  "*f_grid*.");
1021 
1022  // Init Jacobian quantities
1023  Index j_analytical_do = 0;
1024  if (jacobian_do) {
1025  FOR_ANALYTICAL_JACOBIANS_DO(j_analytical_do = 1;)
1026  }
1027  //
1028  const Index nq = j_analytical_do ? jacobian_quantities.nelem() : 0;
1029  ArrayOfTensor3 diy_dpath(nq);
1030  ArrayOfIndex jac_species_i(nq), jac_scat_i(nq), jac_is_t(nq), jac_wind_i(nq);
1031  ArrayOfIndex jac_mag_i(nq), jac_other(nq);
1032 
1033  if (j_analytical_do) {
1034  rtmethods_jacobian_init(jac_species_i,
1035  jac_scat_i,
1036  jac_is_t,
1037  jac_wind_i,
1038  jac_mag_i,
1039  jac_other,
1040  diy_dx,
1041  diy_dpath,
1042  ns,
1043  nf,
1044  np,
1045  nq,
1046  abs_species,
1047  cloudbox_on,
1048  scat_species,
1049  dpnd_field_dx,
1050  jacobian_quantities,
1051  iy_agenda_call1);
1052  }
1053 
1054  // Init iy_aux and fill where possible
1055  const Index naux = iy_aux_vars.nelem();
1056  iy_aux.resize(naux);
1057  //
1058  for (Index i = 0; i < naux; i++) {
1059  iy_aux[i].resize(nf, ns);
1060 
1061  if (iy_aux_vars[i] == "Optical depth") {
1062  } // Filled below
1063  else if (iy_aux_vars[i] == "Radiative background") {
1064  iy_aux[i] = (Numeric)min((Index)2, rbi - 1);
1065  } else {
1066  ostringstream os;
1067  os << "The only allowed strings in *iy_aux_vars* are:\n"
1068  << " \"Radiative background\"\n"
1069  << " \"Optical depth\"\n"
1070  << "but you have selected: \"" << iy_aux_vars[i] << "\"";
1071  throw runtime_error(os.str());
1072  }
1073  }
1074 
1075  // Get atmospheric and radiative variables along the propagation path
1076  //
1077  ppvar_trans_cumulat.resize(np, nf, ns, ns);
1078  Tensor3 J(np, nf, ns);
1079  Tensor4 trans_partial(np, nf, ns, ns);
1080  Tensor5 dtrans_partial_dx_above(np, nq, nf, ns, ns);
1081  Tensor5 dtrans_partial_dx_below(np, nq, nf, ns, ns);
1082  Tensor4 dJ_dx(np, nq, nf, ns);
1083  ArrayOfIndex clear2cloudy;
1084  //
1085  if (np == 1 && rbi == 1) // i.e. ppath is totally outside the atmosphere:
1086  {
1087  ppvar_p.resize(0);
1088  ppvar_t.resize(0);
1089  ppvar_vmr.resize(0, 0);
1090  ppvar_wind.resize(0, 0);
1091  ppvar_mag.resize(0, 0);
1092  ppvar_pnd.resize(0, 0);
1093  ppvar_f.resize(0, 0);
1094  ppvar_iy.resize(0, 0, 0);
1095  } else {
1096  // ppvar_iy
1097  ppvar_iy.resize(nf, ns, np);
1098 
1099  // Basic atmospheric variables
1100  get_ppath_atmvars(ppvar_p,
1101  ppvar_t,
1102  ppvar_nlte,
1103  ppvar_vmr,
1104  ppvar_wind,
1105  ppvar_mag,
1106  ppath,
1107  atmosphere_dim,
1108  p_grid,
1109  t_field,
1110  nlte_field,
1111  vmr_field,
1112  wind_u_field,
1113  wind_v_field,
1114  wind_w_field,
1115  mag_u_field,
1116  mag_v_field,
1117  mag_w_field);
1118 
1119  get_ppath_f(
1120  ppvar_f, ppath, f_grid, atmosphere_dim, rte_alonglos_v, ppvar_wind);
1121 
1122  // here, the cloudbox is on, ie we don't need to check and branch this here
1123  // anymore.
1124  ArrayOfMatrix ppvar_dpnd_dx;
1125  //
1126  get_ppath_cloudvars(clear2cloudy,
1127  ppvar_pnd,
1128  ppvar_dpnd_dx,
1129  ppath,
1130  atmosphere_dim,
1131  cloudbox_limits,
1132  pnd_field,
1133  dpnd_field_dx);
1134 
1135  // Size radiative variables always used
1136  Vector B(nf);
1137  PropagationMatrix K_this(nf, ns), K_past(nf, ns), Kp(nf, ns);
1138  StokesVector a(nf, ns), S(nf, ns), Sp(nf, ns);
1139  ArrayOfIndex lte(np);
1140 
1141  // Init variables only used if analytical jacobians done
1142  Vector dB_dT(0);
1143  ArrayOfPropagationMatrix dK_this_dx(nq), dK_past_dx(nq), dKp_dx(nq);
1144  ArrayOfStokesVector da_dx(nq), dS_dx(nq), dSp_dx(nq);
1145  //
1146  if (j_analytical_do) {
1147  dB_dT.resize(nf);
1148  FOR_ANALYTICAL_JACOBIANS_DO(dK_this_dx[iq] = PropagationMatrix(nf, ns);
1149  dK_past_dx[iq] = PropagationMatrix(nf, ns);
1150  dKp_dx[iq] = PropagationMatrix(nf, ns);
1151  da_dx[iq] = StokesVector(nf, ns);
1152  dS_dx[iq] = StokesVector(nf, ns);
1153  dSp_dx[iq] = StokesVector(nf, ns);)
1154  }
1155 
1156  // Loop ppath points and determine radiative properties
1157  for (Index ip = 0; ip < np; ip++) {
1158  bool temperature_jacobian =
1159  jacobian_do && do_temperature_jacobian(jacobian_quantities);
1160 
1161  get_stepwise_blackbody_radiation(
1162  B, dB_dT, ppvar_f(joker, ip), ppvar_t[ip], temperature_jacobian);
1163 
1164  get_stepwise_clearsky_propmat(ws,
1165  K_this,
1166  S,
1167  lte[ip],
1168  dK_this_dx,
1169  dS_dx,
1170  propmat_clearsky_agenda,
1171  jacobian_quantities,
1172  ppvar_f(joker, ip),
1173  ppvar_mag(joker, ip),
1174  ppath.los(ip, joker),
1175  ppvar_nlte[ip],
1176  ppvar_vmr(joker, ip),
1177  ppvar_t[ip],
1178  ppvar_p[ip],
1179  jac_species_i,
1180  j_analytical_do);
1181 
1182  if (j_analytical_do) {
1183  adapt_stepwise_partial_derivatives(dK_this_dx,
1184  dS_dx,
1185  jacobian_quantities,
1186  ppvar_f(joker, ip),
1187  ppath.los(ip, joker),
1188  ppvar_vmr(joker, ip),
1189  ppvar_t[ip],
1190  ppvar_p[ip],
1191  jac_species_i,
1192  jac_wind_i,
1193  lte[ip],
1194  atmosphere_dim,
1195  j_analytical_do);
1196  }
1197 
1198  if (clear2cloudy[ip] + 1) {
1199  get_stepwise_scattersky_propmat(a,
1200  Kp,
1201  da_dx,
1202  dKp_dx,
1203  jacobian_quantities,
1204  ppvar_pnd(joker, Range(ip, 1)),
1205  ppvar_dpnd_dx,
1206  ip,
1207  scat_data,
1208  ppath.los(ip, joker),
1209  ppvar_t[Range(ip, 1)],
1210  atmosphere_dim,
1211  jacobian_do);
1212  a += K_this;
1213  K_this += Kp;
1214 
1215  if (j_analytical_do) {
1216  FOR_ANALYTICAL_JACOBIANS_DO(da_dx[iq] += dK_this_dx[iq];
1217  dK_this_dx[iq] += dKp_dx[iq];)
1218  }
1219 
1220  Vector aa_grid;
1221  nlinspace(aa_grid, 0, 360, Naa);
1222  //
1223  get_stepwise_scattersky_source(Sp,
1224  dSp_dx,
1225  jacobian_quantities,
1226  ppvar_pnd(joker, ip),
1227  ppvar_dpnd_dx,
1228  ip,
1229  scat_data,
1230  cloudbox_field,
1231  za_grid,
1232  aa_grid,
1233  ppath.los(Range(ip, 1), joker),
1234  ppath.gp_p[ip],
1235  ppvar_t[Range(ip, 1)],
1236  atmosphere_dim,
1237  jacobian_do,
1238  t_interp_order);
1239  S += Sp;
1240 
1241  if (j_analytical_do) {
1242  FOR_ANALYTICAL_JACOBIANS_DO(dS_dx[iq] += dSp_dx[iq];)
1243  }
1244  }
1245 
1246  else // no particles present at this level
1247  {
1248  a = K_this;
1249  if (j_analytical_do) {
1250  FOR_ANALYTICAL_JACOBIANS_DO(da_dx[iq] = dK_this_dx[iq];)
1251  }
1252  }
1253 
1254  get_stepwise_transmission_matrix(
1255  ppvar_trans_cumulat(ip, joker, joker, joker),
1256  trans_partial(ip, joker, joker, joker),
1257  dtrans_partial_dx_above(ip, joker, joker, joker, joker),
1258  dtrans_partial_dx_below(ip, joker, joker, joker, joker),
1259  (ip > 0) ? ppvar_trans_cumulat(ip - 1, joker, joker, joker)
1260  : Tensor3(0, 0, 0),
1261  K_past,
1262  K_this,
1263  dK_past_dx,
1264  dK_this_dx,
1265  (ip > 0) ? ppath.lstep[ip - 1] : Numeric(1.0),
1266  ip == 0);
1267 
1268  get_stepwise_effective_source(J(ip, joker, joker),
1269  dJ_dx(ip, joker, joker, joker),
1270  K_this,
1271  a,
1272  S,
1273  dK_this_dx,
1274  da_dx,
1275  dS_dx,
1276  B,
1277  dB_dT,
1278  jacobian_quantities,
1279  j_analytical_do);
1280 
1281  swap(K_past, K_this);
1282  swap(dK_past_dx, dK_this_dx);
1283  }
1284  }
1285 
1286  // iy_transmission
1287  Tensor3 iy_trans_new;
1288  if (iy_agenda_call1) {
1289  iy_trans_new = ppvar_trans_cumulat(np - 1, joker, joker, joker);
1290  } else {
1291  iy_transmission_mult(iy_trans_new,
1292  iy_transmission,
1293  ppvar_trans_cumulat(np - 1, joker, joker, joker));
1294  }
1295 
1296  // Copy transmission to iy_aux
1297  for (Index i = 0; i < naux; i++) {
1298  if (iy_aux_vars[i] == "Optical depth") {
1299  for (Index iv = 0; iv < nf; iv++) {
1300  iy_aux[i](iv, joker) = -log(ppvar_trans_cumulat(np - 1, iv, 0, 0));
1301  }
1302  }
1303  }
1304 
1305  // Radiative background
1306  get_iy_of_background(ws,
1307  iy,
1308  diy_dx,
1309  iy_trans_new,
1310  iy_id,
1311  jacobian_do,
1312  jacobian_quantities,
1313  ppath,
1314  rte_pos2,
1315  atmosphere_dim,
1316  nlte_field,
1317  cloudbox_on,
1318  stokes_dim,
1319  f_grid,
1320  iy_unit,
1321  surface_props_data,
1322  iy_main_agenda,
1323  iy_space_agenda,
1324  iy_surface_agenda,
1325  iy_cloudbox_agenda,
1326  iy_agenda_call1,
1327  verbosity);
1328  //
1329  ppvar_iy(joker, joker, np - 1) = iy;
1330 
1331  // Radiative transfer calculations
1332  if (np > 1) {
1333  for (Index iv = 0; iv < nf; iv++) {
1334  Vector through_level(ns), from_level(ns);
1335  Vector dfrom_level_dx(ns);
1336  Matrix one_minus_transmission(ns, ns);
1337 
1338  for (Index ip = np - 2; ip >= 0; ip--) {
1339  ConstMatrixView T = trans_partial(ip + 1, iv, joker, joker);
1340 
1341  if (j_analytical_do) {
1342  if (stokes_dim > 1) {
1343  id_mat(one_minus_transmission);
1344  } else {
1345  one_minus_transmission = 1.;
1346  }
1347  one_minus_transmission -= T;
1348  }
1349 
1350  from_level = J(ip, iv, joker);
1351  from_level += J(ip + 1, iv, joker);
1352  from_level *= 0.5;
1353  through_level = iy(iv, joker);
1354  through_level -= from_level;
1355 
1356  if (j_analytical_do) {
1357  FOR_ANALYTICAL_JACOBIANS_DO(
1358  get_diydx(diy_dpath[iq](ip, iv, joker),
1359  diy_dpath[iq](ip + 1, iv, joker),
1360  one_minus_transmission,
1361  ppvar_trans_cumulat(ip, iv, joker, joker),
1362  dtrans_partial_dx_above(ip + 1, iq, iv, joker, joker),
1363  dtrans_partial_dx_below(ip + 1, iq, iv, joker, joker),
1364  through_level,
1365  dJ_dx(ip, iq, iv, joker),
1366  dJ_dx(ip + 1, iq, iv, joker),
1367  stokes_dim);)
1368  }
1369 
1370  // Equation is I1 = T (I0 - 0.5(J_1+J_2)) + 0.5(J_1+J_2)
1371  mult(iy(iv, joker), T, through_level);
1372  iy(iv, joker) += from_level;
1373 
1374  ppvar_iy(iv, joker, ip) = iy(iv, joker);
1375  }
1376  }
1377  }
1378 
1379  // Finalize analytical Jacobians
1380  if (j_analytical_do) {
1381  rtmethods_jacobian_finalisation(ws,
1382  diy_dx,
1383  diy_dpath,
1384  ns,
1385  nf,
1386  np,
1387  atmosphere_dim,
1388  ppath,
1389  ppvar_p,
1390  ppvar_t,
1391  ppvar_vmr,
1392  iy_agenda_call1,
1393  iy_transmission,
1394  water_p_eq_agenda,
1395  jacobian_quantities,
1396  jac_species_i,
1397  jac_is_t);
1398  }
1399 
1400  // Unit conversions
1401  if (iy_agenda_call1) {
1402  rtmethods_unit_conversion(iy,
1403  diy_dx,
1404  ppvar_iy,
1405  ns,
1406  np,
1407  f_grid,
1408  ppath,
1409  jacobian_quantities,
1410  j_analytical_do,
1411  iy_unit);
1412  }
1413 }
1414 
1415 /* Workspace method: Doxygen documentation will be auto-generated */
1416 void iyHybrid2(Workspace& ws,
1417  Matrix& iy,
1418  ArrayOfMatrix& iy_aux,
1419  ArrayOfTensor3& diy_dx,
1420  Vector& ppvar_p,
1421  Vector& ppvar_t,
1422  EnergyLevelMap& ppvar_nlte,
1423  Matrix& ppvar_vmr,
1424  Matrix& ppvar_wind,
1425  Matrix& ppvar_mag,
1426  Matrix& ppvar_pnd,
1427  Matrix& ppvar_f,
1428  Tensor3& ppvar_iy,
1429  Tensor4& ppvar_trans_cumulat,
1430  const Index& iy_id,
1431  const Index& stokes_dim,
1432  const Vector& f_grid,
1433  const Index& atmosphere_dim,
1434  const Vector& p_grid,
1435  const Tensor3& t_field,
1436  const EnergyLevelMap& nlte_field,
1437  const Tensor4& vmr_field,
1438  const ArrayOfArrayOfSpeciesTag& abs_species,
1439  const Tensor3& wind_u_field,
1440  const Tensor3& wind_v_field,
1441  const Tensor3& wind_w_field,
1442  const Tensor3& mag_u_field,
1443  const Tensor3& mag_v_field,
1444  const Tensor3& mag_w_field,
1445  const Index& cloudbox_on,
1446  const ArrayOfIndex& cloudbox_limits,
1447  const Tensor4& pnd_field,
1448  const ArrayOfTensor4& dpnd_field_dx,
1449  const ArrayOfString& scat_species,
1450  const ArrayOfArrayOfSingleScatteringData& scat_data,
1451  const String& iy_unit,
1452  const ArrayOfString& iy_aux_vars,
1453  const Index& jacobian_do,
1454  const ArrayOfRetrievalQuantity& jacobian_quantities,
1455  const Agenda& propmat_clearsky_agenda,
1456  const Agenda& water_p_eq_agenda,
1457  const Agenda& iy_main_agenda,
1458  const Agenda& iy_space_agenda,
1459  const Agenda& iy_surface_agenda,
1460  const Agenda& iy_cloudbox_agenda,
1461  const Index& iy_agenda_call1,
1462  const Tensor3& iy_transmission,
1463  const Ppath& ppath,
1464  const Vector& rte_pos2,
1465  const Numeric& rte_alonglos_v,
1466  const Tensor3& surface_props_data,
1467  const Tensor7& cloudbox_field,
1468  const Vector& za_grid,
1469  const Index& Naa,
1470  const Index& t_interp_order,
1471  const Verbosity& verbosity) {
1472  // If cloudbox off, switch to use clearsky method
1473  if (!cloudbox_on) {
1474  Tensor4 dummy;
1475  iyEmissionStandard(ws,
1476  iy,
1477  iy_aux,
1478  diy_dx,
1479  ppvar_p,
1480  ppvar_t,
1481  ppvar_nlte,
1482  ppvar_vmr,
1483  ppvar_wind,
1484  ppvar_mag,
1485  ppvar_f,
1486  ppvar_iy,
1487  ppvar_trans_cumulat,
1488  dummy,
1489  iy_id,
1490  stokes_dim,
1491  f_grid,
1492  atmosphere_dim,
1493  p_grid,
1494  t_field,
1495  nlte_field,
1496  vmr_field,
1497  abs_species,
1498  wind_u_field,
1499  wind_v_field,
1500  wind_w_field,
1501  mag_u_field,
1502  mag_v_field,
1503  mag_w_field,
1504  cloudbox_on,
1505  iy_unit,
1506  iy_aux_vars,
1507  jacobian_do,
1508  jacobian_quantities,
1509  ppath,
1510  rte_pos2,
1511  propmat_clearsky_agenda,
1512  water_p_eq_agenda,
1513  iy_main_agenda,
1514  iy_space_agenda,
1515  iy_surface_agenda,
1516  iy_cloudbox_agenda,
1517  iy_agenda_call1,
1518  iy_transmission,
1519  rte_alonglos_v,
1520  surface_props_data,
1521  verbosity);
1522  return;
1523  }
1524 
1525  // Some basic sizes
1526  const Index nf = f_grid.nelem();
1527  const Index ns = stokes_dim;
1528  const Index np = ppath.np;
1529 
1530  // Radiative background index
1531  const Index rbi = ppath_what_background(ppath);
1532 
1533  // Throw error if unsupported features are requested
1534  if (atmosphere_dim != 1)
1535  throw runtime_error(
1536  "With cloudbox on, this method handles only 1D calculations.");
1537  if (Naa < 3) throw runtime_error("Naa must be > 2.");
1538  if (jacobian_do)
1539  if (dpnd_field_dx.nelem() != jacobian_quantities.nelem())
1540  throw runtime_error(
1541  "*dpnd_field_dx* not properly initialized:\n"
1542  "Number of elements in dpnd_field_dx must be equal number of jacobian"
1543  " quantities.\n(Note: jacobians have to be defined BEFORE *pnd_field*"
1544  " is calculated/set.");
1545  if (rbi < 1 || rbi > 9)
1546  throw runtime_error(
1547  "ppath.background is invalid. Check your "
1548  "calculation of *ppath*?");
1549  if (rbi == 3 || rbi == 4)
1550  throw runtime_error(
1551  "The propagation path ends inside or at boundary of "
1552  "the cloudbox.\nFor this method, *ppath* must be "
1553  "calculated in this way:\n ppathCalc( cloudbox_on = 0 ).");
1554  // iy_aux_vars checked below
1555  // Checks of i_field
1556  if (cloudbox_field.ncols() != stokes_dim)
1557  throw runtime_error(
1558  "Obtained *cloudbox_field* number of Stokes elements inconsistent with "
1559  "*stokes_dim*.");
1560  if (cloudbox_field.nrows() != 1)
1561  throw runtime_error(
1562  "Obtained *cloudbox_field* has wrong number of azimuth angles.");
1563  if (cloudbox_field.npages() != za_grid.nelem())
1564  throw runtime_error(
1565  "Obtained *cloudbox_field* number of zenith angles inconsistent with "
1566  "*za_grid*.");
1567  if (cloudbox_field.nbooks() != 1)
1568  throw runtime_error(
1569  "Obtained *cloudbox_field* has wrong number of longitude points.");
1570  if (cloudbox_field.nshelves() != 1)
1571  throw runtime_error(
1572  "Obtained *cloudbox_field* has wrong number of latitude points.");
1573  if (cloudbox_field.nvitrines() != cloudbox_limits[1] - cloudbox_limits[0] + 1)
1574  throw runtime_error(
1575  "Obtained *cloudbox_field* number of pressure points inconsistent with "
1576  "*cloudbox_limits*.");
1577  if (cloudbox_field.nlibraries() != nf)
1578  throw runtime_error(
1579  "Obtained *cloudbox_field* number of frequency points inconsistent with "
1580  "*f_grid*.");
1581 
1582  // Init Jacobian quantities
1583  Index j_analytical_do = 0;
1584  if (jacobian_do) FOR_ANALYTICAL_JACOBIANS_DO(j_analytical_do = 1;)
1585  //
1586  const Index nq = j_analytical_do ? jacobian_quantities.nelem() : 0;
1587  ArrayOfTensor3 diy_dpath(nq);
1588  ArrayOfIndex jac_species_i(nq), jac_scat_i(nq), jac_is_t(nq), jac_wind_i(nq);
1589  ArrayOfIndex jac_mag_i(nq), jac_other(nq);
1590 
1591  if (j_analytical_do)
1592  rtmethods_jacobian_init(jac_species_i,
1593  jac_scat_i,
1594  jac_is_t,
1595  jac_wind_i,
1596  jac_mag_i,
1597  jac_other,
1598  diy_dx,
1599  diy_dpath,
1600  ns,
1601  nf,
1602  np,
1603  nq,
1604  abs_species,
1605  cloudbox_on,
1606  scat_species,
1607  dpnd_field_dx,
1608  jacobian_quantities,
1609  iy_agenda_call1);
1610 
1611  // Init iy_aux and fill where possible
1612  const Index naux = iy_aux_vars.nelem();
1613  iy_aux.resize(naux);
1614  //
1615  for (Index i = 0; i < naux; i++) {
1616  iy_aux[i].resize(nf, ns);
1617 
1618  if (iy_aux_vars[i] == "Optical depth") { /*pass*/
1619  } // Filled below
1620  else if (iy_aux_vars[i] == "Radiative background")
1621  iy_aux[i] = (Numeric)min((Index)2, rbi - 1);
1622  else {
1623  ostringstream os;
1624  os << "The only allowed strings in *iy_aux_vars* are:\n"
1625  << " \"Radiative background\"\n"
1626  << " \"Optical depth\"\n"
1627  << "but you have selected: \"" << iy_aux_vars[i] << "\"";
1628  throw runtime_error(os.str());
1629  }
1630  }
1631 
1632  // Get atmospheric and radiative variables along the propagation path
1633  ppvar_trans_cumulat.resize(np, nf, ns, ns);
1634  ppvar_iy.resize(nf, ns, np);
1635 
1636  ArrayOfTransmissionMatrix lyr_tra(np, TransmissionMatrix(nf, ns));
1637  ArrayOfRadiationVector lvl_rad(np, RadiationVector(nf, ns));
1638  ArrayOfArrayOfRadiationVector dlvl_rad(
1639  np, ArrayOfRadiationVector(nq, RadiationVector(nf, ns)));
1640  ArrayOfRadiationVector src_rad(np, RadiationVector(nf, ns));
1641  ArrayOfArrayOfRadiationVector dsrc_rad(
1642  np, ArrayOfRadiationVector(nq, RadiationVector(nf, ns)));
1643 
1644  ArrayOfArrayOfTransmissionMatrix dlyr_tra_above(
1645  np, ArrayOfTransmissionMatrix(nq, TransmissionMatrix(nf, ns)));
1646  ArrayOfArrayOfTransmissionMatrix dlyr_tra_below(
1647  np, ArrayOfTransmissionMatrix(nq, TransmissionMatrix(nf, ns)));
1648 
1649  ArrayOfIndex clear2cloudy;
1650  //
1651  if (np == 1 && rbi == 1) { // i.e. ppath is totally outside the atmosphere:
1652  ppvar_p.resize(0);
1653  ppvar_t.resize(0);
1654  ppvar_vmr.resize(0, 0);
1655  ppvar_wind.resize(0, 0);
1656  ppvar_mag.resize(0, 0);
1657  ppvar_f.resize(0, 0);
1658  ppvar_trans_cumulat = 1;
1659  } else {
1660  // Basic atmospheric variables
1661  get_ppath_atmvars(ppvar_p,
1662  ppvar_t,
1663  ppvar_nlte,
1664  ppvar_vmr,
1665  ppvar_wind,
1666  ppvar_mag,
1667  ppath,
1668  atmosphere_dim,
1669  p_grid,
1670  t_field,
1671  nlte_field,
1672  vmr_field,
1673  wind_u_field,
1674  wind_v_field,
1675  wind_w_field,
1676  mag_u_field,
1677  mag_v_field,
1678  mag_w_field);
1679 
1680  get_ppath_f(
1681  ppvar_f, ppath, f_grid, atmosphere_dim, rte_alonglos_v, ppvar_wind);
1682 
1683  // here, the cloudbox is on, ie we don't need to check and branch this here
1684  // anymore.
1685  ArrayOfMatrix ppvar_dpnd_dx;
1686  //
1687  get_ppath_cloudvars(clear2cloudy,
1688  ppvar_pnd,
1689  ppvar_dpnd_dx,
1690  ppath,
1691  atmosphere_dim,
1692  cloudbox_limits,
1693  pnd_field,
1694  dpnd_field_dx);
1695 
1696  // Size radiative variables always used
1697  Vector B(nf);
1698  PropagationMatrix K_this(nf, ns), K_past(nf, ns), Kp(nf, ns);
1699  StokesVector a(nf, ns), S(nf, ns), Sp(nf, ns);
1700  ArrayOfIndex lte(np);
1701 
1702  // Init variables only used if analytical jacobians done
1703  Vector dB_dT(0);
1704  ArrayOfPropagationMatrix dK_this_dx(nq), dK_past_dx(nq), dKp_dx(nq);
1705  ArrayOfStokesVector da_dx(nq), dS_dx(nq), dSp_dx(nq);
1706 
1707  // HSE variables
1708  Index temperature_derivative_position = -1;
1709  bool do_hse = false;
1710 
1711  if (j_analytical_do) {
1712  dB_dT.resize(nf);
1713  FOR_ANALYTICAL_JACOBIANS_DO(
1714  dK_this_dx[iq] = PropagationMatrix(nf, ns);
1715  dK_past_dx[iq] = PropagationMatrix(nf, ns);
1716  dKp_dx[iq] = PropagationMatrix(nf, ns);
1717  da_dx[iq] = StokesVector(nf, ns);
1718  dS_dx[iq] = StokesVector(nf, ns);
1719  dSp_dx[iq] = StokesVector(nf, ns);
1720  if (jacobian_quantities[iq] == JacPropMatType::Temperature) {
1721  temperature_derivative_position = iq;
1722  do_hse = jacobian_quantities[iq].Subtag() == "HSE on";
1723  })
1724  }
1725  const bool temperature_jacobian =
1726  j_analytical_do and do_temperature_jacobian(jacobian_quantities);
1727 
1728  // Loop ppath points and determine radiative properties
1729  for (Index ip = 0; ip < np; ip++) {
1730  get_stepwise_blackbody_radiation(
1731  B, dB_dT, ppvar_f(joker, ip), ppvar_t[ip], temperature_jacobian);
1732 
1733  get_stepwise_clearsky_propmat(ws,
1734  K_this,
1735  S,
1736  lte[ip],
1737  dK_this_dx,
1738  dS_dx,
1739  propmat_clearsky_agenda,
1740  jacobian_quantities,
1741  ppvar_f(joker, ip),
1742  ppvar_mag(joker, ip),
1743  ppath.los(ip, joker),
1744  ppvar_nlte[ip],
1745  ppvar_vmr(joker, ip),
1746  ppvar_t[ip],
1747  ppvar_p[ip],
1748  jac_species_i,
1749  j_analytical_do);
1750 
1751  if (j_analytical_do)
1752  adapt_stepwise_partial_derivatives(dK_this_dx,
1753  dS_dx,
1754  jacobian_quantities,
1755  ppvar_f(joker, ip),
1756  ppath.los(ip, joker),
1757  ppvar_vmr(joker, ip),
1758  ppvar_t[ip],
1759  ppvar_p[ip],
1760  jac_species_i,
1761  jac_wind_i,
1762  lte[ip],
1763  atmosphere_dim,
1764  j_analytical_do);
1765 
1766  if (clear2cloudy[ip] + 1) {
1767  get_stepwise_scattersky_propmat(a,
1768  Kp,
1769  da_dx,
1770  dKp_dx,
1771  jacobian_quantities,
1772  ppvar_pnd(joker, Range(ip, 1)),
1773  ppvar_dpnd_dx,
1774  ip,
1775  scat_data,
1776  ppath.los(ip, joker),
1777  ppvar_t[Range(ip, 1)],
1778  atmosphere_dim,
1779  jacobian_do);
1780  a += K_this;
1781  K_this += Kp;
1782 
1783  if (j_analytical_do)
1784  FOR_ANALYTICAL_JACOBIANS_DO(da_dx[iq] += dK_this_dx[iq];
1785  dK_this_dx[iq] += dKp_dx[iq];)
1786 
1787  Vector aa_grid;
1788  nlinspace(aa_grid, 0, 360, Naa);
1789  //
1790  get_stepwise_scattersky_source(Sp,
1791  dSp_dx,
1792  jacobian_quantities,
1793  ppvar_pnd(joker, ip),
1794  ppvar_dpnd_dx,
1795  ip,
1796  scat_data,
1797  cloudbox_field,
1798  za_grid,
1799  aa_grid,
1800  ppath.los(Range(ip, 1), joker),
1801  ppath.gp_p[ip],
1802  ppvar_t[Range(ip, 1)],
1803  atmosphere_dim,
1804  jacobian_do,
1805  t_interp_order);
1806  S += Sp;
1807 
1808  if (j_analytical_do)
1809  FOR_ANALYTICAL_JACOBIANS_DO(dS_dx[iq] += dSp_dx[iq];)
1810  } else { // no particles present at this level
1811  a = K_this;
1812  if (j_analytical_do)
1813  FOR_ANALYTICAL_JACOBIANS_DO(da_dx[iq] = dK_this_dx[iq];)
1814  }
1815 
1816  if (ip not_eq 0) {
1817  const Numeric dr_dT_past =
1818  do_hse ? ppath.lstep[ip - 1] / (2.0 * ppvar_t[ip - 1]) : 0;
1819  const Numeric dr_dT_this =
1820  do_hse ? ppath.lstep[ip - 1] / (2.0 * ppvar_t[ip]) : 0;
1821  stepwise_transmission(lyr_tra[ip],
1822  dlyr_tra_above[ip],
1823  dlyr_tra_below[ip],
1824  K_past,
1825  K_this,
1826  dK_past_dx,
1827  dK_this_dx,
1828  ppath.lstep[ip - 1],
1829  dr_dT_past,
1830  dr_dT_this,
1831  temperature_derivative_position);
1832  }
1833 
1834  stepwise_source(src_rad[ip],
1835  dsrc_rad[ip],
1836  K_this,
1837  a,
1838  S,
1839  dK_this_dx,
1840  da_dx,
1841  dS_dx,
1842  B,
1843  dB_dT,
1844  jacobian_quantities,
1845  jacobian_do);
1846 
1847  swap(K_past, K_this);
1848  swap(dK_past_dx, dK_this_dx);
1849  }
1850  }
1851 
1852  const ArrayOfTransmissionMatrix tot_tra =
1853  cumulative_transmission(lyr_tra, CumulativeTransmission::Forward);
1854 
1855  // iy_transmission
1856  Tensor3 iy_trans_new;
1857  if (iy_agenda_call1)
1858  iy_trans_new = tot_tra[np - 1];
1859  else
1860  iy_transmission_mult(iy_trans_new, iy_transmission, tot_tra[np - 1]);
1861 
1862  // Copy transmission to iy_aux
1863  for (Index i = 0; i < naux; i++)
1864  if (iy_aux_vars[i] == "Optical depth")
1865  for (Index iv = 0; iv < nf; iv++)
1866  iy_aux[i](iv, joker) = -log(ppvar_trans_cumulat(np - 1, iv, 0, 0));
1867 
1868  // Radiative background
1869  get_iy_of_background(ws,
1870  iy,
1871  diy_dx,
1872  iy_trans_new,
1873  iy_id,
1874  jacobian_do,
1875  jacobian_quantities,
1876  ppath,
1877  rte_pos2,
1878  atmosphere_dim,
1879  nlte_field,
1880  cloudbox_on,
1881  stokes_dim,
1882  f_grid,
1883  iy_unit,
1884  surface_props_data,
1885  iy_main_agenda,
1886  iy_space_agenda,
1887  iy_surface_agenda,
1888  iy_cloudbox_agenda,
1889  iy_agenda_call1,
1890  verbosity);
1891 
1892  lvl_rad[np - 1] = iy;
1893 
1894  // Radiative transfer calculations
1895  for (Index ip = np - 2; ip >= 0; ip--) {
1896  lvl_rad[ip] = lvl_rad[ip + 1];
1897  update_radiation_vector(lvl_rad[ip],
1898  dlvl_rad[ip],
1899  dlvl_rad[ip + 1],
1900  src_rad[ip],
1901  src_rad[ip + 1],
1902  dsrc_rad[ip],
1903  dsrc_rad[ip + 1],
1904  lyr_tra[ip + 1],
1905  tot_tra[ip],
1906  dlyr_tra_above[ip + 1],
1907  dlyr_tra_below[ip + 1],
1908  RadiativeTransferSolver::Emission);
1909  }
1910 
1911  // Copy back to ARTS external style
1912  iy = lvl_rad[0];
1913  for (Index ip = 0; ip < lvl_rad.nelem(); ip++) {
1914  ppvar_trans_cumulat(ip, joker, joker, joker) = tot_tra[ip];
1915  ppvar_iy(joker, joker, ip) = lvl_rad[ip];
1916  if (j_analytical_do)
1917  FOR_ANALYTICAL_JACOBIANS_DO(diy_dpath[iq](ip, joker, joker) =
1918  dlvl_rad[ip][iq];);
1919  }
1920 
1921  // Finalize analytical Jacobians
1922  if (j_analytical_do)
1923  rtmethods_jacobian_finalisation(ws,
1924  diy_dx,
1925  diy_dpath,
1926  ns,
1927  nf,
1928  np,
1929  atmosphere_dim,
1930  ppath,
1931  ppvar_p,
1932  ppvar_t,
1933  ppvar_vmr,
1934  iy_agenda_call1,
1935  iy_transmission,
1936  water_p_eq_agenda,
1937  jacobian_quantities,
1938  jac_species_i,
1939  jac_is_t);
1940 
1941  // Unit conversions
1942  if (iy_agenda_call1)
1943  rtmethods_unit_conversion(iy,
1944  diy_dx,
1945  ppvar_iy,
1946  ns,
1947  np,
1948  f_grid,
1949  ppath,
1950  jacobian_quantities,
1951  j_analytical_do,
1952  iy_unit);
1953 }
Index
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
get_iy_of_background
void get_iy_of_background(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, ConstTensor3View iy_transmission, const Index &iy_id, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Ppath &ppath, ConstVectorView rte_pos2, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, ConstVectorView f_grid, const String &iy_unit, ConstTensor3View surface_props_data, const Agenda &iy_main_agenda, const Agenda &iy_space_agenda, const Agenda &iy_surface_agenda, const Agenda &iy_cloudbox_agenda, const Index &iy_agenda_call1, const Verbosity &verbosity)
Determines iy of the "background" of a propgation path.
Definition: rte.cc:916
TransmissionMatrix
Class to keep track of Transmission Matrices for Stokes Dim 1-4.
Definition: transmissionmatrix.h:38
ns
#define ns
Definition: legacy_continua.cc:22183
get_ppath_f
void get_ppath_f(Matrix &ppath_f, const Ppath &ppath, ConstVectorView f_grid, const Index &atmosphere_dim, const Numeric &rte_alonglos_v, ConstMatrixView ppath_wind)
Determines the Doppler shifted frequencies along the propagation path.
Definition: rte.cc:1257
Agenda
The Agenda class.
Definition: agenda_class.h:44
get_stepwise_clearsky_propmat
void get_stepwise_clearsky_propmat(Workspace &ws, PropagationMatrix &K, StokesVector &S, Index &lte, ArrayOfPropagationMatrix &dK_dx, ArrayOfStokesVector &dS_dx, const Agenda &propmat_clearsky_agenda, const ArrayOfRetrievalQuantity &jacobian_quantities, ConstVectorView ppath_f_grid, ConstVectorView ppath_magnetic_field, ConstVectorView ppath_line_of_sight, const EnergyLevelMap &ppath_nlte, ConstVectorView ppath_vmrs, const Numeric &ppath_temperature, const Numeric &ppath_pressure, const ArrayOfIndex &jacobian_species, const bool &jacobian_do)
Gets the clearsky propgation matrix and NLTE contributions.
Definition: rte.cc:1322
iy_transmission_mult
void iy_transmission_mult(Tensor3 &iy_trans_total, ConstTensor3View iy_trans_old, ConstTensor3View iy_trans_new)
Multiplicates iy_transmission with transmissions.
Definition: rte.cc:2253
Array::nelem
Index nelem() const
Number of elements.
Definition: array.h:195
iyHybrid
void iyHybrid(Workspace &ws, Matrix &iy, ArrayOfMatrix &iy_aux, ArrayOfTensor3 &diy_dx, Vector &ppvar_p, Vector &ppvar_t, EnergyLevelMap &ppvar_nlte, Matrix &ppvar_vmr, Matrix &ppvar_wind, Matrix &ppvar_mag, Matrix &ppvar_pnd, Matrix &ppvar_f, Tensor3 &ppvar_iy, Tensor4 &ppvar_trans_cumulat, const Index &iy_id, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Tensor3 &t_field, const EnergyLevelMap &nlte_field, const Tensor4 &vmr_field, const ArrayOfArrayOfSpeciesTag &abs_species, const Tensor3 &wind_u_field, const Tensor3 &wind_v_field, const Tensor3 &wind_w_field, const Tensor3 &mag_u_field, const Tensor3 &mag_v_field, const Tensor3 &mag_w_field, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Tensor4 &pnd_field, const ArrayOfTensor4 &dpnd_field_dx, const ArrayOfString &scat_species, const ArrayOfArrayOfSingleScatteringData &scat_data, const String &iy_unit, const ArrayOfString &iy_aux_vars, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &iy_main_agenda, const Agenda &iy_space_agenda, const Agenda &iy_surface_agenda, const Agenda &iy_cloudbox_agenda, const Index &iy_agenda_call1, const Tensor3 &iy_transmission, const Ppath &ppath, const Vector &rte_pos2, const Numeric &rte_alonglos_v, const Tensor3 &surface_props_data, const Tensor7 &cloudbox_field, const Vector &za_grid, const Index &Naa, const Index &t_interp_order, const Verbosity &verbosity)
WORKSPACE METHOD: iyHybrid.
Definition: m_fos.cc:855
Ppath::los
Matrix los
Line-of-sight at each ppath point.
Definition: ppath.h:66
get_stepwise_transmission_matrix
void get_stepwise_transmission_matrix(Tensor3View cumulative_transmission, Tensor3View T, Tensor4View dT_close_dx, Tensor4View dT_far_dx, ConstTensor3View cumulative_transmission_close, const PropagationMatrix &K_close, const PropagationMatrix &K_far, const ArrayOfPropagationMatrix &dK_close_dx, const ArrayOfPropagationMatrix &dK_far_dx, const Numeric &ppath_distance, const bool &first_level, const Numeric &dppath_distance_dT_HSE_close, const Numeric &dppath_distance_dT_HSE_far, const Index &temperature_derivative_position_if_hse_is_active)
Computes layer transmission matrix and cumulative transmission.
Definition: rte.cc:1879
Vector
The Vector class.
Definition: matpackI.h:860
stepwise_transmission
void stepwise_transmission(TransmissionMatrix &T, ArrayOfTransmissionMatrix &dT1, ArrayOfTransmissionMatrix &dT2, const PropagationMatrix &K1, const PropagationMatrix &K2, const ArrayOfPropagationMatrix &dK1, const ArrayOfPropagationMatrix &dK2, const Numeric &r, const Numeric &dr_dtemp1, const Numeric &dr_dtemp2, const Index temp_deriv_pos)
Set the stepwise transmission matrix.
Definition: transmissionmatrix.cc:1353
Tensor4
The Tensor4 class.
Definition: matpackIV.h:421
Range
The range class.
Definition: matpackI.h:160
Ppath::lstep
Vector lstep
The length between ppath points.
Definition: ppath.h:70
ConstTensor7View::nrows
Index nrows() const
Returns the number of rows.
Definition: matpackVII.cc:57
iyHybrid2
void iyHybrid2(Workspace &ws, Matrix &iy, ArrayOfMatrix &iy_aux, ArrayOfTensor3 &diy_dx, Vector &ppvar_p, Vector &ppvar_t, EnergyLevelMap &ppvar_nlte, Matrix &ppvar_vmr, Matrix &ppvar_wind, Matrix &ppvar_mag, Matrix &ppvar_pnd, Matrix &ppvar_f, Tensor3 &ppvar_iy, Tensor4 &ppvar_trans_cumulat, const Index &iy_id, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Tensor3 &t_field, const EnergyLevelMap &nlte_field, const Tensor4 &vmr_field, const ArrayOfArrayOfSpeciesTag &abs_species, const Tensor3 &wind_u_field, const Tensor3 &wind_v_field, const Tensor3 &wind_w_field, const Tensor3 &mag_u_field, const Tensor3 &mag_v_field, const Tensor3 &mag_w_field, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Tensor4 &pnd_field, const ArrayOfTensor4 &dpnd_field_dx, const ArrayOfString &scat_species, const ArrayOfArrayOfSingleScatteringData &scat_data, const String &iy_unit, const ArrayOfString &iy_aux_vars, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &iy_main_agenda, const Agenda &iy_space_agenda, const Agenda &iy_surface_agenda, const Agenda &iy_cloudbox_agenda, const Index &iy_agenda_call1, const Tensor3 &iy_transmission, const Ppath &ppath, const Vector &rte_pos2, const Numeric &rte_alonglos_v, const Tensor3 &surface_props_data, const Tensor7 &cloudbox_field, const Vector &za_grid, const Index &Naa, const Index &t_interp_order, const Verbosity &verbosity)
WORKSPACE METHOD: iyHybrid2.
Definition: m_fos.cc:1416
get_ppath_cloudvars
void get_ppath_cloudvars(ArrayOfIndex &clear2cloudy, Matrix &ppath_pnd, ArrayOfMatrix &ppath_dpnd_dx, const Ppath &ppath, const Index &atmosphere_dim, const ArrayOfIndex &cloudbox_limits, const Tensor4 &pnd_field, const ArrayOfTensor4 &dpnd_field_dx)
Determines the particle fields along a propagation path.
Definition: rte.cc:1153
rtmethods_jacobian_init
void rtmethods_jacobian_init(ArrayOfIndex &jac_species_i, ArrayOfIndex &jac_scat_i, ArrayOfIndex &jac_is_t, ArrayOfIndex &jac_wind_i, ArrayOfIndex &jac_mag_i, ArrayOfIndex &jac_other, ArrayOfTensor3 &diy_dx, ArrayOfTensor3 &diy_dpath, const Index &ns, const Index &nf, const Index &np, const Index &nq, const ArrayOfArrayOfSpeciesTag &abs_species, const Index &cloudbox_on, const ArrayOfString &scat_species, const ArrayOfTensor4 &dpnd_field_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, const Index &iy_agenda_call1, const bool is_active)
This function fixes the initial steps around Jacobian calculations, to be done inside radiative trans...
Definition: rte.cc:2348
Tensor7
The Tensor7 class.
Definition: matpackVII.h:2382
basic_ostringstream
Definition: sstream.h:171
get_ppath_atmvars
void get_ppath_atmvars(Vector &ppath_p, Vector &ppath_t, EnergyLevelMap &ppath_nlte, Matrix &ppath_vmr, Matrix &ppath_wind, Matrix &ppath_mag, const Ppath &ppath, const Index &atmosphere_dim, ConstVectorView p_grid, ConstTensor3View t_field, const EnergyLevelMap &nlte_field, ConstTensor4View vmr_field, ConstTensor3View wind_u_field, ConstTensor3View wind_v_field, ConstTensor3View wind_w_field, ConstTensor3View mag_u_field, ConstTensor3View mag_v_field, ConstTensor3View mag_w_field)
Determines pressure, temperature, VMR, winds and magnetic field for each propgataion path point...
Definition: rte.cc:1034
min
#define min(a, b)
Definition: legacy_continua.cc:20628
DEG2RAD
const Numeric DEG2RAD
StokesVector
Stokes vector is as Propagation matrix but only has 4 possible values.
Definition: propagationmatrix.h:1075
FOR_ANALYTICAL_JACOBIANS_DO
#define FOR_ANALYTICAL_JACOBIANS_DO(what_to_do)
Definition: jacobian.h:405
ConstVectorView::nelem
Index nelem() const
Returns the number of elements.
Definition: matpackI.cc:51
ArrayOfRadiationVector
Array< RadiationVector > ArrayOfRadiationVector
Definition: transmissionmatrix.h:899
ppath_what_background
Index ppath_what_background(const Ppath &ppath)
Returns the case number for the radiative background.
Definition: ppath.cc:1494
my_basic_string< char >
Tensor3
The Tensor3 class.
Definition: matpackIII.h:339
swap
void swap(ComplexVector &v1, ComplexVector &v2)
Swaps two objects.
Definition: complex.cc:731
arts.h
The global header file for ARTS.
_CS_basic_sstream_base::str
_CS_string_type str() const
Definition: sstream.h:491
cumulative_transmission
ArrayOfTransmissionMatrix cumulative_transmission(const ArrayOfTransmissionMatrix &T, const CumulativeTransmission type)
Accumulate the transmission matrix over all layers.
Definition: transmissionmatrix.cc:1527
ArrayOfTransmissionMatrix
Array< TransmissionMatrix > ArrayOfTransmissionMatrix
Definition: transmissionmatrix.h:896
stepwise_source
void stepwise_source(RadiationVector &J, ArrayOfRadiationVector &dJ, const PropagationMatrix &K, const StokesVector &a, const StokesVector &S, const ArrayOfPropagationMatrix &dK, const ArrayOfStokesVector &da, const ArrayOfStokesVector &dS, const ConstVectorView B, const ConstVectorView dB_dT, const ArrayOfRetrievalQuantity &jacobian_quantities, const bool &jacobian_do)
Set the stepwise source.
Definition: transmissionmatrix.cc:1371
math_funcs.h
nlinspace
void nlinspace(Vector &x, const Numeric start, const Numeric stop, const Index n)
nlinspace
Definition: math_funcs.cc:231
get_stepwise_effective_source
void get_stepwise_effective_source(MatrixView J, Tensor3View dJ_dx, const PropagationMatrix &K, const StokesVector &a, const StokesVector &S, const ArrayOfPropagationMatrix &dK_dx, const ArrayOfStokesVector &da_dx, const ArrayOfStokesVector &dS_dx, ConstVectorView B, ConstVectorView dB_dT, const ArrayOfRetrievalQuantity &jacobian_quantities, const bool &jacobian_do)
Gets the effective source at propagation path point.
Definition: rte.cc:1427
rtmethods_unit_conversion
void rtmethods_unit_conversion(Matrix &iy, ArrayOfTensor3 &diy_dx, Tensor3 &ppvar_iy, const Index &ns, const Index &np, const Vector &f_grid, const Ppath &ppath, const ArrayOfRetrievalQuantity &jacobian_quantities, const Index &j_analytical_do, const String &iy_unit)
This function handles the unit conversion to be done at the end of some radiative transfer WSMs...
Definition: rte.cc:2553
iyEmissionStandard
void iyEmissionStandard(Workspace &ws, Matrix &iy, ArrayOfMatrix &iy_aux, ArrayOfTensor3 &diy_dx, Vector &ppvar_p, Vector &ppvar_t, EnergyLevelMap &ppvar_nlte, Matrix &ppvar_vmr, Matrix &ppvar_wind, Matrix &ppvar_mag, Matrix &ppvar_f, Tensor3 &ppvar_iy, Tensor4 &ppvar_trans_cumulat, Tensor4 &ppvar_trans_partial, const Index &iy_id, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Tensor3 &t_field, const EnergyLevelMap &nlte_field, const Tensor4 &vmr_field, const ArrayOfArrayOfSpeciesTag &abs_species, const Tensor3 &wind_u_field, const Tensor3 &wind_v_field, const Tensor3 &wind_w_field, const Tensor3 &mag_u_field, const Tensor3 &mag_v_field, const Tensor3 &mag_w_field, const Index &cloudbox_on, const String &iy_unit, const ArrayOfString &iy_aux_vars, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Ppath &ppath, const Vector &rte_pos2, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &iy_main_agenda, const Agenda &iy_space_agenda, const Agenda &iy_surface_agenda, const Agenda &iy_cloudbox_agenda, const Index &iy_agenda_call1, const Tensor3 &iy_transmission, const Numeric &rte_alonglos_v, const Tensor3 &surface_props_data, const Verbosity &verbosity)
WORKSPACE METHOD: iyEmissionStandard.
Definition: m_rte.cc:563
get_stepwise_scattersky_propmat
void get_stepwise_scattersky_propmat(StokesVector &ap, PropagationMatrix &Kp, ArrayOfStokesVector &dap_dx, ArrayOfPropagationMatrix &dKp_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, ConstMatrixView ppath_1p_pnd, const ArrayOfMatrix &ppath_dpnd_dx, const Index ppath_1p_id, const ArrayOfArrayOfSingleScatteringData &scat_data, ConstVectorView ppath_line_of_sight, ConstVectorView ppath_temperature, const Index &atmosphere_dim, const bool &jacobian_do)
Computes the contribution by scattering at propagation path point.
Definition: rte.cc:1591
joker
const Joker joker
Numeric
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
Matrix
The Matrix class.
Definition: matpackI.h:1193
ConstTensor7View::nlibraries
Index nlibraries() const
Returns the number of libraries.
Definition: matpackVII.cc:42
auto_md.h
get_diydx
void get_diydx(VectorView diy1, VectorView diy2, ConstMatrixView ImT, ConstMatrixView cumulative_transmission, ConstMatrixView dT1, ConstMatrixView dT2, ConstVectorView iYmJ, ConstVectorView dJ1, ConstVectorView dJ2, const Index stokes_dim, const bool transmission_only)
Definition: jacobian.cc:1049
RAD2DEG
const Numeric RAD2DEG
RadiationVector
Radiation Vector for Stokes dimension 1-4.
Definition: transmissionmatrix.h:395
do_temperature_jacobian
bool do_temperature_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants the temperature derivative.
Definition: jacobian.cc:1279
Tensor3::resize
void resize(Index p, Index r, Index c)
Resize function.
Definition: matpackIII.cc:664
mult
void mult(ComplexVectorView y, const ConstComplexMatrixView &M, const ConstComplexVectorView &x)
Matrix-Vector Multiplication.
Definition: complex.cc:1579
PI
const Numeric PI
doit.h
Radiative transfer in cloudbox.
Array
This can be used to make arrays out of anything.
Definition: array.h:40
update_radiation_vector
void update_radiation_vector(RadiationVector &I, ArrayOfRadiationVector &dI1, ArrayOfRadiationVector &dI2, const RadiationVector &J1, const RadiationVector &J2, const ArrayOfRadiationVector &dJ1, const ArrayOfRadiationVector &dJ2, const TransmissionMatrix &T, const TransmissionMatrix &PiT, const ArrayOfTransmissionMatrix &dT1, const ArrayOfTransmissionMatrix &dT2, const RadiativeTransferSolver solver)
Update the Radiation Vector.
Definition: transmissionmatrix.cc:1484
Vector::resize
void resize(Index n)
Resize function.
Definition: matpackI.cc:404
rtmethods_jacobian_finalisation
void rtmethods_jacobian_finalisation(Workspace &ws, ArrayOfTensor3 &diy_dx, ArrayOfTensor3 &diy_dpath, const Index &ns, const Index &nf, const Index &np, const Index &atmosphere_dim, const Ppath &ppath, const Vector &ppvar_p, const Vector &ppvar_t, const Matrix &ppvar_vmr, const Index &iy_agenda_call1, const Tensor3 &iy_transmission, const Agenda &water_p_eq_agenda, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfIndex jac_species_i, const ArrayOfIndex jac_is_t)
This function fixes the last steps to made on the Jacobian in some radiative transfer WSMs...
Definition: rte.cc:2416
get_stepwise_blackbody_radiation
void get_stepwise_blackbody_radiation(VectorView B, VectorView dB_dT, ConstVectorView ppath_f_grid, const Numeric &ppath_temperature, const bool &do_temperature_derivative)
Get the blackbody radiation at propagation path point.
Definition: rte.cc:1307
adapt_stepwise_partial_derivatives
void adapt_stepwise_partial_derivatives(ArrayOfPropagationMatrix &dK_dx, ArrayOfStokesVector &dS_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, ConstVectorView ppath_f_grid, ConstVectorView ppath_line_of_sight, ConstVectorView ppath_vmrs, const Numeric &ppath_temperature, const Numeric &ppath_pressure, const ArrayOfIndex &jacobian_species, const ArrayOfIndex &jacobian_wind, const Index &lte, const Index &atmosphere_dim, const bool &jacobian_do)
Adapts clearsky partial derivatives.
Definition: rte.cc:59
Ppath::np
Index np
Number of points describing the ppath.
Definition: ppath.h:52
ConstTensor7View::nshelves
Index nshelves() const
Returns the number of shelves.
Definition: matpackVII.cc:48
m_xml.h
Workspace methods and template functions for supergeneric XML IO.
ConstMatrixView
A constant view of a Matrix.
Definition: matpackI.h:982
ConstTensor7View::npages
Index npages() const
Returns the number of pages.
Definition: matpackVII.cc:54
Workspace
Workspace class.
Definition: workspace_ng.h:40
get_stepwise_scattersky_source
void get_stepwise_scattersky_source(StokesVector &Sp, ArrayOfStokesVector &dSp_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, ConstVectorView ppath_1p_pnd, const ArrayOfMatrix &ppath_dpnd_dx, const Index ppath_1p_id, const ArrayOfArrayOfSingleScatteringData &scat_data, ConstTensor7View cloudbox_field, ConstVectorView za_grid, ConstVectorView aa_grid, ConstMatrixView ppath_line_of_sight, const GridPos &ppath_pressure, const Vector &temperature, const Index &atmosphere_dim, const bool &jacobian_do, const Index &t_interp_order)
Calculates the stepwise scattering source terms.
Definition: rte.cc:1708
id_mat
void id_mat(MatrixView I)
Identity Matrix.
Definition: lin_alg.cc:2240
ConstTensor7View::nvitrines
Index nvitrines() const
Returns the number of vitrines.
Definition: matpackVII.cc:45
Ppath
The structure to describe a propagation path and releated quantities.
Definition: ppath.h:48
arts_omp.h
Header file for helper functions for OpenMP.
ConstTensor7View::ncols
Index ncols() const
Returns the number of columns.
Definition: matpackVII.cc:60
ConstTensor7View::nbooks
Index nbooks() const
Returns the number of books.
Definition: matpackVII.cc:51
Ppath::gp_p
ArrayOfGridPos gp_p
Index position with respect to the pressure grid.
Definition: ppath.h:82
Tensor4::resize
void resize(Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackIV.cc:1064
Tensor5
The Tensor5 class.
Definition: matpackV.h:506
rte.h
Declaration of functions in rte.cc.
Matrix::resize
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1056