ARTS 2.5.0 (git: 9ee3ac6c)
m_transmitter.cc
Go to the documentation of this file.
1/* Copyright (C) 2012
2 Patrick Eriksson <Patrick.Eriksson@chalmers.se>
3 Stefan Buehler <sbuehler(at)ltu.se>
4
5 This program is free software; you can redistribute it and/or modify it
6 under the terms of the GNU General Public License as published by the
7 Free Software Foundation; either version 2, or (at your option) any
8 later version.
9
10 This program is distributed in the hope that it will be useful,
11 but WITHOUT ANY WARRANTY; without even the implied warranty of
12 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 GNU General Public License for more details.
14
15 You should have received a copy of the GNU General Public License
16 along with this program; if not, write to the Free Software
17 Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
18 USA. */
19
20/*===========================================================================
21 === File description
22 ===========================================================================*/
23
33/*===========================================================================
34 === External declarations
35 ===========================================================================*/
36
37#include "arts.h"
38#include "auto_md.h"
39#include "geodetic.h"
40#include "jacobian.h"
41#include "lin_alg.h"
42#include "logic.h"
43#include "math_funcs.h"
44#include "matpack_complex.h"
45#include "messages.h"
46#include "rte.h"
47#include "sensor.h"
48#include <cmath>
49#include <stdexcept>
50
51extern const Numeric DEG2RAD;
52extern const Numeric PI;
53extern const Numeric RAD2DEG;
54extern const Numeric SPEED_OF_LIGHT;
55
56/* Workspace method: Doxygen documentation will be auto-generated */
57/* Has not been updated since v2.2
58void iyRadioLink(
59 Workspace& ws,
60 Matrix& iy,
61 ArrayOfTensor4& iy_aux,
62 Ppath& ppath,
63 ArrayOfTensor3& diy_dx,
64 const Index& stokes_dim,
65 const Vector& f_grid,
66 const Index& atmosphere_dim,
67 const Vector& p_grid,
68 const Vector& lat_grid,
69 const Vector& lon_grid,
70 const Tensor3& z_field,
71 const Tensor3& t_field,
72 const Tensor4& vmr_field,
73 const ArrayOfArrayOfSpeciesTag& abs_species,
74 const Tensor3& wind_u_field,
75 const Tensor3& wind_v_field,
76 const Tensor3& wind_w_field,
77 const Tensor3& mag_u_field,
78 const Tensor3& mag_v_field,
79 const Tensor3& mag_w_field,
80 const Vector& refellipsoid,
81 const Matrix& z_surface,
82 const Index& cloudbox_on,
83 const ArrayOfIndex& cloudbox_limits,
84 const Tensor4& pnd_field,
85 const ArrayOfArrayOfSingleScatteringData& scat_data,
86 const Matrix& particle_masses,
87 const ArrayOfString& iy_aux_vars,
88 const Index& jacobian_do,
89 const Agenda& ppath_agenda,
90 const Agenda& ppath_step_agenda,
91 const Agenda& propmat_clearsky_agenda,
92 const Agenda& iy_transmitter_agenda,
93 const Index& iy_agenda_call1,
94 const Tensor3& iy_transmittance,
95 const Vector& rte_pos,
96 const Vector& rte_los,
97 const Vector& rte_pos2,
98 const Numeric& rte_alonglos_v,
99 const Numeric& ppath_lmax,
100 const Numeric& ppath_lraytrace,
101 const Index& defocus_method,
102 const Numeric& defocus_shift,
103 const Verbosity& verbosity )
104{
105 // Throw error if unsupported features are requested
106 if( !iy_agenda_call1 )
107 throw runtime_error(
108 "Recursive usage not possible (iy_agenda_call1 must be 1)" );
109 if( !iy_transmittance.empty() )
110 throw runtime_error( "*iy_transmittance* must be empty" );
111 if( jacobian_do )
112 throw runtime_error( "This method does not provide any jacobians and "
113 "*jacobian_do* must be 0." );
114 if( defocus_method < 1 || defocus_method > 2 )
115 throw runtime_error( "Allowed choices for *defocus_method* is 1 and 2." );
116 diy_dx.resize(0);
117
118
119 //- Determine propagation path
120 ppath_agendaExecute( ws, ppath, ppath_lmax, ppath_lraytrace, rte_pos, rte_los,
121 rte_pos2, cloudbox_on, 0, t_field, z_field, vmr_field,
122 f_grid, ppath_agenda );
123
124 //- Set np to zero if ground intersection
125 const Index radback = ppath_what_background(ppath);
126 // np should already be 1 fon non-OK cases, but for extra safety ...
127 if( radback == 0 || radback == 2 )
128 { ppath.np = 1; }
129
130 // Some basic sizes
131 //
132 const Index nf = f_grid.nelem();
133 const Index ns = stokes_dim;
134 const Index np = ppath.np;
135
136
137 // The below copied from iyEmission. Activate for-loop when jacobians
138 // introduced.
139
140 // Set up variable with index of species where we want abs_per_species.
141 // This variable can below be extended in iy_aux part.
142 //
143 ArrayOfIndex iaps(0);
144 //
145 //for( Index i=0; i<jac_species_i.nelem(); i++ )
146 // {
147 // if( jac_species_i[i] >= 0 )
148 // { iaps.push_back( jac_species_i[i] ); }
149 // }
150
151 //=== iy_aux part ===========================================================
152 Index auxPressure = -1,
153 auxTemperature = -1,
154 auxAbsSum = -1,
155 auxPartExt = -1,
156 auxImpactParam = -1,
157 auxFreeSpaceLoss = -1,
158 auxFreeSpaceAtte = -1,
159 auxAtmosphericLoss = -1,
160 auxDefocusingLoss = -1,
161 auxFarRotTotal = -1,
162 auxFarRotSpeed = -1,
163 auxExtraPathDelay = -1,
164 auxBendingAngle = -1;
165 ArrayOfIndex auxAbsSpecies(0), auxAbsIsp(0);
166 ArrayOfIndex auxVmrSpecies(0), auxVmrIsp(0);
167 ArrayOfIndex auxPartCont(0), auxPartContI(0);
168 ArrayOfIndex auxPartField(0), auxPartFieldI(0);
169 //
170 const Index naux = iy_aux_vars.nelem();
171 iy_aux.resize( naux );
172 //
173 for( Index i=0; i<naux; i++ )
174 {
175 if( iy_aux_vars[i] == "Pressure" )
176 { auxPressure = i; iy_aux[i].resize( 1, 1, 1, np ); }
177 else if( iy_aux_vars[i] == "Temperature" )
178 { auxTemperature = i; iy_aux[i].resize( 1, 1, 1, np ); }
179 else if( iy_aux_vars[i].substr(0,13) == "VMR, species " )
180 {
181 Index ispecies;
182 istringstream is(iy_aux_vars[i].substr(13,2));
183 is >> ispecies;
184 if( ispecies < 0 || ispecies>=abs_species.nelem() )
185 {
186 ostringstream os;
187 os << "You have selected VMR of species with index "
188 << ispecies << ".\nThis species does not exist!";
189 throw runtime_error( os.str() );
190 }
191 auxVmrSpecies.push_back(i);
192 auxVmrIsp.push_back(ispecies);
193 iy_aux[i].resize( 1, 1, 1, np );
194 }
195 else if( iy_aux_vars[i] == "Absorption, summed" )
196 { auxAbsSum = i; iy_aux[i].resize( nf, ns, ns, np ); }
197 else if( iy_aux_vars[i] == "Particle extinction, summed" )
198 {
199 auxPartExt = i;
200 iy_aux[i].resize( nf, ns, ns, np );
201 iy_aux[i] = 0;
202 }
203 else if( iy_aux_vars[i].substr(0,20) == "Absorption, species " )
204 {
205 Index ispecies;
206 istringstream is(iy_aux_vars[i].substr(20,2));
207 is >> ispecies;
208 if( ispecies < 0 || ispecies>=abs_species.nelem() )
209 {
210 ostringstream os;
211 os << "You have selected absorption species with index "
212 << ispecies << ".\nThis species does not exist!";
213 throw runtime_error( os.str() );
214 }
215 auxAbsSpecies.push_back(i);
216 const Index ihit = find_first( iaps, ispecies );
217 if( ihit >= 0 )
218 { auxAbsIsp.push_back( ihit ); }
219 else
220 {
221 iaps.push_back(ispecies);
222 auxAbsIsp.push_back( iaps.nelem()-1 );
223 }
224 iy_aux[i].resize( nf, ns, ns, np );
225 }
226 else if( iy_aux_vars[i].substr(0,14) == "Mass content, " )
227 {
228 Index icont;
229 istringstream is(iy_aux_vars[i].substr(14,2));
230 is >> icont;
231 if( icont < 0 || icont>=particle_masses.ncols() )
232 {
233 ostringstream os;
234 os << "You have selected particle mass content category with "
235 << "index " << icont << ".\nThis category is not defined!";
236 throw runtime_error( os.str() );
237 }
238 auxPartCont.push_back(i);
239 auxPartContI.push_back(icont);
240 iy_aux[i].resize( 1, 1, 1, np );
241 }
242 else if( iy_aux_vars[i].substr(0,10) == "PND, type " )
243 {
244 Index ip;
245 istringstream is(iy_aux_vars[i].substr(10,2));
246 is >> ip;
247 if( ip < 0 || ip>=pnd_field.nbooks() )
248 {
249 ostringstream os;
250 os << "You have selected particle number density field with "
251 << "index " << ip << ".\nThis field is not defined!";
252 throw runtime_error( os.str() );
253 }
254 auxPartField.push_back(i);
255 auxPartFieldI.push_back(ip);
256 iy_aux[i].resize( 1, 1, 1, np );
257 }
258 else if( iy_aux_vars[i] == "Impact parameter" )
259 { auxImpactParam = i; iy_aux[i].resize( 1, 1, 1, 1 ); }
260 else if( iy_aux_vars[i] == "Free space loss" )
261 { auxFreeSpaceLoss = i; iy_aux[i].resize( 1, 1, 1, 1 ); }
262 else if( iy_aux_vars[i] == "Free space attenuation" )
263 { auxFreeSpaceAtte = i; iy_aux[i].resize( 1, 1, 1, np ); }
264 else if( iy_aux_vars[i] == "Atmospheric loss" )
265 { auxAtmosphericLoss = i; iy_aux[i].resize( nf, 1, 1, 1 ); }
266 else if( iy_aux_vars[i] == "Defocusing loss" )
267 { auxDefocusingLoss = i; iy_aux[i].resize( 1, 1, 1, 1 ); }
268 else if( iy_aux_vars[i] == "Faraday rotation" )
269 { auxFarRotTotal = i; iy_aux[i].resize( nf, 1, 1, 1 ); iy_aux[i] = 0; }
270 else if( iy_aux_vars[i] == "Faraday speed" )
271 { auxFarRotSpeed = i; iy_aux[i].resize( nf, 1, 1, np ); iy_aux[i] = 0; }
272 else if( iy_aux_vars[i] == "Extra path delay" )
273 { auxExtraPathDelay = i; iy_aux[i].resize( 1, 1, 1, 1 ); }
274 else if( iy_aux_vars[i] == "Bending angle" )
275 { auxBendingAngle = i; iy_aux[i].resize( 1, 1, 1, 1 ); }
276 else
277 {
278 ostringstream os;
279 os << "In *iy_aux_vars* you have included: \"" << iy_aux_vars[i]
280 << "\"\nThis choice is not recognised.";
281 throw runtime_error( os.str() );
282 }
283 }
284
285 // Special stuff to handle Faraday rotation
286 Index ife = -1; // When ready, ife refers abs_per_species
287 if( auxFarRotTotal>=0 || auxFarRotSpeed>=0 )
288 {
289 if( stokes_dim < 3 )
290 throw runtime_error(
291 "To include Faraday rotation, stokes_dim >= 3 is required." );
292
293 // Determine species index of free electrons
294 for( Index sp = 0; sp < abs_species.nelem() && ife < 0; sp++ )
295 {
296 if (abs_species[sp][0].Type() == SpeciesTag::TYPE_FREE_ELECTRONS)
297 { ife = sp; }
298 }
299 // If not found, then aux values already set to zero
300 if( ife < 0 )
301 {
302 auxFarRotTotal = -1;
303 auxFarRotSpeed = -1;
304 }
305 else
306 {
307 const Index ihit = find_first( iaps, ife );
308 if( ihit >= 0 )
309 { ife = ihit; }
310 else
311 {
312 iaps.push_back(ife);
313 ife = iaps.nelem() - 1;
314 }
315 }
316 }
317 //===========================================================================
318
319
320 // Handle cases whn no link was establsihed
321 if( radback == 0 || radback == 2 )
322 {
323 Numeric fillvalue = 0;
324 if( radback == 0 )
325 { fillvalue = NAN; }
326 //
327 iy.resize( nf, stokes_dim );
328 iy = fillvalue;
329 //
330 for( Index i=0; i<naux; i++ )
331 { iy_aux[i] = fillvalue; }
332 //
333 return;
334 }
335
336
337 // Transmitted signal
338 //
339 iy_transmitter_agendaExecute( ws, iy, f_grid,
340 ppath.pos(np-1,Range(0,atmosphere_dim)),
341 ppath.los(np-1,joker), iy_transmitter_agenda );
342 if( iy.ncols() != stokes_dim || iy.nrows() != nf )
343 { throw runtime_error( "The size of *iy* returned from "
344 "*iy_transmitter_agenda* is not correct." ); }
345
346
347 // Get atmospheric and attenuation quantities for each ppath point/step
348 //
349 Vector ppath_p, ppath_t;
350 Matrix ppath_vmr, ppath_pnd, ppath_mag, ppath_wind, ppath_f, ppath_t_nlte;
351 ArrayOfArrayOfPropagationMatrix abs_per_species, dummy_dppath_ext_dx;
352 ArrayOfPropagationMatrix ppath_ext, pnd_ext_mat;
353 Tensor4 trans_partial, trans_cumulat;
354 ArrayOfArrayOfStokesVector dummy_dppath_nlte_dx;
355 ArrayOfStokesVector dummy_ppath_nlte_source;
356 Vector scalar_tau;
357 ArrayOfIndex clear2cloudy, dummy_lte;
358 ArrayOfMatrix dummy_ppath_dpnd_dx;
359 ArrayOfTensor4 dummy_dpnd_field_dx;
360 const Tensor4 t_nlte_field_empty(0,0,0,0);
361 //
362 if( np > 1 )
363 {
364 get_ppath_atmvars( ppath_p, ppath_t, ppath_t_nlte, ppath_vmr,
365 ppath_wind, ppath_mag,
366 ppath, atmosphere_dim, p_grid, t_field, t_nlte_field_empty,
367 vmr_field, wind_u_field, wind_v_field, wind_w_field,
368 mag_u_field, mag_v_field, mag_w_field );
369
370 get_ppath_f( ppath_f, ppath, f_grid, atmosphere_dim,
371 rte_alonglos_v, ppath_wind );
372
373 get_ppath_pmat( ws, ppath_ext, dummy_ppath_nlte_source, dummy_lte, abs_per_species,
374 dummy_dppath_ext_dx, dummy_dppath_nlte_dx,
375 propmat_clearsky_agenda, ArrayOfRetrievalQuantity(0), ppath,
376 ppath_p, ppath_t, ppath_t_nlte, ppath_vmr, ppath_f,
377 ppath_mag, f_grid, stokes_dim, iaps );
378
379 if( !cloudbox_on )
380 {
381 ArrayOfArrayOfIndex extmat_case;
382 get_ppath_trans( trans_partial, extmat_case, trans_cumulat,
383 scalar_tau, ppath, ppath_ext, f_grid, stokes_dim );
384 }
385 else
386 {
387 // Extract basic scattering data
388 Array<ArrayOfArrayOfSingleScatteringData> scat_data_single;
389 ArrayOfArrayOfIndex extmat_case;
390 Tensor3 pnd_abs_vec;
391
392 get_ppath_cloudvars( clear2cloudy, ppath_pnd, dummy_ppath_dpnd_dx,
393 ppath, atmosphere_dim, cloudbox_limits,
394 pnd_field, dummy_dpnd_field_dx );
395 get_ppath_partopt( pnd_abs_vec, pnd_ext_mat, scat_data_single,
396 clear2cloudy, ppath_pnd, ppath, ppath_t,
397 stokes_dim, ppath_f, atmosphere_dim, scat_data,
398 verbosity );
399
400 get_ppath_trans2( trans_partial, extmat_case, trans_cumulat,
401 scalar_tau, ppath, ppath_ext, f_grid, stokes_dim,
402 clear2cloudy, pnd_ext_mat );
403 }
404 }
405
406 // Ppath length variables
407 //
408 Numeric lbg; // Bent geometrical length of ray path
409 Numeric lba; // Bent apparent length of ray path
410 //
411 lbg = ppath.end_lstep;
412 lba = lbg;
413
414 // Do RT calculations
415 //
416 if( np > 1 )
417 {
418 //=== iy_aux part =======================================================
419 // iy_aux for point np-1:
420 // Pressure
421 if( auxPressure >= 0 )
422 { iy_aux[auxPressure](0,0,0,np-1) = ppath_p[np-1]; }
423 // Temperature
424 if( auxTemperature >= 0 )
425 { iy_aux[auxTemperature](0,0,0,np-1) = ppath_t[np-1]; }
426 // VMR
427 for( Index j=0; j<auxVmrSpecies.nelem(); j++ )
428 { iy_aux[auxVmrSpecies[j]](0,0,0,np-1) = ppath_vmr(auxVmrIsp[j],np-1);}
429 // Absorption
430 if( auxAbsSum >= 0 )
431 { for( Index iv=0; iv<nf; iv++ ) {
432 for( Index is1=0; is1<ns; is1++ ){
433 for( Index is2=0; is2<ns; is2++ ){
434 iy_aux[auxAbsSum](iv,is1,is2,np-1) =
435 ppath_ext[np-1](iv,is1,is2); } } } }
436 for( Index j=0; j<auxAbsSpecies.nelem(); j++ )
437 { for( Index iv=0; iv<nf; iv++ ) {
438 for( Index is1=0; is1<stokes_dim; is1++ ){
439 for( Index is2=0; is2<stokes_dim; is2++ ){
440 iy_aux[auxAbsSpecies[j]](iv,is1,is2,np-1) =
441 abs_per_species[np-1][auxAbsIsp[j]](iv,is1,is2); } } } }
442 // Particle properties
443 if( cloudbox_on )
444 {
445 // Extinction
446 if( auxPartExt >= 0 && clear2cloudy[np-1] >= 0 )
447 {
448 const Index ic = clear2cloudy[np-1];
449 for( Index iv=0; iv<nf; iv++ ) {
450 for( Index is1=0; is1<ns; is1++ ){
451 for( Index is2=0; is2<ns; is2++ ){
452 iy_aux[auxPartExt](iv,is1,is2,np-1) =
453 pnd_ext_mat[ic](iv,is1,is2); } } }
454 }
455 // Particle mass content
456 for( Index j=0; j<auxPartCont.nelem(); j++ )
457 { iy_aux[auxPartCont[j]](0,0,0,np-1) = ppath_pnd(joker,np-1) *
458 particle_masses(joker,auxPartContI[j]); }
459 // Particle number density
460 for( Index j=0; j<auxPartField.nelem(); j++ )
461 { iy_aux[auxPartField[j]](0,0,0,np-1) =
462 ppath_pnd(auxPartFieldI[j],np-1); }
463 }
464 // Free space
465 if( auxFreeSpaceAtte >= 0 )
466 { iy_aux[auxFreeSpaceAtte](joker,0,0,np-1) = 2/lbg; }
467 // Faraday speed
468 if( auxFarRotSpeed >= 0 )
469 { for( Index iv=0; iv<nf; iv++ ) {
470 iy_aux[auxFarRotSpeed](iv,0,0,np-1) = 0.5 *
471 abs_per_species[np-1][ife](iv,1,2); } }
472 //=======================================================================
473
474 // Loop ppath steps
475 for( Index ip=np-2; ip>=0; ip-- )
476 {
477 // Lengths
478 lbg += ppath.lstep[ip];
479 lba += ppath.lstep[ip] * (ppath.ngroup[ip]+ppath.ngroup[ip+1]) / 2.0;
480
481 // Atmospheric loss of path step + Faraday rotation
482 if( stokes_dim == 1 )
483 {
484 for( Index iv=0; iv<nf; iv++ )
485 { iy(iv,0) = iy(iv,0) * trans_partial(iv,0,0,ip); }
486 }
487 else
488 {
489 for( Index iv=0; iv<nf; iv++ )
490 {
491 // Unpolarised:
492 if( is_diagonal( trans_partial(iv,joker,joker,ip) ) )
493 {
494 for( Index is=0; is<ns; is++ )
495 { iy(iv,is) = iy(iv,is) * trans_partial(iv,is,is,ip); }
496 }
497 // The general case:
498 else
499 {
500 Vector t1(ns);
501 mult( t1, trans_partial(iv,joker,joker,ip), iy(iv,joker));
502 iy(iv,joker) = t1;
503 }
504 }
505 }
506
507 //=== iy_aux part ===================================================
508 // Pressure
509 if( auxPressure >= 0 )
510 { iy_aux[auxPressure](0,0,0,ip) = ppath_p[ip]; }
511 // Temperature
512 if( auxTemperature >= 0 )
513 { iy_aux[auxTemperature](0,0,0,ip) = ppath_t[ip]; }
514 // VMR
515 for( Index j=0; j<auxVmrSpecies.nelem(); j++ )
516 { iy_aux[auxVmrSpecies[j]](0,0,0,ip) = ppath_vmr(auxVmrIsp[j],ip);}
517 // Absorption
518 if( auxAbsSum >= 0 )
519 { for( Index iv=0; iv<nf; iv++ ) {
520 for( Index is1=0; is1<ns; is1++ ){
521 for( Index is2=0; is2<ns; is2++ ){
522 iy_aux[auxAbsSum](iv,is1,is2,ip) =
523 ppath_ext[ip](iv,is1,is2); } } } }
524 for( Index j=0; j<auxAbsSpecies.nelem(); j++ )
525 { for( Index iv=0; iv<nf; iv++ ) {
526 for( Index is1=0; is1<stokes_dim; is1++ ){
527 for( Index is2=0; is2<stokes_dim; is2++ ){
528 iy_aux[auxAbsSpecies[j]](iv,is1,is2,ip) =
529 abs_per_species[ip][auxAbsIsp[j]](iv,is1,is2); } } } }
530 // Particle properties
531 if( cloudbox_on )
532 {
533 // Extinction
534 if( auxPartExt >= 0 && clear2cloudy[ip] >= 0 )
535 {
536 const Index ic = clear2cloudy[ip];
537 for( Index iv=0; iv<nf; iv++ ) {
538 for( Index is1=0; is1<ns; is1++ ){
539 for( Index is2=0; is2<ns; is2++ ){
540 iy_aux[auxPartExt](iv,is1,is2,ip) =
541 pnd_ext_mat[ic](iv,is1,is2); } } }
542 }
543 // Particle mass content
544 for( Index j=0; j<auxPartCont.nelem(); j++ )
545 { iy_aux[auxPartCont[j]](0,0,0,ip) = ppath_pnd(joker,ip) *
546 particle_masses(joker,auxPartContI[j]); }
547 // Particle number density
548 for( Index j=0; j<auxPartField.nelem(); j++ )
549 { iy_aux[auxPartField[j]](0,0,0,ip) =
550 ppath_pnd(auxPartFieldI[j],ip); }
551 }
552 // Free space loss
553 if( auxFreeSpaceAtte >= 0 )
554 { iy_aux[auxFreeSpaceAtte](joker,0,0,ip) = 2/lbg; }
555 // Faraday rotation, total
556 if( auxFarRotTotal >= 0 )
557 { for( Index iv=0; iv<nf; iv++ ) {
558 iy_aux[auxFarRotTotal](iv,0,0,0) += RAD2DEG * ppath.lstep[ip] *
559 0.25 * ( abs_per_species[ip][ife](iv,1,2)+
560 abs_per_species[ip+1][ife](iv,1,2)); } }
561 // Faraday speed
562 if( auxFarRotSpeed >= 0 )
563 { for( Index iv=0; iv<nf; iv++ ) {
564 iy_aux[auxFarRotSpeed](iv,0,0,ip) = 0.5 *
565 abs_per_species[ip][ife](iv,1,2); } }
566 //===================================================================
567 }
568
569
570 //=== iy_aux part =======================================================
571 if( auxAtmosphericLoss >= 0 )
572 { iy_aux[auxAtmosphericLoss](joker,0,0,0) = iy(joker,0); }
573 if( auxImpactParam >= 0 )
574 {
575 ARTS_ASSERT( ppath.constant >= 0 );
576 iy_aux[auxImpactParam](joker,0,0,0) = ppath.constant;
577 }
578 //=======================================================================
579
580
581 // Remaing length of ppath
582 lbg += ppath.start_lstep;
583 lba += ppath.start_lstep;
584
585
586 // Determine total free space loss
587 Numeric fspl = 1 / ( 4 * PI * lbg*lbg );
588 //
589 if( auxFreeSpaceLoss >= 0 )
590 { iy_aux[auxFreeSpaceLoss] = fspl; }
591
592
593 // Determine defocusing loss
594 Numeric dfl = 1;
595 if( defocus_method == 1 )
596 {
597 defocusing_general( ws, dfl, ppath_step_agenda, atmosphere_dim,
598 p_grid, lat_grid, lon_grid, t_field, z_field,
599 vmr_field, f_grid, refellipsoid,
600 z_surface, ppath, ppath_lmax, ppath_lraytrace,
601 defocus_shift, verbosity );
602 }
603 else if( defocus_method == 2 )
604 { defocusing_sat2sat( ws, dfl, ppath_step_agenda, atmosphere_dim,
605 p_grid, lat_grid, lon_grid, t_field, z_field,
606 vmr_field, f_grid, refellipsoid,
607 z_surface, ppath, ppath_lmax, ppath_lraytrace,
608 defocus_shift, verbosity );
609 }
610 if( auxDefocusingLoss >= 0 )
611 { iy_aux[auxDefocusingLoss] = dfl; }
612
613
614
615 // Include free space and defocusing losses
616 iy *= fspl*dfl;
617
618
619 // Extra path delay
620 if( auxExtraPathDelay >= 0 )
621 {
622 // Radius of rte_pos and rte_pos2
623 const Numeric r1 = ppath.end_pos[0] +
624 pos2refell_r( atmosphere_dim, refellipsoid,
625 lat_grid, lon_grid, ppath.end_pos );
626 const Numeric r2 = ppath.start_pos[0] +
627 pos2refell_r( atmosphere_dim, refellipsoid,
628 lat_grid, lon_grid, ppath.start_pos );
629
630 // Geometrical distance between start and end point
631 Numeric lgd ;
632 if( atmosphere_dim <= 2 )
633 { distance2D( lgd, r1, ppath.end_pos[1], r2, ppath.start_pos[1] ); }
634 else
635 { distance3D( lgd, r1, ppath.end_pos[1], ppath.end_pos[2],
636 r2, ppath.start_pos[1], ppath.start_pos[2] ); }
637 //
638 iy_aux[auxExtraPathDelay] = ( lba - lgd ) / SPEED_OF_LIGHT;
639 }
640
641
642 // Bending angle
643 if( auxBendingAngle >= 0 )
644 {
645 Numeric ba = -999;
646 bending_angle1d( ba, ppath );
647 //
648 iy_aux[auxBendingAngle] = ba;
649 }
650 }
651}
652*/
653
654/* Workspace method: Doxygen documentation will be auto-generated */
656 Matrix& iy,
657 ArrayOfMatrix& iy_aux,
658 ArrayOfTensor3& diy_dx,
659 Vector& ppvar_p,
660 Vector& ppvar_t,
661 EnergyLevelMap& ppvar_nlte,
662 Matrix& ppvar_vmr,
663 Matrix& ppvar_wind,
664 Matrix& ppvar_mag,
665 Matrix& ppvar_pnd,
666 Matrix& ppvar_f,
667 Tensor3& ppvar_iy,
668 Tensor4& ppvar_trans_cumulat,
669 Tensor4& ppvar_trans_partial,
670 const Index& stokes_dim,
671 const Vector& f_grid,
672 const Index& atmosphere_dim,
673 const Vector& p_grid,
674 const Tensor3& t_field,
675 const EnergyLevelMap& nlte_field,
676 const Tensor4& vmr_field,
677 const ArrayOfArrayOfSpeciesTag& abs_species,
678 const Tensor3& wind_u_field,
679 const Tensor3& wind_v_field,
680 const Tensor3& wind_w_field,
681 const Tensor3& mag_u_field,
682 const Tensor3& mag_v_field,
683 const Tensor3& mag_w_field,
684 const Index& cloudbox_on,
685 const ArrayOfIndex& cloudbox_limits,
686 const Tensor4& pnd_field,
687 const ArrayOfTensor4& dpnd_field_dx,
688 const ArrayOfString& scat_species,
689 const ArrayOfArrayOfSingleScatteringData& scat_data,
690 const ArrayOfString& iy_aux_vars,
691 const Index& jacobian_do,
692 const ArrayOfRetrievalQuantity& jacobian_quantities,
693 const Ppath& ppath,
694 const Agenda& propmat_clearsky_agenda,
695 const Agenda& water_p_eq_agenda,
696 const Agenda& iy_transmitter_agenda,
697 const Index& iy_agenda_call1,
698 const Tensor3& iy_transmittance,
699 const Numeric& rte_alonglos_v,
700 const Verbosity&) {
701 // Init Jacobian quantities?
702 const Index j_analytical_do = jacobian_do ? do_analytical_jacobian<1>(jacobian_quantities) : 0;
703
704 // Some basic sizes
705 const Index nf = f_grid.nelem();
706 const Index ns = stokes_dim;
707 const Index np = ppath.np;
708 const Index nq = j_analytical_do ? jacobian_quantities.nelem() : 0;
709
710 // Radiative background index
711 const Index rbi = ppath_what_background(ppath);
712
713 // Checks of input
714 // Throw error if unsupported features are requested
715 if (!iy_agenda_call1)
716 throw runtime_error(
717 "Recursive usage not possible (iy_agenda_call1 must be 1)");
718 if (!iy_transmittance.empty())
719 throw runtime_error("*iy_transmittance* must be empty");
720 if (rbi < 1 || rbi > 9)
721 throw runtime_error(
722 "ppath.background is invalid. Check your "
723 "calculation of *ppath*?");
724 if (jacobian_do) {
725 if (dpnd_field_dx.nelem() != jacobian_quantities.nelem())
726 throw runtime_error(
727 "*dpnd_field_dx* not properly initialized:\n"
728 "Number of elements in dpnd_field_dx must be equal number of jacobian"
729 " quantities.\n(Note: jacobians have to be defined BEFORE *pnd_field*"
730 " is calculated/set.");
731 }
732 // iy_aux_vars checked below
733
734 // Transmitted signal
736 iy,
737 f_grid,
738 ppath.pos(np - 1, Range(0, atmosphere_dim)),
739 ppath.los(np - 1, joker),
740 iy_transmitter_agenda);
741 if (iy.ncols() != ns || iy.nrows() != nf) {
742 ostringstream os;
743 os << "The size of *iy* returned from *iy_transmitter_agenda* is\n"
744 << "not correct:\n"
745 << " expected size = [" << nf << "," << stokes_dim << "]\n"
746 << " size of iy = [" << iy.nrows() << "," << iy.ncols() << "]\n";
747 throw runtime_error(os.str());
748 }
749
750 // Set diy_dpath if we are doing are doing jacobian calculations
751 ArrayOfTensor3 diy_dpath = j_analytical_do ? get_standard_diy_dpath(jacobian_quantities, np, nf, ns, false) : ArrayOfTensor3(0);
752
753 // Set the species pointers if we are doing jacobian
754 const ArrayOfIndex jac_species_i = j_analytical_do ? get_pointers_for_analytical_species(jacobian_quantities, abs_species) : ArrayOfIndex(0);
755
756 // Start diy_dx out if we are doing the first run and are doing jacobian calculations
757 if (j_analytical_do and iy_agenda_call1) diy_dx = get_standard_starting_diy_dx(jacobian_quantities, np, nf, ns, false);
758
759 // Checks that the scattering species are treated correctly if their derivatives are needed (we can here discard the Array)
760 if (j_analytical_do and iy_agenda_call1) get_pointers_for_scat_species(jacobian_quantities, scat_species, cloudbox_on);
761
762 // Init iy_aux and fill where possible
763 const Index naux = iy_aux_vars.nelem();
764 iy_aux.resize(naux);
765 //
766 Index auxOptDepth = -1;
767 //
768 for (Index i = 0; i < naux; i++) {
769 iy_aux[i].resize(nf, ns);
770 iy_aux[i] = 0;
771
772 if (iy_aux_vars[i] == "Radiative background")
773 iy_aux[i](joker, 0) = (Numeric)min((Index)2, rbi - 1);
774 else if (iy_aux_vars[i] == "Optical depth")
775 auxOptDepth = i;
776 else {
777 ostringstream os;
778 os << "The only allowed strings in *iy_aux_vars* are:\n"
779 << " \"Radiative background\"\n"
780 << " \"Optical depth\"\n"
781 << "but you have selected: \"" << iy_aux_vars[i] << "\"";
782 throw runtime_error(os.str());
783 }
784 }
785
786 // Get atmospheric and radiative variables along the propagation path
787 //
788 ppvar_trans_cumulat.resize(np, nf, ns, ns);
789 ppvar_trans_partial.resize(np, nf, ns, ns);
790
791 ArrayOfRadiationVector lvl_rad(np, RadiationVector(nf, ns));
794
800
801 ArrayOfIndex clear2cloudy;
802 //
803 if (np == 1 && rbi == 1) { // i.e. ppath is totally outside the atmosphere:
804 ppvar_p.resize(0);
805 ppvar_t.resize(0);
806 ppvar_vmr.resize(0, 0);
807 ppvar_wind.resize(0, 0);
808 ppvar_mag.resize(0, 0);
809 ppvar_pnd.resize(0, 0);
810 ppvar_f.resize(0, 0);
811 ppvar_iy.resize(0, 0, 0);
812 ppvar_trans_cumulat = 0;
813 ppvar_trans_partial = 0;
814 for (Index iv = 0; iv < nf; iv++) {
815 for (Index is = 0; is < ns; is++) {
816 ppvar_trans_cumulat(0,iv,is,is) = 1;
817 ppvar_trans_partial(0,iv,is,is) = 1;
818 }
819 }
820
821 } else {
822 // ppvar_iy
823 ppvar_iy.resize(nf, ns, np);
824 ppvar_iy(joker, joker, np - 1) = iy;
825
826 // Basic atmospheric variables
827 get_ppath_atmvars(ppvar_p,
828 ppvar_t,
829 ppvar_nlte,
830 ppvar_vmr,
831 ppvar_wind,
832 ppvar_mag,
833 ppath,
834 atmosphere_dim,
835 p_grid,
836 t_field,
837 nlte_field,
838 vmr_field,
839 wind_u_field,
840 wind_v_field,
841 wind_w_field,
842 mag_u_field,
843 mag_v_field,
844 mag_w_field);
845
847 ppvar_f, ppath, f_grid, atmosphere_dim, rte_alonglos_v, ppvar_wind);
848
849 // pnd_field
850 ArrayOfMatrix ppvar_dpnd_dx(0);
851 //
852 if (cloudbox_on)
853 get_ppath_cloudvars(clear2cloudy,
854 ppvar_pnd,
855 ppvar_dpnd_dx,
856 ppath,
857 atmosphere_dim,
858 cloudbox_limits,
859 pnd_field,
860 dpnd_field_dx);
861 else {
862 clear2cloudy.resize(np);
863 for (Index ip = 0; ip < np; ip++) clear2cloudy[ip] = -1;
864 }
865
866 // Size radiative variables always used
867 PropagationMatrix K_this(nf, ns), K_past(nf, ns), Kp(nf, ns);
868 StokesVector a(nf, ns), S(nf, ns), Sp(nf, ns);
869 ArrayOfIndex lte(np);
870
871 // Init variables only used if analytical jacobians done
872 Vector dB_dT(0);
873 ArrayOfPropagationMatrix dK_this_dx(nq), dK_past_dx(nq), dKp_dx(nq);
874 ArrayOfStokesVector da_dx(nq), dS_dx(nq), dSp_dx(nq);
875 //
876 Index temperature_derivative_position = -1;
877 bool do_hse = false;
878
879 if (j_analytical_do) {
880 dB_dT.resize(nf);
882 dK_past_dx[iq] = PropagationMatrix(nf, ns);
883 dKp_dx[iq] = PropagationMatrix(nf, ns);
884 da_dx[iq] = StokesVector(nf, ns);
885 dS_dx[iq] = StokesVector(nf, ns);
886 dSp_dx[iq] = StokesVector(nf, ns);
887 if (jacobian_quantities[iq] == Jacobian::Atm::Temperature) {
888 temperature_derivative_position = iq;
889 do_hse = jacobian_quantities[iq].Subtag() ==
890 "HSE on";
891 })
892 }
893
894 // Loop ppath points and determine radiative properties
895 for (Index ip = 0; ip < np; ip++) {
897 K_this,
898 S,
899 lte[ip],
900 dK_this_dx,
901 dS_dx,
902 propmat_clearsky_agenda,
903 jacobian_quantities,
904 ppvar_f(joker, ip),
905 ppvar_mag(joker, ip),
906 ppath.los(ip, joker),
907 ppvar_nlte[ip],
908 ppvar_vmr(joker, ip),
909 ppvar_t[ip],
910 ppvar_p[ip],
911 j_analytical_do);
912
913 if (j_analytical_do) {
915 dS_dx,
916 jacobian_quantities,
917 ppvar_f(joker, ip),
918 ppath.los(ip, joker),
919 lte[ip],
920 atmosphere_dim,
921 j_analytical_do);
922 }
923
924 if (clear2cloudy[ip] + 1) {
926 Kp,
927 da_dx,
928 dKp_dx,
929 jacobian_quantities,
930 ppvar_pnd(joker, Range(ip, 1)),
931 ppvar_dpnd_dx,
932 ip,
933 scat_data,
934 ppath.los(ip, joker),
935 ppvar_t[Range(ip, 1)],
936 atmosphere_dim,
937 jacobian_do);
938 K_this += Kp;
939
940 if (j_analytical_do) {
941 FOR_ANALYTICAL_JACOBIANS_DO(dK_this_dx[iq] += dKp_dx[iq];)
942 }
943 }
944
945 // Transmission
946 if (ip not_eq 0) {
947 const Numeric dr_dT_past =
948 do_hse ? ppath.lstep[ip - 1] / (2.0 * ppvar_t[ip - 1]) : 0;
949 const Numeric dr_dT_this =
950 do_hse ? ppath.lstep[ip - 1] / (2.0 * ppvar_t[ip]) : 0;
951 stepwise_transmission(lyr_tra[ip],
952 dlyr_tra_above[ip],
953 dlyr_tra_below[ip],
954 K_past,
955 K_this,
956 dK_past_dx,
957 dK_this_dx,
958 ppath.lstep[ip - 1],
959 dr_dT_past,
960 dr_dT_this,
961 temperature_derivative_position);
962 }
963
964 swap(K_past, K_this);
965 swap(dK_past_dx, dK_this_dx);
966 }
967 }
968
969 const ArrayOfTransmissionMatrix tot_tra =
971
972 // iy_aux: Optical depth
973 if (auxOptDepth >= 0) {
974 for (Index iv = 0; iv < nf; iv++)
975 iy_aux[auxOptDepth](iv, 0) = -std::log(tot_tra[np - 1](iv, 0, 0));
976 }
977
978 lvl_rad[np - 1] = iy;
979
980 // Radiative transfer calculations
981 for (Index ip = np - 2; ip >= 0; ip--) {
982 lvl_rad[ip] = lvl_rad[ip + 1];
983 update_radiation_vector(lvl_rad[ip],
984 dlvl_rad[ip],
985 dlvl_rad[ip + 1],
990 lyr_tra[ip + 1],
991 tot_tra[ip],
992 dlyr_tra_above[ip + 1],
993 dlyr_tra_below[ip + 1],
998 Numeric(),
999 Vector(),
1000 Vector(),
1001 0,
1002 0,
1004 }
1005
1006 // Copy back to ARTS external style
1007 iy = lvl_rad[0];
1008 //
1009 if (np > 1) {
1010 for (Index ip = 0; ip < lvl_rad.nelem(); ip++) {
1011 ppvar_trans_cumulat(ip, joker, joker, joker) = tot_tra[ip];
1012 ppvar_trans_partial(ip, joker, joker, joker) = lyr_tra[ip];
1013 ppvar_iy(joker, joker, ip) = lvl_rad[ip];
1014 if (j_analytical_do)
1015 FOR_ANALYTICAL_JACOBIANS_DO(diy_dpath[iq](ip, joker, joker) =
1016 dlvl_rad[ip][iq];);
1017 }
1018 }
1019
1020 // Finalize analytical Jacobians
1021 if (j_analytical_do) {
1023 diy_dx,
1024 diy_dpath,
1025 ns,
1026 nf,
1027 np,
1028 atmosphere_dim,
1029 ppath,
1030 ppvar_p,
1031 ppvar_t,
1032 ppvar_vmr,
1033 iy_agenda_call1,
1034 iy_transmittance,
1035 water_p_eq_agenda,
1036 jacobian_quantities,
1037 jac_species_i);
1038 }
1039}
1040
1041/* Workspace method: Doxygen documentation will be auto-generated */
1043 const Index& stokes_dim,
1044 const Vector& f_grid,
1045 const ArrayOfIndex& instrument_pol,
1046 const Verbosity&) {
1047 const Index nf = f_grid.nelem();
1048
1049 if (instrument_pol.nelem() != nf)
1050 throw runtime_error(
1051 "The length of *f_grid* and the number of elements "
1052 "in *instrument_pol* must be equal.");
1053
1054 iy.resize(nf, stokes_dim);
1055
1056 for (Index i = 0; i < nf; i++) {
1057 stokes2pol(iy(i, joker), stokes_dim, instrument_pol[i], 1);
1058 }
1059}
1060
1061/* Workspace method: Doxygen documentation will be auto-generated */
1063 const Index& stokes_dim,
1064 const Vector& f_grid,
1065 const ArrayOfIndex& instrument_pol,
1066 const Verbosity&) {
1067 const Index nf = f_grid.nelem();
1068
1069 if (instrument_pol.nelem() != 1)
1070 throw runtime_error(
1071 "The number of elements in *instrument_pol* must be 1.");
1072
1073 iy.resize(nf, stokes_dim);
1074
1075 stokes2pol(iy(0, joker), stokes_dim, instrument_pol[0], 1);
1076
1077 for (Index i = 1; i < nf; i++) {
1078 iy(i, joker) = iy(0, joker);
1079 }
1080}
Array< Index > ArrayOfIndex
An array of Index.
Definition: array.h:39
Array< Tensor3 > ArrayOfTensor3
An array of Tensor3.
Definition: array.h:73
The global header file for ARTS.
void iy_transmitter_agendaExecute(Workspace &ws, Matrix &iy, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Agenda &input_agenda)
Definition: auto_md.cc:23863
The Agenda class.
Definition: agenda_class.h:44
This can be used to make arrays out of anything.
Definition: array.h:107
Index nelem() const ARTS_NOEXCEPT
Number of elements.
Definition: array.h:195
Index nrows() const ARTS_NOEXCEPT
Returns the number of rows.
Definition: matpackI.cc:433
Index ncols() const ARTS_NOEXCEPT
Returns the number of columns.
Definition: matpackI.cc:436
bool empty() const
Check if variable is empty.
Definition: matpackIII.cc:36
Index nelem() const ARTS_NOEXCEPT
Returns the number of elements.
Definition: matpackI.cc:48
The Matrix class.
Definition: matpackI.h:1225
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1056
Radiation Vector for Stokes dimension 1-4.
The range class.
Definition: matpackI.h:165
Stokes vector is as Propagation matrix but only has 4 possible values.
The Tensor3 class.
Definition: matpackIII.h:339
void resize(Index p, Index r, Index c)
Resize function.
Definition: matpackIII.cc:658
The Tensor4 class.
Definition: matpackIV.h:421
void resize(Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackIV.cc:1058
Class to keep track of Transmission Matrices for Stokes Dim 1-4.
The Vector class.
Definition: matpackI.h:876
void resize(Index n)
Resize function.
Definition: matpackI.cc:408
Workspace class.
Definition: workspace_ng.h:40
ArrayOfIndex get_pointers_for_analytical_species(const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfArrayOfSpeciesTag &abs_species)
Help function for analytical jacobian calculations.
Definition: jacobian.cc:564
ArrayOfTensor3 get_standard_starting_diy_dx(const ArrayOfRetrievalQuantity &jacobian_quantities, Index np, Index nf, Index ns, bool active)
Help function for analytical jacobian calculations.
Definition: jacobian.cc:605
ArrayOfTensor3 get_standard_diy_dpath(const ArrayOfRetrievalQuantity &jacobian_quantities, Index np, Index nf, Index ns, bool active)
Help function for analytical jacobian calculations.
Definition: jacobian.cc:596
ArrayOfIndex get_pointers_for_scat_species(const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfString &scat_species, const bool cloudbox_on)
Help function for analytical jacobian calculations.
Definition: jacobian.cc:618
Routines for setting up the jacobian.
#define FOR_ANALYTICAL_JACOBIANS_DO(what_to_do)
Definition: jacobian.h:563
#define ns
#define min(a, b)
Linear algebra functions.
Header file for logic.cc.
void iy_transmitterSinglePol(Matrix &iy, const Index &stokes_dim, const Vector &f_grid, const ArrayOfIndex &instrument_pol, const Verbosity &)
WORKSPACE METHOD: iy_transmitterSinglePol.
const Numeric RAD2DEG
void iyTransmissionStandard(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, Tensor4 &ppvar_trans_partial, 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 ArrayOfString &iy_aux_vars, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Ppath &ppath, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &iy_transmitter_agenda, const Index &iy_agenda_call1, const Tensor3 &iy_transmittance, const Numeric &rte_alonglos_v, const Verbosity &)
WORKSPACE METHOD: iyTransmissionStandard.
const Numeric PI
const Numeric SPEED_OF_LIGHT
const Numeric DEG2RAD
void iy_transmitterMultiplePol(Matrix &iy, const Index &stokes_dim, const Vector &f_grid, const ArrayOfIndex &instrument_pol, const Verbosity &)
WORKSPACE METHOD: iy_transmitterMultiplePol.
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
void swap(ComplexVector &v1, ComplexVector &v2)
Swaps two objects.
const Joker joker
Declarations having to do with the four output streams.
Temperature
Definition: jacobian.h:57
invlib::Vector< ArtsVector > Vector
invlib wrapper type for ARTS vectors.
Definition: oem.h:31
Index ppath_what_background(const Ppath &ppath)
Returns the case number for the radiative background.
Definition: ppath.cc:1476
#define a
Array< PropagationMatrix > ArrayOfPropagationMatrix
void get_stepwise_scattersky_propmat(StokesVector &ap, PropagationMatrix &Kp, ArrayOfStokesVector &dap_dx, ArrayOfPropagationMatrix &dKp_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, const ConstMatrixView &ppath_1p_pnd, const ArrayOfMatrix &ppath_dpnd_dx, const Index ppath_1p_id, const ArrayOfArrayOfSingleScatteringData &scat_data, const ConstVectorView &ppath_line_of_sight, const ConstVectorView &ppath_temperature, const Index &atmosphere_dim, const bool &jacobian_do)
Computes the contribution by scattering at propagation path point.
Definition: rte.cc:1316
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, const ConstVectorView &ppath_f_grid, const ConstVectorView &ppath_magnetic_field, const ConstVectorView &ppath_line_of_sight, const EnergyLevelMap &ppath_nlte, const ConstVectorView &ppath_vmrs, const Numeric &ppath_temperature, const Numeric &ppath_pressure, const bool &jacobian_do)
Gets the clearsky propgation matrix and NLTE contributions.
Definition: rte.cc:1131
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, const ConstVectorView &p_grid, const ConstTensor3View &t_field, const EnergyLevelMap &nlte_field, const ConstTensor4View &vmr_field, const ConstTensor3View &wind_u_field, const ConstTensor3View &wind_v_field, const ConstTensor3View &wind_w_field, const ConstTensor3View &mag_u_field, const ConstTensor3View &mag_v_field, const ConstTensor3View &mag_w_field)
Determines pressure, temperature, VMR, winds and magnetic field for each propgataion path point.
Definition: rte.cc:843
void get_ppath_f(Matrix &ppath_f, const Ppath &ppath, const ConstVectorView &f_grid, const Index &atmosphere_dim, const Numeric &rte_alonglos_v, const ConstMatrixView &ppath_wind)
Determines the Doppler shifted frequencies along the propagation path.
Definition: rte.cc:1066
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_transmittance, const Agenda &water_p_eq_agenda, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfIndex &jac_species_i)
This function fixes the last steps to made on the Jacobian in some radiative transfer WSMs.
Definition: rte.cc:2068
void adapt_stepwise_partial_derivatives(ArrayOfPropagationMatrix &dK_dx, ArrayOfStokesVector &dS_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, const ConstVectorView &ppath_f_grid, const ConstVectorView &ppath_line_of_sight, const Index &lte, const Index &atmosphere_dim, const bool &jacobian_do)
Adapts clearsky partial derivatives.
Definition: rte.cc:55
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:962
Declaration of functions in rte.cc.
void stokes2pol(VectorView w, const Index &stokes_dim, const Index &ipol_1based, const Numeric nv)
stokes2pol
Definition: sensor.cc:1057
Sensor modelling functions.
The structure to describe a propagation path and releated quantities.
Definition: ppath.h:48
Matrix los
Line-of-sight at each ppath point.
Definition: ppath.h:66
Index np
Number of points describing the ppath.
Definition: ppath.h:52
Matrix pos
The distance between start pos and the last position in pos.
Definition: ppath.h:64
Vector lstep
The length between ppath points.
Definition: ppath.h:70
ArrayOfTransmissionMatrix cumulative_transmission(const ArrayOfTransmissionMatrix &T, const CumulativeTransmission type)
Accumulate the transmission matrix over all layers.
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 PropagationMatrix &K1, const PropagationMatrix &K2, const ArrayOfPropagationMatrix &dK1, const ArrayOfPropagationMatrix &dK2, const Numeric r, const Vector &dr1, const Vector &dr2, const Index ia, const Index iz, const RadiativeTransferSolver solver)
Update the Radiation Vector.
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.
Array< RadiationVector > ArrayOfRadiationVector
Array< TransmissionMatrix > ArrayOfTransmissionMatrix