ARTS 2.5.11 (git: 6827797f)
m_surface.cc
Go to the documentation of this file.
1/*===========================================================================
2 === File description
3 ===========================================================================*/
4
16/*===========================================================================
17 === External declarations
18 ===========================================================================*/
19
20#include <cmath>
21#include "array.h"
22#include "arts.h"
23#include "arts_constants.h"
24#include "auto_md.h"
25#include "check_input.h"
26#include "matpack_complex.h"
27#include "fastem.h"
28#include "geodetic.h"
29#include "interpolation.h"
30#include "math_funcs.h"
31#include "matpack_math.h"
32#include "messages.h"
33#include "physics_funcs.h"
34#include "ppath.h"
35#include "rte.h"
36#include "special_interp.h"
37#include "surface.h"
38#include "tessem.h"
39#include "arts_conversions.h"
40#include "gas_scattering.h"
41
42using Constant::pi;
44
45inline constexpr Numeric EARTH_RADIUS=Constant::earth_radius;
46inline constexpr Numeric DEG2RAD = Conversion::deg2rad(1);
47
48/*===========================================================================
49 === The functions (in alphabetical order)
50 ===========================================================================*/
51
52/* Workspace method: Doxygen documentation will be auto-generated */
53void FastemStandAlone(Matrix& emissivity,
54 Matrix& reflectivity,
55 const Vector& f_grid,
56 const Numeric& surface_skin_t,
57 const Numeric& za,
58 const Numeric& salinity,
59 const Numeric& wind_speed,
60 const Numeric& rel_aa,
61 const Vector& transmittance,
62 const Index& fastem_version,
63 const Verbosity&) {
64 const Index nf = f_grid.nelem();
65
66 chk_if_in_range("zenith angle", za, 90, 180);
67 chk_if_in_range_exclude("surface skin temperature", surface_skin_t, 260, 373);
68 chk_if_in_range_exclude_high("salinity", salinity, 0, 1);
69 chk_if_in_range_exclude_high("wind speed", wind_speed, 0, 100);
70 chk_if_in_range("azimuth angle", rel_aa, -180, 180);
71 chk_vector_length("transmittance", "f_grid", transmittance, f_grid);
72 if (fastem_version < 3 || fastem_version > 6)
73 throw std::runtime_error(
74 "Invalid fastem version: 3 <= fastem_version <= 6");
75
76 emissivity.resize(nf, 4);
77 reflectivity.resize(nf, 4);
78
79 const Numeric t = max(surface_skin_t, Numeric(270));
80
81 for (Index i = 0; i < nf; i++) {
82 if (f_grid[i] > 250e9)
83 throw std::runtime_error("Only frequency <= 250 GHz are allowed");
84 chk_if_in_range("transmittance", transmittance[i], 0, 1);
85
86 Vector e, r;
87
88 fastem(e,
89 r,
90 f_grid[i],
91 za,
92 t,
93 salinity,
94 wind_speed,
95 transmittance[i],
96 rel_aa,
97 fastem_version);
98
99 emissivity(i, joker) = e;
100 reflectivity(i, joker) = r;
101 }
102
103 // FASTEM does not work close to the horizon (at least v6). Make sure values
104 // are inside [0,1]. Then seems best to make sure that e+r=1.
105 for (Index i = 0; i < nf; i++) {
106 for (Index s = 0; s < 2; s++) {
107 if (emissivity(i, s) > 1) {
108 emissivity(i, s) = 1;
109 reflectivity(i, s) = 0;
110 }
111 if (emissivity(i, s) < 0) {
112 emissivity(i, s) = 0;
113 reflectivity(i, s) = 1;
114 }
115 if (reflectivity(i, s) > 1) {
116 emissivity(i, s) = 0;
117 reflectivity(i, s) = 1;
118 }
119 if (reflectivity(i, s) < 0) {
120 emissivity(i, s) = 1;
121 reflectivity(i, s) = 0;
122 }
123 }
124 }
125}
126
127/* Workspace method: Doxygen documentation will be auto-generated */
128void InterpGriddedField2ToPosition(Numeric& outvalue,
129 const Index& atmosphere_dim,
130 const Vector& lat_grid,
131 const Vector& lat_true,
132 const Vector& lon_true,
133 const Vector& rtp_pos,
134 const GriddedField2& gfield2,
135 const Verbosity&) {
136 // Set expected order of grids
137 Index gfield_latID = 0;
138 Index gfield_lonID = 1;
139
140 // Basic checks and sizes
141 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
142 chk_latlon_true(atmosphere_dim, lat_grid, lat_true, lon_true);
143 chk_rte_pos(atmosphere_dim, rtp_pos);
144 gfield2.checksize_strict();
145 //
146 chk_griddedfield_gridname(gfield2, gfield_latID, "Latitude");
147 chk_griddedfield_gridname(gfield2, gfield_lonID, "Longitude");
148 //
149 const Index nlat = gfield2.data.nrows();
150 const Index nlon = gfield2.data.ncols();
151 //
152 if (nlat < 2 || nlon < 2) {
153 ostringstream os;
154 os << "The data in *gfield2* must span a geographical region. That is,\n"
155 << "the latitude and longitude grids must have a length >= 2.";
156 throw runtime_error(os.str());
157 }
158
159 const Vector& GFlat = gfield2.get_numeric_grid(gfield_latID);
160 const Vector& GFlon = gfield2.get_numeric_grid(gfield_lonID);
161
162 // Determine true geographical position
163 Vector lat(1), lon(1);
165 lat[0], lon[0], atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
166
167 // Ensure correct coverage of lon grid
168 Vector lon_shifted;
169 lon_shiftgrid(lon_shifted, GFlon, lon[0]);
170
171 // Check if lat/lon we need are actually covered
172 chk_if_in_range("rtp_pos.lat", lat[0], GFlat[0], GFlat[nlat - 1]);
173 chk_if_in_range("rtp_pos.lon", lon[0], lon_shifted[0], lon_shifted[nlon - 1]);
174
175 // Interpolate in lat and lon
176 //
177 GridPos gp_lat, gp_lon;
178 gridpos(gp_lat, GFlat, lat[0]);
179 gridpos(gp_lon, lon_shifted, lon[0]);
180 Vector itw(4);
181 interpweights(itw, gp_lat, gp_lon);
182 outvalue = interp(itw, gfield2.data, gp_lat, gp_lon);
183}
184
185
186/* Workspace method: Doxygen documentation will be auto-generated */
187void InterpSurfaceFieldToPosition(Numeric& outvalue,
188 const Index& atmosphere_dim,
189 const Vector& lat_grid,
190 const Vector& lon_grid,
191 const Vector& rtp_pos,
192 const Matrix& z_surface,
193 const Matrix& field,
194 const Verbosity& verbosity) {
195 // Input checks (dummy p_grid)
196 chk_atm_grids(atmosphere_dim, uniform_grid(2, 2, -1), lat_grid, lon_grid);
198 "input argument *field*", field, atmosphere_dim, lat_grid, lon_grid);
199 chk_rte_pos(atmosphere_dim, rtp_pos);
200 //
201 const Numeric zmax = max(z_surface);
202 const Numeric zmin = min(z_surface);
203 const Numeric dzok = 1;
204 if (rtp_pos[0] < zmin - dzok || rtp_pos[0] > zmax + dzok) {
205 ostringstream os;
206 os << "The given position does not match *z_surface*.\nThe altitude in "
207 << "*rtp_pos* is " << rtp_pos[0] / 1e3 << " km.\n"
208 << "The altitude range covered by *z_surface* is [" << zmin / 1e3 << ","
209 << zmax / 1e3 << "] km.\n"
210 << "One possible mistake is to mix up *rtp_pos* and *rte_los*.";
211 throw runtime_error(os.str());
212 }
213
214 if (atmosphere_dim == 1) {
215 outvalue = field(0, 0);
216 } else {
217 chk_interpolation_grids("Latitude interpolation", lat_grid, rtp_pos[1]);
218 GridPos gp_lat, gp_lon;
219 gridpos(gp_lat, lat_grid, rtp_pos[1]);
220 if (atmosphere_dim == 3) {
221 chk_interpolation_grids("Longitude interpolation", lon_grid, rtp_pos[2]);
222 gridpos(gp_lon, lon_grid, rtp_pos[2]);
223 }
224 //
225 outvalue = interp_atmsurface_by_gp(atmosphere_dim, field, gp_lat, gp_lon);
226 }
227
228 // Interpolate
230 out3 << " Result = " << outvalue << "\n";
231}
232
233/* Workspace method: Doxygen documentation will be auto-generated */
235 Matrix& iy,
236 ArrayOfTensor3& diy_dx,
237 const Tensor3& iy_transmittance,
238 const Index& iy_id,
239 const Index& jacobian_do,
240 const Index& atmosphere_dim,
241 const EnergyLevelMap& nlte_field,
242 const Index& cloudbox_on,
243 const Index& stokes_dim,
244 const Vector& f_grid,
245 const Vector& rtp_pos,
246 const Vector& rtp_los,
247 const Vector& rte_pos2,
248 const String& iy_unit,
249 const Agenda& iy_main_agenda,
250 const Numeric& surface_skin_t,
251 const Numeric& salinity,
252 const Numeric& wind_speed,
253 const Numeric& wind_direction,
254 const Index& fastem_version,
255 const Verbosity& verbosity) {
256 // Most obvious input checks are performed in specular_losCalc and surfaceFastem
257
258 // Obtain radiance and transmission for specular direction
259
260 // Determine specular direction
261 Vector specular_los, surface_normal;
262 specular_losCalcNoTopography(specular_los,
263 surface_normal,
264 rtp_pos,
265 rtp_los,
266 atmosphere_dim,
267 verbosity);
268
269 // Use iy_aux to get optical depth for downwelling radiation.
270 ArrayOfString iy_aux_vars(1);
271 iy_aux_vars[0] = "Optical depth";
272
273 // Calculate iy for downwelling radiation
274 // Note that iy_transmittance used here lacks surface R. Fixed below.
275 //
276 const Index nf = f_grid.nelem();
277 Vector transmittance(nf);
278 Vector geo_pos;
279 ArrayOfMatrix iy_aux;
280 Ppath ppath;
281 //
283 iy,
284 iy_aux,
285 ppath,
286 diy_dx,
287 geo_pos,
288 0,
289 iy_transmittance,
290 iy_aux_vars,
291 iy_id,
292 iy_unit,
293 cloudbox_on,
294 jacobian_do,
295 f_grid,
296 nlte_field,
297 rtp_pos,
298 specular_los,
299 rte_pos2,
300 iy_main_agenda);
301
302 // Convert tau to transmissions
303 for (Index i = 0; i < nf; i++) {
304 transmittance[i] = exp(-iy_aux[0](i, 0));
305 }
306
307 // Call Fastem by surface_RTprop version
308 //
309 Matrix surface_los;
310 Tensor4 surface_rmatrix;
311 Matrix surface_emission;
312 //
313 surfaceFastem(surface_los,
314 surface_rmatrix,
315 surface_emission,
316 atmosphere_dim,
317 stokes_dim,
318 f_grid,
319 rtp_pos,
320 rtp_los,
321 surface_skin_t,
322 salinity,
323 wind_speed,
324 wind_direction,
325 transmittance,
326 fastem_version,
327 verbosity);
328
329 // Add up
330 //
331 Tensor3 I(1, nf, stokes_dim);
332 I(0, joker, joker) = iy;
333 Matrix sensor_los_dummy(1, 1, 0);
334 //
335 surface_calc(iy, I, sensor_los_dummy, surface_rmatrix, surface_emission);
336
337 // Adjust diy_dx, if necessary.
338 // For vector cases this is a slight approximation, as the order of the
339 // different transmission and reflectivities matters.
340 if (iy_transmittance.npages()) {
341 for (Index q = 0; q < diy_dx.nelem(); q++) {
342 for (Index p = 0; p < diy_dx[q].npages(); p++) {
343 for (Index i = 0; i < nf; i++) {
344 Vector x{diy_dx[q](p, i, joker)};
345 mult(diy_dx[q](p, i, joker), surface_rmatrix(0, i, joker, joker), x);
346 }
347 }
348 }
349 }
350}
351
352/* Workspace method: Doxygen documentation will be auto-generated */
354 Matrix& iy,
355 ArrayOfTensor3& diy_dx,
356 ArrayOfTensor4& dsurface_rmatrix_dx,
357 ArrayOfMatrix& dsurface_emission_dx,
358 const Tensor3& iy_transmittance,
359 const Index& iy_id,
360 const Index& jacobian_do,
361 const Index& suns_do,
362 const Index& atmosphere_dim,
363 const EnergyLevelMap& nlte_field,
364 const Index& cloudbox_on,
365 const Index& stokes_dim,
366 const Vector& f_grid,
367 const Vector& lat_grid,
368 const Vector& lon_grid,
369 const Matrix& z_surface,
370 const Vector& refellipsoid,
371 const Vector& rtp_pos,
372 const Vector& rtp_los,
373 const Vector& rte_pos2,
374 const String& iy_unit,
375 const Tensor3& surface_reflectivity,
376 const Tensor3& surface_props_data,
377 const ArrayOfString& surface_props_names,
378 const ArrayOfString& dsurface_names,
379 const ArrayOfRetrievalQuantity& jacobian_quantities,
380 const Agenda& iy_main_agenda,
381 const Verbosity& verbosity) {
382
383 // Input checks
384 ARTS_USER_ERROR_IF(atmosphere_dim==2, "This method does not work for 2d atmospheres.")
385 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
386 chk_rte_pos(atmosphere_dim, rtp_pos);
387 chk_rte_los(atmosphere_dim, rtp_los);
388 chk_size("iy",iy,f_grid.nelem(),stokes_dim);
389
390 // Check surface_data
391 surface_props_check(atmosphere_dim,
392 lat_grid,
393 lon_grid,
394 surface_props_data,
395 surface_props_names);
396
397 // Interplation grid positions and weights
398 ArrayOfGridPos gp_lat(1), gp_lon(1);
399 Matrix itw;
401 gp_lat[0], gp_lon[0], atmosphere_dim, lat_grid, lon_grid, rtp_pos);
402 interp_atmsurface_gp2itw(itw, atmosphere_dim, gp_lat, gp_lon);
403
404 // Skin temperature
405 Vector surface_skin_t(1);
406 surface_props_interp(surface_skin_t,
407 "Skin temperature",
408 atmosphere_dim,
409 gp_lat,
410 gp_lon,
411 itw,
412 surface_props_data,
413 surface_props_names);
414
415 chk_not_negative("surface_skin_t", surface_skin_t[0]);
416
417
418
419 //Allocate
420 Tensor3 iy_trans_new;
421 Matrix iy_incoming;
422 Vector specular_los, surface_normal;
423 ArrayOfMatrix iy_aux;
424 Ppath ppath;
425 ArrayOfTensor3 diy_dx_dumb;
426 Matrix surface_los;
427 Tensor4 surface_rmatrix;
428 Matrix surface_emission;
429 Tensor3 I(1,f_grid.nelem(),stokes_dim);
430 if (iy_transmittance.npages()) {
431 iy_trans_new=iy_transmittance;
432 }
433
434 //get specular line of sight
435 specular_losCalc(specular_los,
436 surface_normal,
437 rtp_pos,
438 rtp_los,
439 atmosphere_dim,
440 lat_grid,
441 lon_grid,
442 refellipsoid,
443 z_surface,
444 0,
445 verbosity);
446
447 // Calculate incoming radiation directions, surface reflection matrix and
448 // emission vector
449 surfaceFlatReflectivity(surface_los,
450 surface_rmatrix,
451 surface_emission,
452 f_grid,
453 stokes_dim,
454 atmosphere_dim,
455 rtp_pos,
456 rtp_los,
457 specular_los,
458 surface_skin_t[0],
459 surface_reflectivity,
460 verbosity);
461
462
463
464 // Jacobian part
465 if (jacobian_do) {
466 dsurface_check(surface_props_names,
467 dsurface_names,
468 dsurface_rmatrix_dx,
469 dsurface_emission_dx);
470
471 Index irq;
472
473 // Skin temperature
474 irq = find_first(dsurface_names, String("Skin temperature"));
475 if (irq >= 0) {
476 const Numeric dd = 0.1;
477 Matrix surface_los2;
478 surfaceFlatReflectivity(surface_los2,
479 dsurface_rmatrix_dx[irq],
480 dsurface_emission_dx[irq],
481 f_grid,
482 stokes_dim,
483 atmosphere_dim,
484 rtp_pos,
485 rtp_los,
486 specular_los,
487 surface_skin_t[0]+dd,
488 surface_reflectivity,
489 verbosity);
490 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
491 dsurface_rmatrix_dx[irq] /= dd;
492 dsurface_emission_dx[irq] -= surface_emission;
493 dsurface_emission_dx[irq] /= dd;
494 }
495 }
496
498 iy,
499 diy_dx,
500 surface_los,
501 surface_rmatrix,
502 surface_emission,
503 dsurface_names,
504 dsurface_rmatrix_dx,
505 dsurface_emission_dx,
506 iy_transmittance,
507 iy_id,
508 jacobian_do,
509 suns_do,
510 jacobian_quantities,
511 atmosphere_dim,
512 nlte_field,
513 cloudbox_on,
514 stokes_dim,
515 f_grid,
516 rtp_pos,
517 rtp_los,
518 rte_pos2,
519 iy_unit,
520 iy_main_agenda,
521 verbosity);
522
523
524}
525
526/* Workspace method: Doxygen documentation will be auto-generated */
528 Workspace& ws,
529 Matrix& iy,
530 const Vector& rtp_pos,
531 const Vector& rtp_los,
532 const Index& stokes_dim,
533 const Vector& f_grid,
534 const Index& atmosphere_dim,
535 const Vector& p_grid,
536 const Vector& lat_grid,
537 const Vector& lon_grid,
538 const Tensor3& z_field,
539 const Tensor3& t_field,
540 const EnergyLevelMap& nlte_field,
541 const Tensor4& vmr_field,
542 const ArrayOfArrayOfSpeciesTag& abs_species,
543 const Tensor3& wind_u_field,
544 const Tensor3& wind_v_field,
545 const Tensor3& wind_w_field,
546 const Tensor3& mag_u_field,
547 const Tensor3& mag_v_field,
548 const Tensor3& mag_w_field,
549 const Matrix& z_surface,
550 const Tensor3& surface_reflectivity,
551 const Vector& refellipsoid,
552 const Tensor4& pnd_field,
553 const ArrayOfTensor4& dpnd_field_dx,
554 const ArrayOfString& scat_species,
555 const ArrayOfArrayOfSingleScatteringData& scat_data,
556 const Numeric& ppath_lmax,
557 const Numeric& ppath_lraytrace,
558 const Index& ppath_inside_cloudbox_do,
559 const Index& cloudbox_on,
560 const ArrayOfIndex& cloudbox_limits,
561 const Index& suns_do,
562 const Index& gas_scattering_do,
563 const Index& jacobian_do,
564 const ArrayOfRetrievalQuantity& jacobian_quantities,
565 const ArrayOfSun& suns,
566 const Numeric& rte_alonglos_v,
567 const String& iy_unit,
568 const Agenda& propmat_clearsky_agenda,
569 const Agenda& water_p_eq_agenda,
570 const Agenda& gas_scattering_agenda,
571 const Agenda& ppath_step_agenda,
572 const Verbosity& verbosity) {
573 //Check for correct unit
574 ARTS_USER_ERROR_IF(iy_unit != "1" && suns_do,
575 "If suns are present only iy_unit=\"1\" can be used.");
576
577 chk_size("iy", iy, f_grid.nelem(), stokes_dim);
578
579 if (suns_do) {
580 Matrix iy_incoming;
581 Index suns_visible;
582 Vector specular_los;
583
584 //Get incoming direct radiation, if no sun is in line of sight, then
585 //iy_incoming is zero.
587 iy_incoming,
588 suns_visible,
589 specular_los,
590 rtp_pos,
591 rtp_los,
592 stokes_dim,
593 f_grid,
594 atmosphere_dim,
595 p_grid,
596 lat_grid,
597 lon_grid,
598 z_field,
599 t_field,
600 nlte_field,
601 vmr_field,
602 abs_species,
603 wind_u_field,
604 wind_v_field,
605 wind_w_field,
606 mag_u_field,
607 mag_v_field,
608 mag_w_field,
609 z_surface,
610 refellipsoid,
611 pnd_field,
612 dpnd_field_dx,
613 scat_species,
614 scat_data,
615 ppath_lmax,
616 ppath_lraytrace,
617 ppath_inside_cloudbox_do,
618 cloudbox_on,
619 cloudbox_limits,
620 gas_scattering_do,
621 jacobian_do,
622 jacobian_quantities,
623 suns,
624 rte_alonglos_v,
625 propmat_clearsky_agenda,
626 water_p_eq_agenda,
627 gas_scattering_agenda,
628 ppath_step_agenda,
629 verbosity);
630
631 if (suns_visible) {
632 Matrix surface_los;
633 Tensor4 surface_rmatrix;
634 Matrix surface_emission;
635 Numeric surface_skin_t_dummy = 1.;
636
637 //As we only consider the reflection of the direct radiation we do not consider
638 //any type of surface emission, therefore we set the surface skin temperature
639 //to a dummy value and later on set the surface emission vector (matrix) to zero.
640 surfaceFlatReflectivity(surface_los,
641 surface_rmatrix,
642 surface_emission,
643 f_grid,
644 stokes_dim,
645 atmosphere_dim,
646 rtp_pos,
647 rtp_los,
648 specular_los,
649 surface_skin_t_dummy,
650 surface_reflectivity,
651 verbosity);
652
653 surface_emission *= 0.;
654
655 Tensor3 I(1, f_grid.nelem(), stokes_dim);
656 I(0, joker, joker) = iy_incoming;
657
658 surface_calc(iy, I, surface_los, surface_rmatrix, surface_emission);
659
660 //TODO: add reflectivity jacobian
661 }
662 }
663}
664
665/* Workspace method: Doxygen documentation will be auto-generated */
667 Matrix& iy,
668 ArrayOfTensor3& diy_dx,
669 ArrayOfTensor4& dsurface_rmatrix_dx,
670 ArrayOfMatrix& dsurface_emission_dx,
671 const Tensor3& iy_transmittance,
672 const Index& iy_id,
673 const Index& jacobian_do,
674 const Index& suns_do,
675 const Index& atmosphere_dim,
676 const EnergyLevelMap& nlte_field,
677 const Index& cloudbox_on,
678 const Index& stokes_dim,
679 const Vector& f_grid,
680 const Vector& lat_grid,
681 const Vector& lon_grid,
682 const Matrix& z_surface,
683 const Vector& refellipsoid,
684 const Vector& rtp_pos,
685 const Vector& rtp_los,
686 const Vector& rte_pos2,
687 const String& iy_unit,
688 const GriddedField3& surface_complex_refr_index,
689 const Tensor3& surface_props_data,
690 const ArrayOfString& surface_props_names,
691 const ArrayOfString& dsurface_names,
692 const ArrayOfRetrievalQuantity& jacobian_quantities,
693 const Agenda& iy_main_agenda,
694 const Verbosity& verbosity) {
695
696 // Input checks
697 ARTS_USER_ERROR_IF(atmosphere_dim==2, "This method does not work for 2d atmospheres.")
698 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
699 chk_rte_pos(atmosphere_dim, rtp_pos);
700 chk_rte_los(atmosphere_dim, rtp_los);
701 chk_size("iy",iy,f_grid.nelem(),stokes_dim);
702
703 // Check surface_data
704 surface_props_check(atmosphere_dim,
705 lat_grid,
706 lon_grid,
707 surface_props_data,
708 surface_props_names);
709
710 // Interplation grid positions and weights
711 ArrayOfGridPos gp_lat(1), gp_lon(1);
712 Matrix itw;
714 gp_lat[0], gp_lon[0], atmosphere_dim, lat_grid, lon_grid, rtp_pos);
715 interp_atmsurface_gp2itw(itw, atmosphere_dim, gp_lat, gp_lon);
716
717 // Skin temperature
718 Vector surface_skin_t(1);
719 surface_props_interp(surface_skin_t,
720 "Skin temperature",
721 atmosphere_dim,
722 gp_lat,
723 gp_lon,
724 itw,
725 surface_props_data,
726 surface_props_names);
727
728 chk_not_negative("surface_skin_t", surface_skin_t[0]);
729
730
731 //Allocate
732 Tensor3 iy_trans_new;
733 Matrix iy_incoming;
734 Vector specular_los, surface_normal;
735 ArrayOfMatrix iy_aux;
736 Ppath ppath;
737 ArrayOfTensor3 diy_dx_dumb;
738 Matrix surface_los;
739 Tensor4 surface_rmatrix;
740 Matrix surface_emission;
741 Tensor3 I(1,f_grid.nelem(),stokes_dim);
742 if (iy_transmittance.npages()) {
743 iy_trans_new=iy_transmittance;
744 }
745
746 //get specular line of sight
747 specular_losCalc(specular_los,
748 surface_normal,
749 rtp_pos,
750 rtp_los,
751 atmosphere_dim,
752 lat_grid,
753 lon_grid,
754 refellipsoid,
755 z_surface,
756 0,
757 verbosity);
758
759 // Calculate incoming radiation directions, surface reflection matrix and
760 // emission vector
761 surfaceFlatRefractiveIndex(surface_los,
762 surface_rmatrix,
763 surface_emission,
764 f_grid,
765 stokes_dim,
766 atmosphere_dim,
767 rtp_pos,
768 rtp_los,
769 specular_los,
770 surface_skin_t[0],
771 surface_complex_refr_index,
772 verbosity);
773
774
775
776 // Jacobian part
777 if (jacobian_do) {
778 dsurface_check(surface_props_names,
779 dsurface_names,
780 dsurface_rmatrix_dx,
781 dsurface_emission_dx);
782
783 Index irq;
784
785 // Skin temperature
786 irq = find_first(dsurface_names, String("Skin temperature"));
787 if (irq >= 0) {
788 const Numeric dd = 0.1;
789 Matrix surface_los2;
790 surfaceFlatRefractiveIndex(surface_los2,
791 dsurface_rmatrix_dx[irq],
792 dsurface_emission_dx[irq],
793 f_grid,
794 stokes_dim,
795 atmosphere_dim,
796 rtp_pos,
797 rtp_los,
798 specular_los,
799 surface_skin_t[0]+dd,
800 surface_complex_refr_index,
801 verbosity);
802 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
803 dsurface_rmatrix_dx[irq] /= dd;
804 dsurface_emission_dx[irq] -= surface_emission;
805 dsurface_emission_dx[irq] /= dd;
806 }
807 }
808
810 iy,
811 diy_dx,
812 surface_los,
813 surface_rmatrix,
814 surface_emission,
815 dsurface_names,
816 dsurface_rmatrix_dx,
817 dsurface_emission_dx,
818 iy_transmittance,
819 iy_id,
820 jacobian_do,
821 suns_do,
822 jacobian_quantities,
823 atmosphere_dim,
824 nlte_field,
825 cloudbox_on,
826 stokes_dim,
827 f_grid,
828 rtp_pos,
829 rtp_los,
830 rte_pos2,
831 iy_unit,
832 iy_main_agenda,
833 verbosity);
834
835
836}
837
838/* Workspace method: Doxygen documentation will be auto-generated */
840 Workspace& ws,
841 Matrix& iy,
842 const Vector& rtp_pos,
843 const Vector& rtp_los,
844 const Index& stokes_dim,
845 const Vector& f_grid,
846 const Index& atmosphere_dim,
847 const Vector& p_grid,
848 const Vector& lat_grid,
849 const Vector& lon_grid,
850 const Tensor3& z_field,
851 const Tensor3& t_field,
852 const EnergyLevelMap& nlte_field,
853 const Tensor4& vmr_field,
854 const ArrayOfArrayOfSpeciesTag& abs_species,
855 const Tensor3& wind_u_field,
856 const Tensor3& wind_v_field,
857 const Tensor3& wind_w_field,
858 const Tensor3& mag_u_field,
859 const Tensor3& mag_v_field,
860 const Tensor3& mag_w_field,
861 const Matrix& z_surface,
862 const GriddedField3& surface_complex_refr_index,
863 const Vector& refellipsoid,
864 const Tensor4& pnd_field,
865 const ArrayOfTensor4& dpnd_field_dx,
866 const ArrayOfString& scat_species,
867 const ArrayOfArrayOfSingleScatteringData& scat_data,
868 const Numeric& ppath_lmax,
869 const Numeric& ppath_lraytrace,
870 const Index& ppath_inside_cloudbox_do,
871 const Index& cloudbox_on,
872 const ArrayOfIndex& cloudbox_limits,
873 const Index& suns_do,
874 const Index& gas_scattering_do,
875 const Index& jacobian_do,
876 const ArrayOfRetrievalQuantity& jacobian_quantities,
877 const ArrayOfSun& suns,
878 const Numeric& rte_alonglos_v,
879 const String& iy_unit,
880 const Agenda& propmat_clearsky_agenda,
881 const Agenda& water_p_eq_agenda,
882 const Agenda& gas_scattering_agenda,
883 const Agenda& ppath_step_agenda,
884 const Verbosity& verbosity) {
885 //Check for correct unit
886 ARTS_USER_ERROR_IF(iy_unit != "1" && suns_do,
887 "If suns are present only iy_unit=\"1\" can be used.");
888
889 chk_size("iy", iy, f_grid.nelem(), stokes_dim);
890
891 if (suns_do) {
892 Matrix iy_incoming;
893 Index suns_visible;
894 Vector specular_los;
895
896 //Get incoming direct radiation, if no sun is in line of sight, then
897 //iy_incoming is zero.
899 iy_incoming,
900 suns_visible,
901 specular_los,
902 rtp_pos,
903 rtp_los,
904 stokes_dim,
905 f_grid,
906 atmosphere_dim,
907 p_grid,
908 lat_grid,
909 lon_grid,
910 z_field,
911 t_field,
912 nlte_field,
913 vmr_field,
914 abs_species,
915 wind_u_field,
916 wind_v_field,
917 wind_w_field,
918 mag_u_field,
919 mag_v_field,
920 mag_w_field,
921 z_surface,
922 refellipsoid,
923 pnd_field,
924 dpnd_field_dx,
925 scat_species,
926 scat_data,
927 ppath_lmax,
928 ppath_lraytrace,
929 ppath_inside_cloudbox_do,
930 cloudbox_on,
931 cloudbox_limits,
932 gas_scattering_do,
933 jacobian_do,
934 jacobian_quantities,
935 suns,
936 rte_alonglos_v,
937 propmat_clearsky_agenda,
938 water_p_eq_agenda,
939 gas_scattering_agenda,
940 ppath_step_agenda,
941 verbosity);
942
943 if (suns_visible) {
944 Matrix surface_los;
945 Tensor4 surface_rmatrix;
946 Matrix surface_emission;
947 Numeric surface_skin_t_dummy = 1.;
948
949 //As we only consider the reflection of the direct radiation we do not consider
950 //any type of surface emission, therefore we set the surface skin temperature
951 //to a dummy value and later on set the surface emission vector (matrix) to zero.
952 surfaceFlatRefractiveIndex(surface_los,
953 surface_rmatrix,
954 surface_emission,
955 f_grid,
956 stokes_dim,
957 atmosphere_dim,
958 rtp_pos,
959 rtp_los,
960 specular_los,
961 surface_skin_t_dummy,
962 surface_complex_refr_index,
963 verbosity);
964
965 surface_emission *= 0.;
966
967 Tensor3 I(1, f_grid.nelem(), stokes_dim);
968 I(0, joker, joker) = iy_incoming;
969
970 surface_calc(iy, I, surface_los, surface_rmatrix, surface_emission);
971
972 //TODO: add reflectivity jacobian
973 }
974 }
975}
976
977/* Workspace method: Doxygen documentation will be auto-generated */
978void iySurfaceInit(Matrix& iy,
979 const Vector& f_grid,
980 const Index& stokes_dim,
981 const Verbosity&) {
982 iy.resize(f_grid.nelem(), stokes_dim);
983 iy = 0.;
984}
985
986/* Workspace method: Doxygen documentation will be auto-generated */
988 Matrix& iy,
989 ArrayOfTensor3& diy_dx,
990 const Tensor3& iy_transmittance,
991 const Index& iy_id,
992 const Index& jacobian_do,
993 const Index& suns_do,
994 const Index& atmosphere_dim,
995 const EnergyLevelMap& nlte_field,
996 const Index& cloudbox_on,
997 const Index& stokes_dim,
998 const Vector& f_grid,
999 const Vector& lat_grid,
1000 const Vector& lon_grid,
1001 const Matrix& z_surface,
1002 const Vector& refellipsoid,
1003 const Vector& rtp_pos,
1004 const Vector& rtp_los,
1005 const Vector& rte_pos2,
1006 const String& iy_unit,
1007 const Vector& surface_scalar_reflectivity,
1008 const Tensor3& surface_props_data,
1009 const ArrayOfString& surface_props_names,
1010 const ArrayOfString& dsurface_names,
1011 const ArrayOfRetrievalQuantity& jacobian_quantities,
1012 const Agenda& iy_main_agenda,
1013 const Index& N_za,
1014 const Index& N_aa,
1015 const Verbosity& verbosity) {
1016 // Input checks
1017 ARTS_USER_ERROR_IF(surface_scalar_reflectivity.nelem() != f_grid.nelem() &&
1018 surface_scalar_reflectivity.nelem() != 1,
1019 "The number of elements in *surface_scalar_reflectivity* should\n",
1020 "match length of *f_grid* or be 1.",
1021 "\n length of *f_grid* : ", f_grid.nelem(),
1022 "\n length of *surface_scalar_reflectivity* : ",
1023 surface_scalar_reflectivity.nelem());
1024 ARTS_USER_ERROR_IF(min(surface_scalar_reflectivity) < 0 ||
1025 max(surface_scalar_reflectivity) > 1,
1026 "All values in *surface_scalar_reflectivity* must be inside [0,1].");
1027 ARTS_USER_ERROR_IF(atmosphere_dim==2, "This method does not work for 2d atmospheres.");
1028 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
1029 chk_rte_pos(atmosphere_dim, rtp_pos);
1030 chk_rte_los(atmosphere_dim, rtp_los);
1031 chk_size("iy",iy,f_grid.nelem(),stokes_dim);
1032
1033 // Check surface_data
1034 surface_props_check(atmosphere_dim,
1035 lat_grid,
1036 lon_grid,
1037 surface_props_data,
1038 surface_props_names);
1039
1040 // Interplation grid positions and weights
1041 ArrayOfGridPos gp_lat(1), gp_lon(1);
1042 Matrix itw;
1044 gp_lat[0], gp_lon[0], atmosphere_dim, lat_grid, lon_grid, rtp_pos);
1045 interp_atmsurface_gp2itw(itw, atmosphere_dim, gp_lat, gp_lon);
1046
1047 // Skin temperature
1048 Vector surface_skin_t(1);
1049 surface_props_interp(surface_skin_t,
1050 "Skin temperature",
1051 atmosphere_dim,
1052 gp_lat,
1053 gp_lon,
1054 itw,
1055 surface_props_data,
1056 surface_props_names);
1057
1058 chk_not_negative("surface_skin_t", surface_skin_t[0]);
1059
1060
1061
1062
1063 const Index nf = f_grid.nelem();
1064
1065 Vector za_grid, aa_grid, za_grid_weights;
1066
1067 //get angular grids. here we use gauss quadrature. As the function is not intended for
1068 //a hemisphere but for full sphere, we have to tweak it a little bit and use only
1069 //the first half of the angles.
1071 aa_grid,
1072 za_grid_weights,
1073 N_za * 2,
1074 N_aa,
1075 "double_gauss",
1076 verbosity);
1077
1078 Tensor3 iy_trans_new;
1079
1080
1081 ArrayOfString iy_aux_var(0);
1082 if (suns_do) {
1083 iy_aux_var.resize(1);
1084 iy_aux_var[0] = "Direct radiation";
1085 }
1086
1087 Matrix iy_temp;
1088 Matrix Flx(nf, 1, 0.);
1089
1090 if (iy_transmittance.npages()) {
1091 iy_trans_new=iy_transmittance;
1092 }
1093
1094
1095 ArrayOfTensor3 diy_dx_dumb;
1096 iy.resize(nf, stokes_dim);
1097 iy = 0;
1098
1099 Vector los;
1100
1101 //No azimuth dependency
1102 if (N_aa == 1) {
1103 //Set azimuth angle index
1104 Index i_aa = 0;
1105
1106 for (Index i_za = 0; i_za < N_za; i_za++) {
1107
1108 if (atmosphere_dim==1){
1109 los.resize(1);
1110 }
1111 else{
1112 los.resize(2);
1113 // For clearsky, the skylight is in general approximately independent of the azimuth direction.
1114 // So, we consider only the line of sight plane.
1115 los[1] = 0.;
1116 }
1117 los[0] = za_grid[i_za];
1118
1119
1120 // Calculate downwelling radiation for a los
1121 ArrayOfMatrix iy_aux;
1122 Ppath ppath;
1123 Vector geo_pos_dumb;
1124 Index iy_id_new = iy_id + i_za + i_aa * N_za + 1;
1126 iy_temp,
1127 iy_aux,
1128 ppath,
1129 diy_dx_dumb,
1130 geo_pos_dumb,
1131 0,
1132 iy_trans_new,
1133 iy_aux_var,
1134 iy_id_new,
1135 iy_unit,
1136 cloudbox_on,
1137 jacobian_do,
1138 f_grid,
1139 nlte_field,
1140 rtp_pos,
1141 los,
1142 rte_pos2,
1143 iy_main_agenda);
1144
1145 //For the case that a sun is present and the los is towards a sun, we
1146 //subtract the direct radiation, so that only the diffuse radiation is considered here.
1147 //If sun is within los iy_aux[0] is the incoming and attenuated direct (sun)
1148 //radiation at the surface. Otherwise, it is zero.
1149 if (suns_do) {
1150 iy_temp -= iy_aux[0];
1151 }
1152 iy_temp *= abs(cos(deg2rad(los[0])) * sin(deg2rad(los[0]))) *
1153 za_grid_weights[i_za];
1154 Flx(joker, 0) += iy_temp(joker, 0);
1155 }
1156 //Azimuthal integration just gives a factor of 2*pi
1157 Flx *= 2 * pi;
1158
1159 //Azimuth dependency
1160 } else {
1161 // before we start with the actual RT calculation, we calculate the integration
1162 // weights for a trapezoidal integration over the azimuth.
1163 Vector aa_grid_weights;
1164 Vector aa_grid_rad = aa_grid;
1165
1166 //For the integration we need the aa grid in rad
1167 for (Index i = 0; i < aa_grid_rad.nelem(); i++) {
1168 aa_grid_rad[i] = deg2rad(aa_grid_rad[i]);
1169 }
1170
1171 calculate_int_weights_arbitrary_grid(aa_grid_weights, aa_grid_rad);
1172
1173 //Allocate
1174 los.resize(2);
1175
1176 for (Index i_aa = 0; i_aa < N_aa; i_aa++) {
1177 for (Index i_za = 0; i_za < N_za; i_za++) {
1178
1179
1180 los[0] = za_grid[i_za];
1181 los[1] = aa_grid[i_aa];
1182
1183 // Calculate downwelling radiation for LOS ilos
1184 ArrayOfMatrix iy_aux;
1185 Ppath ppath;
1186 Vector geo_pos_dumb;
1187 Index iy_id_new = iy_id + i_za + i_aa * N_za + 1;
1189 iy_temp,
1190 iy_aux,
1191 ppath,
1192 diy_dx_dumb,
1193 geo_pos_dumb,
1194 0,
1195 iy_trans_new,
1196 iy_aux_var,
1197 iy_id_new,
1198 iy_unit,
1199 cloudbox_on,
1200 jacobian_do,
1201 f_grid,
1202 nlte_field,
1203 rtp_pos,
1204 los,
1205 rte_pos2,
1206 iy_main_agenda);
1207
1208 //For the case that a sun is present and the los is towards a sun, we
1209 //subtract the direct radiation, so that only the diffuse radiation is considered here.
1210 //If sun is within los iy_aux[0] is the incoming and attenuated direct (sun)
1211 //radiation at the surface. Otherwise it is zero.
1212 if (suns_do) {
1213 iy_temp -= iy_aux[0];
1214 }
1215 iy_temp *= abs(cos(deg2rad(los[0])) * sin(deg2rad(los[0]))) *
1216 za_grid_weights[i_za] * aa_grid_weights[i_aa];
1217 Flx(joker, 0) += iy_temp(joker, 0);
1218 }
1219 }
1220 }
1221
1222 Vector specular_los, surface_normal;
1223 specular_losCalc(specular_los,
1224 surface_normal,
1225 rtp_pos,
1226 rtp_los,
1227 atmosphere_dim,
1228 lat_grid,
1229 lon_grid,
1230 refellipsoid,
1231 z_surface,
1232 0,
1233 verbosity);
1234
1235 //Surface emission
1236 Vector b(nf);
1237 planck(b, f_grid, surface_skin_t[0]);
1238
1239 //Upwelling radiation
1240 Numeric r;
1241 for (Index i_freq = 0; i_freq < f_grid.nelem(); i_freq++) {
1242 if (i_freq == 0 || surface_scalar_reflectivity.nelem() > 1) {
1243 r = surface_scalar_reflectivity[i_freq];
1244 }
1245 iy(i_freq, 0) = r *
1246 cos(deg2rad(specular_los[0])) / pi * Flx(i_freq, 0) +
1247 (1 - r) * b[i_freq];
1248 }
1249
1250 // Surface Jacobians
1251 if (jacobian_do && dsurface_names.nelem()) {
1252
1253 Index irq;
1254 irq = find_first(dsurface_names, String("Skin temperature"));
1255
1256 if (irq >= 0){
1257 // Find index among jacobian quantities
1258 Index ihit = -1;
1259 for (Index j = 0; j < jacobian_quantities.nelem() && ihit < 0; j++) {
1260 if (dsurface_names[irq] == jacobian_quantities[j].Subtag()) {
1261 ihit = j;
1262 }
1263 }
1264
1265 // Derivative of surface skin temperature, as observed at the surface
1266 Matrix diy_dpos0(f_grid.nelem(), stokes_dim, 0.), diy_dpos;
1267 dplanck_dt(diy_dpos0(joker, 0), f_grid, surface_skin_t[0]);
1268
1269 Vector emissivity=surface_scalar_reflectivity;
1270 emissivity*=-1;
1271 emissivity+=1;
1272 diy_dpos0*=ExhaustiveMatrixView{emissivity};
1273
1274 // Weight with transmission to sensor
1275 iy_transmittance_mult(diy_dpos, iy_transmittance, diy_dpos0);
1276
1277 // Put into diy_dx
1278 diy_from_pos_to_rgrids(diy_dx[ihit],
1279 jacobian_quantities[ihit],
1280 diy_dpos,
1281 atmosphere_dim,
1282 rtp_pos);
1283
1284 }
1285 }
1286}
1287
1288/* Workspace method: Doxygen documentation will be auto-generated */
1290 Workspace& ws,
1291 Matrix& iy,
1292 const Vector& rtp_pos,
1293 const Index& stokes_dim,
1294 const Vector& f_grid,
1295 const Index& atmosphere_dim,
1296 const Vector& p_grid,
1297 const Vector& lat_grid,
1298 const Vector& lon_grid,
1299 const Tensor3& z_field,
1300 const Tensor3& t_field,
1301 const EnergyLevelMap& nlte_field,
1302 const Tensor4& vmr_field,
1303 const ArrayOfArrayOfSpeciesTag& abs_species,
1304 const Tensor3& wind_u_field,
1305 const Tensor3& wind_v_field,
1306 const Tensor3& wind_w_field,
1307 const Tensor3& mag_u_field,
1308 const Tensor3& mag_v_field,
1309 const Tensor3& mag_w_field,
1310 const Matrix& z_surface,
1311 const Vector& surface_scalar_reflectivity,
1312 const Vector& refellipsoid,
1313 const Tensor4& pnd_field,
1314 const ArrayOfTensor4& dpnd_field_dx,
1315 const ArrayOfString& scat_species,
1316 const ArrayOfArrayOfSingleScatteringData& scat_data,
1317 const Numeric& ppath_lmax,
1318 const Numeric& ppath_lraytrace,
1319 const Index& cloudbox_on,
1320 const ArrayOfIndex& cloudbox_limits,
1321 const Index& suns_do,
1322 const Index& gas_scattering_do,
1323 const Index& jacobian_do,
1324 const ArrayOfRetrievalQuantity& jacobian_quantities,
1325 const ArrayOfSun& suns,
1326 const Numeric& rte_alonglos_v,
1327 const String& iy_unit,
1328 const Agenda& propmat_clearsky_agenda,
1329 const Agenda& water_p_eq_agenda,
1330 const Agenda& gas_scattering_agenda,
1331 const Agenda& ppath_step_agenda,
1332 const Verbosity& verbosity) {
1333 //Allocate
1334 ArrayOfVector sun_rte_los(suns.nelem());
1335 ArrayOfMatrix transmitted_sunlight;
1336 ArrayOfArrayOfTensor3 dtransmitted_sunlight;
1337 ArrayOfPpath sun_ppaths(suns.nelem());
1338 ArrayOfIndex suns_visible(suns.nelem());
1339
1340 //Check for correct unit
1341 ARTS_USER_ERROR_IF(iy_unit != "1" && suns_do,
1342 "If suns are present only iy_unit=\"1\" can be used.");
1343 //Check size of iy
1344 chk_size("iy",iy,f_grid.nelem(),stokes_dim);
1345
1346 ARTS_USER_ERROR_IF(surface_scalar_reflectivity.nelem() != f_grid.nelem() &&
1347 surface_scalar_reflectivity.nelem() != 1,
1348 "The number of elements in *surface_scalar_reflectivity* should\n",
1349 "match length of *f_grid* or be 1.",
1350 "\n length of *f_grid* : ", f_grid.nelem(),
1351 "\n length of *surface_scalar_reflectivity* : ",
1352 surface_scalar_reflectivity.nelem());
1353 ARTS_USER_ERROR_IF(min(surface_scalar_reflectivity) < 0 ||
1354 max(surface_scalar_reflectivity) > 1,
1355 "All values in *surface_scalar_reflectivity* must be inside [0,1].");
1356 // Input checks
1357 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1358 chk_rte_pos(atmosphere_dim, rtp_pos);
1359
1360
1361 //do something only if there is a sun
1362 if (suns_do) {
1363 //get sun ppaths
1364 get_sun_ppaths(ws,
1365 sun_ppaths,
1366 suns_visible,
1367 sun_rte_los,
1368 rtp_pos,
1369 suns,
1370 f_grid,
1371 atmosphere_dim,
1372 p_grid,
1373 lat_grid,
1374 lon_grid,
1375 z_field,
1376 z_surface,
1377 refellipsoid,
1378 ppath_lmax,
1379 ppath_lraytrace,
1380 ppath_step_agenda,
1381 verbosity);
1382
1384 transmitted_sunlight,
1385 dtransmitted_sunlight,
1386 stokes_dim,
1387 f_grid,
1388 atmosphere_dim,
1389 p_grid,
1390 lat_grid,
1391 lon_grid,
1392 t_field,
1393 nlte_field,
1394 vmr_field,
1395 abs_species,
1396 wind_u_field,
1397 wind_v_field,
1398 wind_w_field,
1399 mag_u_field,
1400 mag_v_field,
1401 mag_w_field,
1402 cloudbox_on,
1403 cloudbox_limits,
1404 gas_scattering_do,
1405 1,
1406 sun_ppaths,
1407 suns,
1408 suns_visible,
1409 refellipsoid,
1410 pnd_field,
1411 dpnd_field_dx,
1412 scat_species,
1413 scat_data,
1414 jacobian_do,
1415 jacobian_quantities,
1416 propmat_clearsky_agenda,
1417 water_p_eq_agenda,
1418 gas_scattering_agenda,
1419 rte_alonglos_v,
1420 verbosity);
1421
1422 //loop over suns
1423 Vector specular_los, surface_normal;
1424
1425 for (Index i_sun = 0; i_sun < suns.nelem(); i_sun++) {
1426 if (suns_visible[i_sun]) {
1427 //sun_rte_los is the direction toward the sun, but we need from the sun
1428
1429 Vector incoming_los;
1430 mirror_los(incoming_los,sun_rte_los[i_sun], atmosphere_dim);
1431
1432 specular_losCalc(specular_los,
1433 surface_normal,
1434 rtp_pos,
1435 incoming_los,
1436 atmosphere_dim,
1437 lat_grid,
1438 lon_grid,
1439 refellipsoid,
1440 z_surface,
1441 0,
1442 verbosity);
1443
1444 // Only the first component of transmitted_sunlight is relevant.
1445 // Comment taken from surfaceLambertianSimple.
1446 // This follows VDISORT that refers to: K. L. Coulson, Polarization and Intensity of Light in
1447 // the Atmosphere (1989), page 229
1448 // (Thanks to Michael Kahnert for providing this information!)
1449 // Update: Is the above for a later edition? We have not found a copy of
1450 // that edition. In a 1988 version of the book, the relevant page seems
1451 // to be 232.
1452
1453 Vector iy_surface_direct(f_grid.nelem());
1454 Numeric r;
1455 for (Index i_freq = 0; i_freq < f_grid.nelem(); i_freq++) {
1456 if (i_freq == 0 || surface_scalar_reflectivity.nelem() > 1) {
1457 r = surface_scalar_reflectivity[i_freq];
1458 }
1459 iy_surface_direct[i_freq] = r / pi *
1460 transmitted_sunlight[i_sun](i_freq, 0) *
1461 cos(deg2rad(specular_los[0]));
1462 }
1463
1464 //Add the surface contribution of the direct part
1465 iy(joker, 0) += iy_surface_direct;
1466
1467
1468 }
1469 }
1470 }
1471}
1472
1473/* Workspace method: Doxygen documentation will be auto-generated */
1475 Matrix& iy,
1476 ArrayOfTensor3& diy_dx,
1477 Numeric& surface_skin_t,
1478 Matrix& surface_los,
1479 Tensor4& surface_rmatrix,
1480 Matrix& surface_emission,
1481 const Tensor3& iy_transmittance,
1482 const Index& iy_id,
1483 const Index& jacobian_do,
1484 const Index& suns_do,
1485 const Index& atmosphere_dim,
1486 const EnergyLevelMap& nlte_field,
1487 const Index& cloudbox_on,
1488 const Index& stokes_dim,
1489 const Vector& f_grid,
1490 const Vector& rtp_pos,
1491 const Vector& rtp_los,
1492 const Vector& rte_pos2,
1493 const String& iy_unit,
1494 const Agenda& iy_main_agenda,
1495 const Agenda& surface_rtprop_agenda,
1496 const Verbosity&) {
1497 // Input checks
1498 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1499 chk_rte_pos(atmosphere_dim, rtp_pos);
1500 chk_rte_los(atmosphere_dim, rtp_los);
1501
1502 // Call *surface_rtprop_agenda*
1504 surface_skin_t,
1505 surface_emission,
1506 surface_los,
1507 surface_rmatrix,
1508 f_grid,
1509 rtp_pos,
1510 rtp_los,
1511 surface_rtprop_agenda);
1512
1513 // Check output of *surface_rtprop_agenda*
1514 const Index nlos = surface_los.nrows();
1515 const Index nf = f_grid.nelem();
1516 //
1517 if (nlos) // if 0, blackbody ground and not all checks are needed
1518 {
1519 if (surface_los.ncols() != rtp_los.nelem())
1520 throw runtime_error("Number of columns in *surface_los* is not correct.");
1521 if (nlos != surface_rmatrix.nbooks())
1522 throw runtime_error(
1523 "Mismatch in size of *surface_los* and *surface_rmatrix*.");
1524 if (surface_rmatrix.npages() != nf)
1525 throw runtime_error(
1526 "Mismatch in size of *surface_rmatrix* and *f_grid*.");
1527 if (surface_rmatrix.nrows() != stokes_dim ||
1528 surface_rmatrix.ncols() != stokes_dim)
1529 throw runtime_error(
1530 "Mismatch between size of *surface_rmatrix* and *stokes_dim*.");
1531 }
1532 if (surface_emission.ncols() != stokes_dim)
1533 throw runtime_error(
1534 "Mismatch between size of *surface_emission* and *stokes_dim*.");
1535 if (surface_emission.nrows() != nf)
1536 throw runtime_error("Mismatch in size of *surface_emission* and f_grid*.");
1537
1538 // Variable to hold down-welling radiation
1539 Tensor3 I(nlos, nf, stokes_dim);
1540
1541 ArrayOfString iy_aux_var(0);
1542 if (suns_do) iy_aux_var.emplace_back("Direct radiation");
1543
1544 // Loop *surface_los*-es. If no such LOS, we are ready.
1545 if (nlos > 0) {
1546 for (Index ilos = 0; ilos < nlos; ilos++) {
1547 Vector los{surface_los(ilos, joker)};
1548
1549 // Include surface reflection matrix in *iy_transmittance*
1550 // If iy_transmittance is empty, this is interpreted as the
1551 // variable is not needed.
1552 //
1553 Tensor3 iy_trans_new;
1554 //
1555 if (iy_transmittance.npages()) {
1556 iy_transmittance_mult(iy_trans_new,
1557 iy_transmittance,
1558 surface_rmatrix(ilos, joker, joker, joker));
1559 }
1560
1561 // Calculate downwelling radiation for LOS ilos
1562 //
1563 {
1564 ArrayOfMatrix iy_aux;
1565 Ppath ppath;
1566 Vector geo_pos;
1567 Index iy_id_new = iy_id + ilos + 1;
1569 iy,
1570 iy_aux,
1571 ppath,
1572 diy_dx,
1573 geo_pos,
1574 0,
1575 iy_trans_new,
1576 iy_aux_var,
1577 iy_id_new,
1578 iy_unit,
1579 cloudbox_on,
1580 jacobian_do,
1581 f_grid,
1582 nlte_field,
1583 rtp_pos,
1584 los,
1585 rte_pos2,
1586 iy_main_agenda);
1587
1588 //For the case that a sun is present and the los is towards a sun, we
1589 //subtract the direct radiation, so that only the diffuse radiation is considered here.
1590 //If sun is within los iy_aux[0] is the incoming and attenuated direct (sun)
1591 //radiation at the surface. Otherwise it is zero.
1592 if (suns_do){
1593 iy-=iy_aux[0];
1594 }
1595
1596 }
1597
1598 if (iy.ncols() != stokes_dim || iy.nrows() != nf) {
1599 ostringstream os;
1600 os << "The size of *iy* returned from *" << iy_main_agenda.name()
1601 << "* is\n"
1602 << "not correct:\n"
1603 << " expected size = [" << nf << "," << stokes_dim << "]\n"
1604 << " size of iy = [" << iy.nrows() << "," << iy.ncols() << "]\n";
1605 throw runtime_error(os.str());
1606 }
1607
1608 I(ilos, joker, joker) = iy;
1609 }
1610 }
1611
1612 // Add up
1613 surface_calc(iy, I, surface_los, surface_rmatrix, surface_emission);
1614}
1615
1616/* Workspace method: Doxygen documentation will be auto-generated */
1618 Matrix& iy,
1619 ArrayOfTensor3& diy_dx,
1620 const Matrix& surface_los,
1621 const Tensor4& surface_rmatrix,
1622 const Matrix& surface_emission,
1623 const ArrayOfString& dsurface_names,
1624 const ArrayOfTensor4& dsurface_rmatrix_dx,
1625 const ArrayOfMatrix& dsurface_emission_dx,
1626 const Tensor3& iy_transmittance,
1627 const Index& iy_id,
1628 const Index& jacobian_do,
1629 const Index& suns_do,
1630 const ArrayOfRetrievalQuantity& jacobian_quantities,
1631 const Index& atmosphere_dim,
1632 const EnergyLevelMap& nlte_field,
1633 const Index& cloudbox_on,
1634 const Index& stokes_dim,
1635 const Vector& f_grid,
1636 const Vector& rtp_pos,
1637 const Vector& rtp_los,
1638 const Vector& rte_pos2,
1639 const String& iy_unit,
1640 const Agenda& iy_main_agenda,
1641 const Verbosity&) {
1642 // Input checks
1643 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1644 chk_rte_pos(atmosphere_dim, rtp_pos);
1645 chk_rte_los(atmosphere_dim, rtp_los);
1646
1647 // Check provided surface rtprop variables
1648 const Index nlos = surface_los.nrows();
1649 const Index nf = f_grid.nelem();
1650 //
1651 if (nlos) // if 0, blackbody ground and not all checks are needed
1652 {
1653 if (surface_los.ncols() != rtp_los.nelem())
1654 throw runtime_error("Number of columns in *surface_los* is not correct.");
1655 if (nlos != surface_rmatrix.nbooks())
1656 throw runtime_error(
1657 "Mismatch in size of *surface_los* and *surface_rmatrix*.");
1658 if (surface_rmatrix.npages() != nf)
1659 throw runtime_error(
1660 "Mismatch in size of *surface_rmatrix* and *f_grid*.");
1661 if (surface_rmatrix.nrows() != stokes_dim ||
1662 surface_rmatrix.ncols() != stokes_dim)
1663 throw runtime_error(
1664 "Mismatch between size of *surface_rmatrix* and *stokes_dim*.");
1665 }
1666 if (surface_emission.ncols() != stokes_dim)
1667 throw runtime_error(
1668 "Mismatch between size of *surface_emission* and *stokes_dim*.");
1669 if (surface_emission.nrows() != nf)
1670 throw runtime_error("Mismatch in size of *surface_emission* and f_grid*.");
1671
1672 // Variable to hold down-welling radiation
1673 Tensor3 I(nlos, nf, stokes_dim);
1674
1675 ArrayOfString iy_aux_var(0);
1676 if (suns_do) iy_aux_var.emplace_back("Direct radiation");
1677
1678 // Loop *surface_los*-es.
1679 if (nlos > 0) {
1680 for (Index ilos = 0; ilos < nlos; ilos++) {
1681 Vector los{surface_los(ilos, joker)};
1682
1683 // Include surface reflection matrix in *iy_transmittance*
1684 // If iy_transmittance is empty, this is interpreted as the
1685 // variable is not needed.
1686 //
1687 Tensor3 iy_trans_new;
1688 //
1689 if (iy_transmittance.npages()) {
1690 iy_transmittance_mult(iy_trans_new,
1691 iy_transmittance,
1692 surface_rmatrix(ilos, joker, joker, joker));
1693 }
1694
1695 // Calculate downwelling radiation for LOS ilos
1696 //
1697 {
1698 ArrayOfMatrix iy_aux;
1699 Ppath ppath;
1700 Vector geo_pos;
1702 iy,
1703 iy_aux,
1704 ppath,
1705 diy_dx,
1706 geo_pos,
1707 0,
1708 iy_trans_new,
1709 iy_aux_var,
1710 iy_id,
1711 iy_unit,
1712 cloudbox_on,
1713 jacobian_do,
1714 f_grid,
1715 nlte_field,
1716 rtp_pos,
1717 los,
1718 rte_pos2,
1719 iy_main_agenda);
1720
1721 //For the case that a sun is present and the los is towards a sun, we
1722 //subtract the direct radiation, so that only the diffuse radiation is considered here.
1723 //If sun is within los iy_aux[0] is the incoming and attenuated direct (sun)
1724 //radiation at the surface. Otherwise it is zero.
1725 if (suns_do){
1726 iy-=iy_aux[0];
1727 }
1728
1729 }
1730
1731 if (iy.ncols() != stokes_dim || iy.nrows() != nf) {
1732 ostringstream os;
1733 os << "The size of *iy* returned from *" << iy_main_agenda.name()
1734 << "* is\n"
1735 << "not correct:\n"
1736 << " expected size = [" << nf << "," << stokes_dim << "]\n"
1737 << " size of iy = [" << iy.nrows() << "," << iy.ncols() << "]\n";
1738 throw runtime_error(os.str());
1739 }
1740
1741 I(ilos, joker, joker) = iy;
1742 }
1743 }
1744
1745 // Add up
1746 surface_calc(iy, I, surface_los, surface_rmatrix, surface_emission);
1747
1748 // Surface Jacobians
1749 if (jacobian_do && dsurface_names.nelem()) {
1750 // Loop dsurface_names
1751 for (Index i = 0; i < dsurface_names.nelem(); i++) {
1752 // Error if derivatives not calculated
1753 // Or should we accept this?
1754 if (dsurface_emission_dx[i].empty() || dsurface_rmatrix_dx[i].empty()) {
1755 ostringstream os;
1756 os << "The derivatives for surface quantity: " << dsurface_names[i]
1757 << "\nwere not calculated by *iy_surface_agenda*.\n"
1758 << "That is, *dsurface_emission_dx* and/or *dsurface_rmatrix_dx*\n"
1759 << "are empty.";
1760 throw runtime_error(os.str());
1761 } else {
1762 // Find index among jacobian quantities
1763 Index ihit = -1;
1764 for (Index j = 0; j < jacobian_quantities.nelem() && ihit < 0; j++) {
1765 if (dsurface_names[i] == jacobian_quantities[j].Subtag()) {
1766 ihit = j;
1767 }
1768 }
1769 ARTS_ASSERT(ihit >= 0);
1770 // Derivative, as observed at the surface
1771 Matrix diy_dpos0, diy_dpos;
1772 surface_calc(diy_dpos0,
1773 I,
1774 surface_los,
1775 dsurface_rmatrix_dx[i],
1776 dsurface_emission_dx[i]);
1777 // Weight with transmission to sensor
1778 iy_transmittance_mult(diy_dpos, iy_transmittance, diy_dpos0);
1779 // Put into diy_dx
1780 diy_from_pos_to_rgrids(diy_dx[ihit],
1781 jacobian_quantities[ihit],
1782 diy_dpos,
1783 atmosphere_dim,
1784 rtp_pos);
1785 }
1786 }
1787 }
1788}
1789
1790/* Workspace method: Doxygen documentation will be auto-generated */
1791void specular_losCalcNoTopography(Vector& specular_los,
1792 Vector& surface_normal,
1793 const Vector& rtp_pos,
1794 const Vector& rtp_los,
1795 const Index& atmosphere_dim,
1796 const Verbosity&) {
1797 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1798 chk_rte_pos(atmosphere_dim, rtp_pos);
1799 chk_rte_los(atmosphere_dim, rtp_los);
1800
1801 surface_normal.resize(max(Index(1), atmosphere_dim - 1));
1802 specular_los.resize(max(Index(1), atmosphere_dim - 1));
1803
1804 if (atmosphere_dim == 1) {
1805 surface_normal[0] = 0;
1806 if (rtp_los[0] < 90) {
1807 throw runtime_error(
1808 "Invalid zenith angle. The zenith angle corresponds "
1809 "to observe the surface from below.");
1810 }
1811 specular_los[0] = 180 - rtp_los[0];
1812 }
1813
1814 else if (atmosphere_dim == 2) {
1815 specular_los[0] = sign(rtp_los[0]) * 180 - rtp_los[0];
1816 surface_normal[0] = 0;
1817 }
1818
1819 else if (atmosphere_dim == 3) {
1820 specular_los[0] = 180 - rtp_los[0];
1821 specular_los[1] = rtp_los[1];
1822 surface_normal[0] = 0;
1823 surface_normal[1] = 0;
1824 }
1825}
1826
1827/* Workspace method: Doxygen documentation will be auto-generated */
1828void specular_losCalc(Vector& specular_los,
1829 Vector& surface_normal,
1830 const Vector& rtp_pos,
1831 const Vector& rtp_los,
1832 const Index& atmosphere_dim,
1833 const Vector& lat_grid,
1834 const Vector& lon_grid,
1835 const Vector& refellipsoid,
1836 const Matrix& z_surface,
1837 const Index& ignore_surface_slope,
1838 const Verbosity& verbosity) {
1839 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1840 chk_rte_pos(atmosphere_dim, rtp_pos);
1841 chk_rte_los(atmosphere_dim, rtp_los);
1842 chk_if_in_range("ignore_surface_slope", ignore_surface_slope, 0, 1);
1843
1844 // Use special function if there is no slope, or it is ignored
1845 if (atmosphere_dim == 1 || ignore_surface_slope) {
1846 specular_losCalcNoTopography(specular_los,
1847 surface_normal,
1848 rtp_pos,
1849 rtp_los,
1850 atmosphere_dim,
1851 verbosity);
1852 return;
1853 }
1854
1855 surface_normal.resize(max(Index(1), atmosphere_dim - 1));
1856 specular_los.resize(max(Index(1), atmosphere_dim - 1));
1857
1858 if (atmosphere_dim == 2) {
1859 chk_interpolation_grids("Latitude interpolation", lat_grid, rtp_pos[1]);
1860 GridPos gp_lat;
1861 gridpos(gp_lat, lat_grid, rtp_pos[1]);
1862 Numeric c1; // Radial slope of the surface
1864 c1, lat_grid, refellipsoid, z_surface(joker, 0), gp_lat, rtp_los[0]);
1865 Vector itw(2);
1866 interpweights(itw, gp_lat);
1867 const Numeric rv_surface = refell2d(refellipsoid, lat_grid, gp_lat) +
1868 interp(itw, z_surface(joker, 0), gp_lat);
1869 surface_normal[0] = -plevel_angletilt(rv_surface, c1);
1870 if (abs(rtp_los[0] - surface_normal[0]) < 90) {
1871 throw runtime_error(
1872 "Invalid zenith angle. The zenith angle corresponds "
1873 "to observe the surface from below.");
1874 }
1875 specular_los[0] =
1876 sign(rtp_los[0]) * 180 - rtp_los[0] + 2 * surface_normal[0];
1877 }
1878
1879 else {
1880 // Calculate surface normal in South-North direction
1881 chk_interpolation_grids("Latitude interpolation", lat_grid, rtp_pos[1]);
1882 chk_interpolation_grids("Longitude interpolation", lon_grid, rtp_pos[2]);
1883 GridPos gp_lat, gp_lon;
1884 gridpos(gp_lat, lat_grid, rtp_pos[1]);
1885 gridpos(gp_lon, lon_grid, rtp_pos[2]);
1886 Numeric c1, c2;
1888 c1, c2, lat_grid, lon_grid, refellipsoid, z_surface, gp_lat, gp_lon, 0);
1889 Vector itw(4);
1890 interpweights(itw, gp_lat, gp_lon);
1891 const Numeric rv_surface = refell2d(refellipsoid, lat_grid, gp_lat) +
1892 interp(itw, z_surface, gp_lat, gp_lon);
1893 const Numeric zaSN = 90 - plevel_angletilt(rv_surface, c1);
1894 // The same for East-West
1895 plevel_slope_3d(c1,
1896 c2,
1897 lat_grid,
1898 lon_grid,
1899 refellipsoid,
1900 z_surface,
1901 gp_lat,
1902 gp_lon,
1903 90);
1904 const Numeric zaEW = 90 - plevel_angletilt(rv_surface, c1);
1905 // Convert to Cartesian, and determine normal by cross-product
1906 Vector tangentSN(3), tangentEW(3), normal(3);
1907 zaaa2cart(tangentSN[0], tangentSN[1], tangentSN[2], zaSN, 0);
1908 zaaa2cart(tangentEW[0], tangentEW[1], tangentEW[2], zaEW, 90);
1909 cross3(normal, tangentSN, tangentEW);
1910 // Convert rtp_los to cartesian and flip direction
1911 Vector di(3);
1912 zaaa2cart(di[0], di[1], di[2], rtp_los[0], rtp_los[1]);
1913 di *= -1;
1914 // Set LOS normal vector
1915 cart2zaaa(
1916 surface_normal[0], surface_normal[1], normal[0], normal[1], normal[2]);
1917 if (abs(rtp_los[0] - surface_normal[0]) < 90) {
1918 throw runtime_error(
1919 "Invalid zenith angle. The zenith angle corresponds "
1920 "to observe the surface from below.");
1921 }
1922 // Specular direction is 2(dn*di)dn-di, where dn is the normal vector
1923 Vector speccart(3);
1924 const Numeric fac = 2 * (normal * di);
1925 for (Index i = 0; i < 3; i++) {
1926 speccart[i] = fac * normal[i] - di[i];
1927 }
1928 cart2zaaa(specular_los[0],
1929 specular_los[1],
1930 speccart[0],
1931 speccart[1],
1932 speccart[2]);
1933 }
1934}
1935
1936/* Workspace method: Doxygen documentation will be auto-generated */
1937void surfaceBlackbody(Matrix& surface_los,
1938 Tensor4& surface_rmatrix,
1939 Matrix& surface_emission,
1940 const Index& atmosphere_dim,
1941 const Vector& f_grid,
1942 const Index& stokes_dim,
1943 const Vector& rtp_pos,
1944 const Vector& rtp_los,
1945 const Numeric& surface_skin_t,
1946 const Verbosity& verbosity) {
1947 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
1948 chk_rte_pos(atmosphere_dim, rtp_pos);
1949 chk_rte_los(atmosphere_dim, rtp_los);
1950 chk_not_negative("surface_skin_t", surface_skin_t);
1951
1953 out2 << " Sets variables to model a blackbody surface with a temperature "
1954 << " of " << surface_skin_t << " K.\n";
1955
1956 surface_los.resize(0, 0);
1957 surface_rmatrix.resize(0, 0, 0, 0);
1958
1959 const Index nf = f_grid.nelem();
1960
1961 Vector b(nf);
1962 planck(b, f_grid, surface_skin_t);
1963
1964 surface_emission.resize(nf, stokes_dim);
1965 surface_emission = 0.0;
1966
1967 for (Index iv = 0; iv < nf; iv++) {
1968 surface_emission(iv, 0) = b[iv];
1969 for (Index is = 1; is < stokes_dim; is++) {
1970 surface_emission(iv, is) = 0;
1971 }
1972 }
1973}
1974
1975/* Workspace method: Doxygen documentation will be auto-generated */
1976void surfaceFastem(Matrix& surface_los,
1977 Tensor4& surface_rmatrix,
1978 Matrix& surface_emission,
1979 const Index& atmosphere_dim,
1980 const Index& stokes_dim,
1981 const Vector& f_grid,
1982 const Vector& rtp_pos,
1983 const Vector& rtp_los,
1984 const Numeric& surface_skin_t,
1985 const Numeric& salinity,
1986 const Numeric& wind_speed,
1987 const Numeric& wind_direction,
1988 const Vector& transmittance,
1989 const Index& fastem_version,
1990 const Verbosity& verbosity) {
1991 // Input checks
1992 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1993 chk_rte_pos(atmosphere_dim, rtp_pos);
1994 chk_rte_los(atmosphere_dim, rtp_los);
1995 chk_if_in_range("wind_direction", wind_direction, -180, 180);
1996
1997 const Index nf = f_grid.nelem();
1998
1999 // Determine specular direction
2000 Vector specular_los, surface_normal;
2001 specular_losCalcNoTopography(specular_los,
2002 surface_normal,
2003 rtp_pos,
2004 rtp_los,
2005 atmosphere_dim,
2006 verbosity);
2007
2008 // Determine relative azimuth
2009 //
2010 // According to email from James Hocking, UkMet:
2011 // All azimuth angles are defined as being measured clockwise from north
2012 // (i.e. if the projection onto the Earth's surface of the view path lies due
2013 // north the azimuth angle is 0 and if it lies due east the azimuth angle is
2014 // 90 degrees). The relative azimuth is the wind direction (azimuth) angle
2015 // minus the satellite azimuth angle.
2016 Numeric rel_azimuth = wind_direction; // Always true for 1D
2017 if (atmosphere_dim == 2 && rtp_los[0] < 0) {
2018 rel_azimuth -= 180;
2019 resolve_lon(rel_azimuth, -180, 180);
2020 } else if (atmosphere_dim == 3) {
2021 rel_azimuth -= rtp_los[1];
2022 resolve_lon(rel_azimuth, -180, 180);
2023 }
2024
2025 // Call FASTEM
2026 Matrix emissivity, reflectivity;
2027 FastemStandAlone(emissivity,
2028 reflectivity,
2029 f_grid,
2030 surface_skin_t,
2031 abs(rtp_los[0]),
2032 salinity,
2033 wind_speed,
2034 rel_azimuth,
2035 transmittance,
2036 fastem_version,
2037 verbosity);
2038
2039 // Set surface_los
2040 surface_los.resize(1, specular_los.nelem());
2041 surface_los(0, joker) = specular_los;
2042
2043 // Surface emission
2044 //
2045 Vector b(nf);
2046 planck(b, f_grid, surface_skin_t);
2047 //
2048 surface_emission.resize(nf, stokes_dim);
2049 for (Index i = 0; i < nf; i++) {
2050 // I
2051 surface_emission(i, 0) = b[i] * 0.5 * (emissivity(i, 0) + emissivity(i, 1));
2052 // Q
2053 if (stokes_dim >= 2) {
2054 surface_emission(i, 1) =
2055 b[i] * 0.5 * (emissivity(i, 0) - emissivity(i, 1));
2056 }
2057 // U and V
2058 for (Index j = 2; j < stokes_dim; j++) {
2059 surface_emission(i, j) = b[i] * emissivity(i, j);
2060 }
2061 }
2062
2063 // Surface reflectivity matrix
2064 //
2065 surface_rmatrix.resize(1, nf, stokes_dim, stokes_dim);
2066 surface_rmatrix = 0.0;
2067 for (Index iv = 0; iv < nf; iv++) {
2068 surface_rmatrix(0, iv, 0, 0) =
2069 0.5 * (reflectivity(iv, 0) + reflectivity(iv, 1));
2070 if (stokes_dim >= 2) {
2071 surface_rmatrix(0, iv, 0, 1) =
2072 0.5 * (reflectivity(iv, 0) - reflectivity(iv, 1));
2073 ;
2074 surface_rmatrix(0, iv, 1, 0) = surface_rmatrix(0, iv, 0, 1);
2075 surface_rmatrix(0, iv, 1, 1) = surface_rmatrix(0, iv, 0, 0);
2076
2077 for (Index i = 2; i < stokes_dim; i++) {
2078 surface_rmatrix(0, iv, i, i) = surface_rmatrix(0, iv, 0, 0);
2079 }
2080 }
2081 }
2082}
2083
2084/* Workspace method: Doxygen documentation will be auto-generated */
2085void surfaceMapToLinearPolarisation(Matrix& surface_emission,
2086 Tensor4& surface_rmatrix,
2087 const Index& stokes_dim,
2088 const Numeric& pol_angle,
2089 const Verbosity& ) {
2090 ARTS_USER_ERROR_IF (stokes_dim != 1,
2091 "You should only use this method where the main calculations are "
2092 "done with *stokes_dim* set to 1.");
2093
2094 Index local_stokes = surface_emission.ncols();
2095 ARTS_USER_ERROR_IF (local_stokes < 2,
2096 "This method requires that the input surface proporties match a Stokes "
2097 "dimension of 2-4. Incoming *surface_emission* matches stokes_dim=1.");
2098
2099 const Index nf = surface_emission.nrows();
2100 const Index nlos = surface_rmatrix.nbooks();
2101 Matrix se = surface_emission;
2102 Tensor4 sr = surface_rmatrix;
2103 surface_emission.resize(nf, 1);
2104 surface_rmatrix.resize(nlos, nf, 1, 1);
2105
2106 if (local_stokes == 2) {
2107 const Numeric alpha = DEG2RAD * pol_angle;
2108 const Numeric c2 = pow( cos(alpha), 2 );
2109 const Numeric s2 = pow( sin(alpha), 2 );
2110 for (Index f=0; f<nf; ++f) {
2111 const Numeric bv = se(f,0) + se(f,1); // Note that we have to multiply with 2
2112 const Numeric bh = se(f,0) - se(f,1); // as we place a single pol as I
2113 surface_emission(f,0) = c2*bv + s2*bh;
2114 for (Index l=0; l<nlos; ++l) {
2115 const Numeric rv = sr(l,f,0,0) + sr(l,f,1,0);
2116 const Numeric rh = sr(l,f,0,0) - sr(l,f,1,0);
2117 surface_rmatrix(l,f,0,0) = c2*rv + s2*rh;
2118 }
2119 }
2120 } else {
2121 Matrix Cm, Cs, Lp, Lm;
2122 mueller_stokes2modif(Cm, local_stokes);
2123 mueller_modif2stokes(Cs, local_stokes);
2124 mueller_rotation(Lp, local_stokes, pol_angle);
2125 mueller_rotation(Lm, local_stokes, -pol_angle);
2126 Matrix Mleft(local_stokes,local_stokes), Mright(local_stokes,local_stokes);
2127 mult(Mleft, Cm, Lp);
2128 mult(Mright, Lm, Cs);
2129 //
2130 Vector Vr(local_stokes);
2131 Matrix Tmp(local_stokes,local_stokes), Mr(local_stokes,local_stokes);
2132 //
2133 for (Index f=0; f<nf; ++f) {
2134 mult(Vr, Mleft, se(f,joker) ); // Note that we have to multiply with 2
2135 surface_emission(f,0) = 2.0 * Vr[0]; // as we place a single pol as I
2136 for (Index l=0; l<nlos; ++l) {
2137 mult(Tmp, sr(l,f,joker,joker), Mright);
2138 mult(Mr, Mleft, Tmp);
2139 surface_rmatrix(l,f,0,0) = Tmp(0,0);
2140 }
2141 }
2142 }
2143}
2144
2145/* Workspace method: Doxygen documentation will be auto-generated */
2146void surfaceTelsem(Matrix& surface_los,
2147 Tensor4& surface_rmatrix,
2148 Matrix& surface_emission,
2149 const Index& atmosphere_dim,
2150 const Index& stokes_dim,
2151 const Vector& f_grid,
2152 const Vector& lat_grid,
2153 const Vector& lat_true,
2154 const Vector& lon_true,
2155 const Vector& rtp_pos,
2156 const Vector& rtp_los,
2157 const Vector& specular_los,
2158 const Numeric& surface_skin_t,
2159 const TelsemAtlas& atlas,
2160 const Numeric& r_min,
2161 const Numeric& r_max,
2162 const Numeric& d_max,
2163 const Verbosity& verbosity) {
2164 // Input checks
2165 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2166 chk_rte_pos(atmosphere_dim, rtp_pos);
2167 chk_rte_los(atmosphere_dim, rtp_los);
2168 chk_rte_los(atmosphere_dim, specular_los);
2170 "surface skin temperature", surface_skin_t, 190.0, 373.0);
2171
2172 const Index nf = f_grid.nelem();
2173 Matrix surface_rv_rh(nf, 2);
2174
2175 // Lookup cell in atlas.
2176
2177 Numeric lat, lon;
2179 lat, lon, atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
2180 chk_if_in_range("Latitude input to TELSEM2", lat, -90.0, 90.0);
2181
2182 lon = fmod(lon, 360.0);
2183 if (lon < 0) {
2184 lon += 360.0;
2185 }
2186 chk_if_in_range("Longitude input to TELSEM2", lon, 0.0, 360.0);
2187
2188 Index cellnumber = atlas.calc_cellnum(lat, lon);
2189 // Check if cell is in atlas.
2190 if (!atlas.contains(cellnumber)) {
2191 if (d_max <= 0.0) {
2192 throw std::runtime_error(
2193 "Given coordinates are not contained in "
2194 " TELSEM atlas. To enable nearest neighbour"
2195 "interpolation set *d_max* to a positive "
2196 "value.");
2197 } else {
2198 cellnumber = atlas.calc_cellnum_nearest_neighbor(lat, lon);
2199 Numeric lat_nn, lon_nn;
2200 std::tie(lat_nn, lon_nn) = atlas.get_coordinates(cellnumber);
2201 Numeric d = sphdist(lat, lon, lat_nn, lon_nn);
2202 if (d > d_max) {
2203 std::ostringstream out{};
2204 out << "Distance of nearest neighbour exceeds provided limit (";
2205 out << d << " > " << d_max << ").";
2206 throw std::runtime_error(out.str());
2207 }
2208 }
2209 }
2210
2211 Index class1 = atlas.get_class1(cellnumber);
2212 Index class2 = atlas.get_class2(cellnumber);
2213
2214 Vector emis_h = atlas.get_emis_h(cellnumber);
2215 Vector emis_v = atlas.get_emis_v(cellnumber);
2216 Numeric e_v, e_h;
2217
2218 // Incidence angle
2219 const Numeric theta = calc_incang(rtp_los, specular_los);
2220 ARTS_ASSERT(theta <= 90);
2221
2222 for (Index i = 0; i < nf; ++i) {
2223 if (f_grid[i] < 5e9)
2224 throw std::runtime_error("Only frequency >= 5 GHz are allowed");
2225 if (f_grid[i] > 900e9)
2226 throw std::runtime_error("Only frequency <= 900 GHz are allowed");
2227
2228 // For frequencies 700 <= f <= 900 GHz use 700 GHz to avoid extrapolation
2229 // outside of TELSEM specifications.
2230 Numeric f = std::min(f_grid[i], 700e9) * 1e-9;
2231 std::tie(e_v, e_h) =
2232 atlas.emis_interp(theta, f, class1, class2, emis_v, emis_h);
2233
2234 surface_rv_rh(i, 0) = min(max(1.0 - e_v, r_min), r_max);
2235 surface_rv_rh(i, 1) = min(max(1.0 - e_h, r_min), r_max);
2236 }
2237
2238 surfaceFlatRvRh(surface_los,
2239 surface_rmatrix,
2240 surface_emission,
2241 f_grid,
2242 stokes_dim,
2243 atmosphere_dim,
2244 rtp_pos,
2245 rtp_los,
2246 specular_los,
2247 surface_skin_t,
2248 surface_rv_rh,
2249 verbosity);
2250}
2251
2252/* Workspace method: Doxygen documentation will be auto-generated */
2253void surfaceTessem(Matrix& surface_los,
2254 Tensor4& surface_rmatrix,
2255 Matrix& surface_emission,
2256 const Index& atmosphere_dim,
2257 const Index& stokes_dim,
2258 const Vector& f_grid,
2259 const Vector& rtp_pos,
2260 const Vector& rtp_los,
2261 const Numeric& surface_skin_t,
2262 const TessemNN& net_h,
2263 const TessemNN& net_v,
2264 const Numeric& salinity,
2265 const Numeric& wind_speed,
2266 const Verbosity& verbosity) {
2267 // Input checks
2268 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2269 chk_rte_pos(atmosphere_dim, rtp_pos);
2270 chk_rte_los(atmosphere_dim, rtp_los);
2272 "surface skin temperature", surface_skin_t, 260.0, 373.0);
2273 chk_if_in_range_exclude_high("salinity", salinity, 0, 1);
2274 chk_if_in_range_exclude_high("wind speed", wind_speed, 0, 100);
2275
2276 // Determine specular direction
2277 Vector specular_los, surface_normal;
2278 specular_losCalcNoTopography(specular_los,
2279 surface_normal,
2280 rtp_pos,
2281 rtp_los,
2282 atmosphere_dim,
2283 verbosity);
2284
2285 // TESSEM in and out
2286 //
2287 Vector out(2);
2288 VectorView e_h = out[Range(0, 1)];
2289 VectorView e_v = out[Range(1, 1)];
2290 //
2291 Vector in(5);
2292 in[1] = 180.0 - abs(rtp_los[0]);
2293 in[2] = wind_speed;
2294 in[3] = surface_skin_t;
2295 in[4] = salinity;
2296
2297 // Get Rv and Rh
2298 //
2299 const Index nf = f_grid.nelem();
2300 Matrix surface_rv_rh(nf, 2);
2301 //
2302 for (Index i = 0; i < nf; ++i) {
2303 if (f_grid[i] < 5e9)
2304 throw std::runtime_error("Only frequency >= 5 GHz are allowed");
2305 if (f_grid[i] > 900e9)
2306 throw std::runtime_error("Only frequency <= 900 GHz are allowed");
2307
2308 in[0] = f_grid[i];
2309
2310 tessem_prop_nn(e_h, net_h, in);
2311 tessem_prop_nn(e_v, net_v, in);
2312
2313 surface_rv_rh(i, 0) = min(max(1 - e_v[0], (Numeric)0), (Numeric)1);
2314 surface_rv_rh(i, 1) = min(max(1 - e_h[0], (Numeric)0), (Numeric)1);
2315 }
2316
2317 surfaceFlatRvRh(surface_los,
2318 surface_rmatrix,
2319 surface_emission,
2320 f_grid,
2321 stokes_dim,
2322 atmosphere_dim,
2323 rtp_pos,
2324 rtp_los,
2325 specular_los,
2326 surface_skin_t,
2327 surface_rv_rh,
2328 verbosity);
2329}
2330
2331/* Workspace method: Doxygen documentation will be auto-generated */
2332void surfaceFlatRefractiveIndex(Matrix& surface_los,
2333 Tensor4& surface_rmatrix,
2334 Matrix& surface_emission,
2335 const Vector& f_grid,
2336 const Index& stokes_dim,
2337 const Index& atmosphere_dim,
2338 const Vector& rtp_pos,
2339 const Vector& rtp_los,
2340 const Vector& specular_los,
2341 const Numeric& surface_skin_t,
2342 const GriddedField3& surface_complex_refr_index,
2343 const Verbosity& verbosity) {
2346
2347 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2348 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
2349 chk_rte_pos(atmosphere_dim, rtp_pos);
2350 chk_rte_los(atmosphere_dim, rtp_los);
2351 chk_rte_los(atmosphere_dim, specular_los);
2352 chk_not_negative("surface_skin_t", surface_skin_t);
2353
2354 // Interpolate *surface_complex_refr_index*
2355 //
2356 const Index nf = f_grid.nelem();
2357 //
2358 Matrix n_real(nf, 1), n_imag(nf, 1);
2359 //
2360 complex_n_interp(n_real,
2361 n_imag,
2362 surface_complex_refr_index,
2363 "surface_complex_refr_index",
2364 f_grid,
2365 Vector(1, surface_skin_t));
2366
2367 out2 << " Sets variables to model a flat surface\n";
2368 out3 << " surface temperature: " << surface_skin_t << " K.\n";
2369
2370 surface_los.resize(1, specular_los.nelem());
2371 surface_los(0, joker) = specular_los;
2372
2373 surface_emission.resize(nf, stokes_dim);
2374 surface_rmatrix.resize(1, nf, stokes_dim, stokes_dim);
2375
2376 // Incidence angle
2377 const Numeric incang = calc_incang(rtp_los, specular_los);
2378 ARTS_ASSERT(incang <= 90);
2379
2380 // Complex (amplitude) reflection coefficients
2381 Complex Rv, Rh;
2382
2383 for (Index iv = 0; iv < nf; iv++) {
2384 // Set n2 (refractive index of surface medium)
2385 Complex n2(n_real(iv, 0), n_imag(iv, 0));
2386
2387 // Amplitude reflection coefficients
2388 fresnel(Rv, Rh, Numeric(1.0), n2, incang);
2389
2390 // Fill reflection matrix and emission vector
2391 surface_specular_R_and_b(surface_rmatrix(0, iv, joker, joker),
2392 surface_emission(iv, joker),
2393 Rv,
2394 Rh,
2395 f_grid[iv],
2396 stokes_dim,
2397 surface_skin_t);
2398 }
2399}
2400
2401/* Workspace method: Doxygen documentation will be auto-generated */
2402void surfaceFlatReflectivity(Matrix& surface_los,
2403 Tensor4& surface_rmatrix,
2404 Matrix& surface_emission,
2405 const Vector& f_grid,
2406 const Index& stokes_dim,
2407 const Index& atmosphere_dim,
2408 const Vector& rtp_pos,
2409 const Vector& rtp_los,
2410 const Vector& specular_los,
2411 const Numeric& surface_skin_t,
2412 const Tensor3& surface_reflectivity,
2413 const Verbosity& verbosity) {
2415
2416 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2417 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
2418 chk_rte_pos(atmosphere_dim, rtp_pos);
2419 chk_rte_los(atmosphere_dim, rtp_los);
2420 chk_rte_los(atmosphere_dim, specular_los);
2421 chk_not_negative("surface_skin_t", surface_skin_t);
2422
2423 const Index nf = f_grid.nelem();
2424
2425 if (surface_reflectivity.nrows() != stokes_dim &&
2426 surface_reflectivity.ncols() != stokes_dim) {
2427 ostringstream os;
2428 os << "The number of rows and columnss in *surface_reflectivity* must\n"
2429 << "match *stokes_dim*."
2430 << "\n stokes_dim : " << stokes_dim
2431 << "\n number of rows in *surface_reflectivity* : "
2432 << surface_reflectivity.nrows()
2433 << "\n number of columns in *surface_reflectivity* : "
2434 << surface_reflectivity.ncols() << "\n";
2435 throw runtime_error(os.str());
2436 }
2437
2438 if (surface_reflectivity.npages() != nf &&
2439 surface_reflectivity.npages() != 1) {
2440 ostringstream os;
2441 os << "The number of pages in *surface_reflectivity* should\n"
2442 << "match length of *f_grid* or be 1."
2443 << "\n length of *f_grid* : " << nf
2444 << "\n dimension of *surface_reflectivity* : "
2445 << surface_reflectivity.npages() << "\n";
2446 throw runtime_error(os.str());
2447 }
2448
2449 out2 << " Sets variables to model a flat surface\n";
2450
2451 surface_los.resize(1, specular_los.nelem());
2452 surface_los(0, joker) = specular_los;
2453
2454 surface_emission.resize(nf, stokes_dim);
2455 surface_rmatrix.resize(1, nf, stokes_dim, stokes_dim);
2456
2457 Matrix R, IR(stokes_dim, stokes_dim);
2458
2459 Vector b(nf);
2460 planck(b, f_grid, surface_skin_t);
2461
2462 Vector B(stokes_dim, 0);
2463
2464 for (Index iv = 0; iv < nf; iv++) {
2465 if (iv == 0 || surface_reflectivity.npages() > 1) {
2466 R = surface_reflectivity(iv, joker, joker);
2467 for (Index i = 0; i < stokes_dim; i++) {
2468 for (Index j = 0; j < stokes_dim; j++) {
2469 if (i == j) {
2470 IR(i, j) = 1 - R(i, j);
2471 } else {
2472 IR(i, j) = -R(i, j);
2473 }
2474 }
2475 }
2476 }
2477
2478 surface_rmatrix(0, iv, joker, joker) = R;
2479
2480 B[0] = b[iv];
2481 mult(surface_emission(iv, joker), IR, B);
2482 }
2483}
2484
2485/* Workspace method: Doxygen documentation will be auto-generated */
2486void surfaceFlatRvRh(Matrix& surface_los,
2487 Tensor4& surface_rmatrix,
2488 Matrix& surface_emission,
2489 const Vector& f_grid,
2490 const Index& stokes_dim,
2491 const Index& atmosphere_dim,
2492 const Vector& rtp_pos,
2493 const Vector& rtp_los,
2494 const Vector& specular_los,
2495 const Numeric& surface_skin_t,
2496 const Matrix& surface_rv_rh,
2497 const Verbosity&) {
2498 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2499 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
2500 chk_rte_pos(atmosphere_dim, rtp_pos);
2501 chk_rte_los(atmosphere_dim, rtp_los);
2502 chk_rte_los(atmosphere_dim, specular_los);
2503 chk_not_negative("surface_skin_t", surface_skin_t);
2504
2505 const Index nf = f_grid.nelem();
2506
2507 if (surface_rv_rh.ncols() != 2) {
2508 ostringstream os;
2509 os << "The number of columns in *surface_rv_rh* must be two,\n"
2510 << "but the actual number of columns is " << surface_rv_rh.ncols()
2511 << "\n";
2512 throw runtime_error(os.str());
2513 }
2514
2515 if (surface_rv_rh.nrows() != nf && surface_rv_rh.nrows() != 1) {
2516 ostringstream os;
2517 os << "The number of rows in *surface_rv_rh* should\n"
2518 << "match length of *f_grid* or be 1."
2519 << "\n length of *f_grid* : " << nf
2520 << "\n rows in *surface_rv_rh* : " << surface_rv_rh.nrows() << "\n";
2521 throw runtime_error(os.str());
2522 }
2523
2524 if (min(surface_rv_rh) < 0 || max(surface_rv_rh) > 1) {
2525 throw runtime_error("All values in *surface_rv_rh* must be inside [0,1].");
2526 }
2527
2528 surface_los.resize(1, specular_los.nelem());
2529 surface_los(0, joker) = specular_los;
2530
2531 surface_emission.resize(nf, stokes_dim);
2532 surface_rmatrix.resize(1, nf, stokes_dim, stokes_dim);
2533
2534 surface_emission = 0;
2535 surface_rmatrix = 0;
2536
2537 Vector b(nf);
2538 planck(b, f_grid, surface_skin_t);
2539
2540 Numeric rmean = 0.0, rdiff = 0.0;
2541
2542 for (Index iv = 0; iv < nf; iv++) {
2543 if (iv == 0 || surface_rv_rh.nrows() > 1) {
2544 rmean = 0.5 * (surface_rv_rh(iv, 0) + surface_rv_rh(iv, 1));
2545 rdiff = 0.5 * (surface_rv_rh(iv, 0) - surface_rv_rh(iv, 1));
2546 }
2547
2548 surface_emission(iv, 0) = (1.0 - rmean) * b[iv];
2549 surface_rmatrix(0, iv, 0, 0) = rmean;
2550
2551 if (stokes_dim > 1) {
2552 surface_emission(iv, 1) = -rdiff * b[iv];
2553
2554 surface_rmatrix(0, iv, 0, 1) = rdiff;
2555 surface_rmatrix(0, iv, 1, 0) = rdiff;
2556 surface_rmatrix(0, iv, 1, 1) = rmean;
2557
2558 for (Index i = 2; i < stokes_dim; i++) {
2559 surface_rmatrix(0, iv, i, i) = rmean;
2560 }
2561 }
2562 }
2563}
2564
2565/* Workspace method: Doxygen documentation will be auto-generated */
2566void surfaceFlatScalarReflectivity(Matrix& surface_los,
2567 Tensor4& surface_rmatrix,
2568 Matrix& surface_emission,
2569 const Vector& f_grid,
2570 const Index& stokes_dim,
2571 const Index& atmosphere_dim,
2572 const Vector& rtp_pos,
2573 const Vector& rtp_los,
2574 const Vector& specular_los,
2575 const Numeric& surface_skin_t,
2576 const Vector& surface_scalar_reflectivity,
2577 const Verbosity&) {
2578 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2579 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
2580 chk_rte_pos(atmosphere_dim, rtp_pos);
2581 chk_rte_los(atmosphere_dim, rtp_los);
2582 chk_rte_los(atmosphere_dim, specular_los);
2583 chk_not_negative("surface_skin_t", surface_skin_t);
2584
2585 const Index nf = f_grid.nelem();
2586
2587 if (surface_scalar_reflectivity.nelem() != nf &&
2588 surface_scalar_reflectivity.nelem() != 1) {
2589 ostringstream os;
2590 os << "The number of elements in *surface_scalar_reflectivity* should\n"
2591 << "match length of *f_grid* or be 1."
2592 << "\n length of *f_grid* : " << nf
2593 << "\n length of *surface_scalar_reflectivity* : "
2594 << surface_scalar_reflectivity.nelem() << "\n";
2595 throw runtime_error(os.str());
2596 }
2597
2598 if (min(surface_scalar_reflectivity) < 0 ||
2599 max(surface_scalar_reflectivity) > 1) {
2600 throw runtime_error(
2601 "All values in *surface_scalar_reflectivity* must be inside [0,1].");
2602 }
2603
2604 surface_los.resize(1, specular_los.nelem());
2605 surface_los(0, joker) = specular_los;
2606
2607 surface_emission.resize(nf, stokes_dim);
2608 surface_rmatrix.resize(1, nf, stokes_dim, stokes_dim);
2609
2610 surface_emission = 0;
2611 surface_rmatrix = 0;
2612
2613 Vector b(nf);
2614 planck(b, f_grid, surface_skin_t);
2615
2616 Numeric r = 0.0;
2617
2618 for (Index iv = 0; iv < nf; iv++) {
2619 if (iv == 0 || surface_scalar_reflectivity.nelem() > 1) {
2620 r = surface_scalar_reflectivity[iv];
2621 }
2622
2623 surface_emission(iv, 0) = (1.0 - r) * b[iv];
2624 surface_rmatrix(0, iv, 0, 0) = r;
2625 for (Index i = 1; i < stokes_dim; i++) {
2626 surface_rmatrix(0, iv, i, i) = r;
2627 }
2628 }
2629}
2630
2631/* Workspace method: Doxygen documentation will be auto-generated */
2632void surfaceLambertianSimple(Matrix& surface_los,
2633 Tensor4& surface_rmatrix,
2634 Matrix& surface_emission,
2635 const Vector& f_grid,
2636 const Index& stokes_dim,
2637 const Index& atmosphere_dim,
2638 const Vector& rtp_pos,
2639 const Vector& rtp_los,
2640 const Vector& surface_normal,
2641 const Numeric& surface_skin_t,
2642 const Vector& surface_scalar_reflectivity,
2643 const Index& lambertian_nza,
2644 const Numeric& za_pos,
2645 const Verbosity&) {
2646 const Index nf = f_grid.nelem();
2647
2648 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2649 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
2650 chk_rte_pos(atmosphere_dim, rtp_pos);
2651 chk_rte_los(atmosphere_dim, rtp_los);
2652 chk_not_negative("surface_skin_t", surface_skin_t);
2653 chk_if_in_range("za_pos", za_pos, 0, 1);
2654
2655 if (surface_scalar_reflectivity.nelem() != nf &&
2656 surface_scalar_reflectivity.nelem() != 1) {
2657 ostringstream os;
2658 os << "The number of elements in *surface_scalar_reflectivity* should\n"
2659 << "match length of *f_grid* or be 1."
2660 << "\n length of *f_grid* : " << nf
2661 << "\n length of *surface_scalar_reflectivity* : "
2662 << surface_scalar_reflectivity.nelem() << "\n";
2663 throw runtime_error(os.str());
2664 }
2665
2666 if (min(surface_scalar_reflectivity) < 0 ||
2667 max(surface_scalar_reflectivity) > 1) {
2668 throw runtime_error(
2669 "All values in *surface_scalar_reflectivity* must be inside [0,1].");
2670 }
2671
2672 // Allocate and init everything to zero
2673 //
2674 surface_los.resize(lambertian_nza, rtp_los.nelem());
2675 surface_rmatrix.resize(lambertian_nza, nf, stokes_dim, stokes_dim);
2676 surface_emission.resize(nf, stokes_dim);
2677 //
2678 surface_los = 0.0;
2679 surface_rmatrix = 0.0;
2680 surface_emission = 0.0;
2681
2682 // Help variables
2683 //
2684 const Numeric dza = (90.0 - abs(surface_normal[0])) / (Numeric)lambertian_nza;
2685 const Vector za_lims=uniform_grid(0.0, lambertian_nza + 1, dza);
2686
2687 // surface_los
2688 for (Index ip = 0; ip < lambertian_nza; ip++) {
2689 surface_los(ip, 0) = za_lims[ip] + za_pos * dza;
2690 if (atmosphere_dim == 2) {
2691 if (rtp_los[0] < 0) {
2692 surface_los(ip, 0) *= -1.0;
2693 }
2694
2695 } else if (atmosphere_dim == 3) {
2696 surface_los(ip, 1) = rtp_los[1];
2697 }
2698 }
2699
2700 Vector b(nf);
2701 planck(b, f_grid, surface_skin_t);
2702
2703 // Loop frequencies and set remaining values
2704 //
2705 Numeric r = 0.0;
2706 //
2707 for (Index iv = 0; iv < nf; iv++) {
2708 // Get reflectivity
2709 if (iv == 0 || surface_scalar_reflectivity.nelem() > 1) {
2710 r = surface_scalar_reflectivity[iv];
2711 }
2712
2713 // surface_rmatrix:
2714 // Only element (0,0) is set to be non-zero. This follows VDISORT
2715 // that refers to: K. L. Coulson, Polarization and Intensity of Light in
2716 // the Atmosphere (1989), page 229
2717 // (Thanks to Michael Kahnert for providing this information!)
2718 // Update: Is the above for a later edition? We have not found a copy of
2719 // that edition. In a 1988 version of the book, the relevant page seems
2720 // to be 232.
2721 for (Index ip = 0; ip < lambertian_nza; ip++) {
2722 const Numeric w =
2723 r * 0.5 *
2724 (cos(2 * DEG2RAD * za_lims[ip]) - cos(2 * DEG2RAD * za_lims[ip + 1]));
2725 surface_rmatrix(ip, iv, 0, 0) = w;
2726 }
2727
2728 // surface_emission
2729 surface_emission(iv, 0) = (1 - r) * b[iv];
2730 }
2731}
2732
2733/* Workspace method: Doxygen documentation will be auto-generated */
2735 GriddedField3& surface_complex_refr_index,
2736 const Index& atmosphere_dim,
2737 const Vector& lat_grid,
2738 const Vector& lat_true,
2739 const Vector& lon_true,
2740 const Vector& rtp_pos,
2741 const GriddedField5& complex_n_field,
2742 const Verbosity&) {
2743 // Set expected order of grids
2744 Index gfield_fID = 0;
2745 Index gfield_tID = 1;
2746 Index gfield_compID = 2;
2747 Index gfield_latID = 3;
2748 Index gfield_lonID = 4;
2749
2750 // Basic checks and sizes
2751 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2752 chk_latlon_true(atmosphere_dim, lat_grid, lat_true, lon_true);
2753 chk_rte_pos(atmosphere_dim, rtp_pos);
2754 complex_n_field.checksize_strict();
2755 //
2756 chk_griddedfield_gridname(complex_n_field, gfield_fID, "Frequency");
2757 chk_griddedfield_gridname(complex_n_field, gfield_tID, "Temperature");
2758 chk_griddedfield_gridname(complex_n_field, gfield_compID, "Complex");
2759 chk_griddedfield_gridname(complex_n_field, gfield_latID, "Latitude");
2760 chk_griddedfield_gridname(complex_n_field, gfield_lonID, "Longitude");
2761 //
2762 const Index nf = complex_n_field.data.nshelves();
2763 const Index nt = complex_n_field.data.nbooks();
2764 const Index nn = complex_n_field.data.npages();
2765 const Index nlat = complex_n_field.data.nrows();
2766 const Index nlon = complex_n_field.data.ncols();
2767 //
2768 if (nlat < 2 || nlon < 2) {
2769 ostringstream os;
2770 os << "The data in *complex_refr_index_field* must span a geographical "
2771 << "region. That is,\nthe latitude and longitude grids must have a "
2772 << "length >= 2.";
2773 throw runtime_error(os.str());
2774 }
2775 //
2776 if (nn != 2) {
2777 ostringstream os;
2778 os << "The data in *complex_refr_index_field* must have exactly two "
2779 << "pages. One page each\nfor the real and imaginary part of the "
2780 << "complex refractive index.";
2781 throw runtime_error(os.str());
2782 }
2783
2784 const Vector& GFlat = complex_n_field.get_numeric_grid(gfield_latID);
2785 const Vector& GFlon = complex_n_field.get_numeric_grid(gfield_lonID);
2786
2787 // Determine true geographical position
2788 Vector lat(1), lon(1);
2790 lat[0], lon[0], atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
2791
2792 // Ensure correct coverage of lon grid
2793 Vector lon_shifted;
2794 lon_shiftgrid(lon_shifted, GFlon, lon[0]);
2795
2796 // Check if lat/lon we need are actually covered
2797 chk_if_in_range("rtp_pos.lat", lat[0], GFlat[0], GFlat[nlat - 1]);
2798 chk_if_in_range("rtp_pos.lon", lon[0], lon_shifted[0], lon_shifted[nlon - 1]);
2799
2800 // Size and fills grids of *surface_complex_refr_index*
2801 surface_complex_refr_index.resize(nf, nt, 2);
2802 surface_complex_refr_index.set_grid_name(0, "Frequency");
2803 surface_complex_refr_index.set_grid(
2804 0, complex_n_field.get_numeric_grid(gfield_fID));
2805 surface_complex_refr_index.set_grid_name(1, "Temperature");
2806 surface_complex_refr_index.set_grid(
2807 1, complex_n_field.get_numeric_grid(gfield_tID));
2808 surface_complex_refr_index.set_grid_name(2, "Complex");
2809 surface_complex_refr_index.set_grid(2, {"real", "imaginary"});
2810
2811 // Interpolate in lat and lon
2812 //
2813 GridPos gp_lat, gp_lon;
2814 gridpos(gp_lat, GFlat, lat[0]);
2815 gridpos(gp_lon, lon_shifted, lon[0]);
2816 Vector itw(4);
2817 interpweights(itw, gp_lat, gp_lon);
2818 //
2819 for (Index iv = 0; iv < nf; iv++) {
2820 for (Index it = 0; it < nt; it++) {
2821 surface_complex_refr_index.data(iv, it, 0) = interp(
2822 itw, complex_n_field.data(iv, it, 0, joker, joker), gp_lat, gp_lon);
2823 surface_complex_refr_index.data(iv, it, 1) = interp(
2824 itw, complex_n_field.data(iv, it, 1, joker, joker), gp_lat, gp_lon);
2825 }
2826 }
2827}
2828
2829/* Workspace method: Doxygen documentation will be auto-generated */
2830void surface_reflectivityFromGriddedField6(Tensor3& surface_reflectivity,
2831 const Index& stokes_dim,
2832 const Vector& f_grid,
2833 const Index& atmosphere_dim,
2834 const Vector& lat_grid,
2835 const Vector& lat_true,
2836 const Vector& lon_true,
2837 const Vector& rtp_pos,
2838 const Vector& rtp_los,
2839 const GriddedField6& r_field,
2840 const Verbosity&) {
2841 // Basic checks and sizes
2842 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2843 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
2844 chk_latlon_true(atmosphere_dim, lat_grid, lat_true, lon_true);
2845 chk_rte_pos(atmosphere_dim, rtp_pos);
2846 chk_rte_los(atmosphere_dim, rtp_los);
2847 r_field.checksize_strict();
2848 chk_griddedfield_gridname(r_field, 0, "Frequency");
2849 chk_griddedfield_gridname(r_field, 1, "Stokes element");
2850 chk_griddedfield_gridname(r_field, 2, "Stokes element");
2851 chk_griddedfield_gridname(r_field, 3, "Incidence angle");
2852 chk_griddedfield_gridname(r_field, 4, "Latitude");
2853 chk_griddedfield_gridname(r_field, 5, "Longitude");
2854 //
2855 const Index nf_in = r_field.data.nvitrines();
2856 const Index ns2 = r_field.data.nshelves();
2857 const Index ns1 = r_field.data.nbooks();
2858 const Index nza = r_field.data.npages();
2859 const Index nlat = r_field.data.nrows();
2860 const Index nlon = r_field.data.ncols();
2861 //
2862 if (nlat < 2 || nlon < 2) {
2863 ostringstream os;
2864 os << "The data in *r_field* must span a geographical region. That is,\n"
2865 << "the latitude and longitude grids must have a length >= 2.";
2866 throw runtime_error(os.str());
2867 }
2868 //
2869 if (nza < 2) {
2870 ostringstream os;
2871 os << "The data in *r_field* must span a range of zenith angles. That\n"
2872 << "is the zenith angle grid must have a length >= 2.";
2873 throw runtime_error(os.str());
2874 }
2875 if (ns1 < stokes_dim || ns2 < stokes_dim || ns1 > 4 || ns2 > 4) {
2876 ostringstream os;
2877 os << "The \"Stokes dimensions\" must have a size that is >= "
2878 << "*stokes_dim* (but not exceeding 4).";
2879 throw runtime_error(os.str());
2880 }
2881
2882 // Determine true geographical position
2883 Vector lat(1), lon(1);
2885 lat[0], lon[0], atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
2886
2887 // Ensure correct coverage of lon grid
2888 Vector lon_shifted;
2889 lon_shiftgrid(lon_shifted, r_field.get_numeric_grid(5), lon[0]);
2890
2891 // Interpolate in lat and lon
2892 Tensor4 r_f_za(nf_in, stokes_dim, stokes_dim, nza);
2893 {
2895 "Latitude interpolation", r_field.get_numeric_grid(4), lat[0]);
2896 chk_interpolation_grids("Longitude interpolation", lon_shifted, lon[0]);
2897 GridPos gp_lat, gp_lon;
2898 gridpos(gp_lat, r_field.get_numeric_grid(4), lat[0]);
2899 gridpos(gp_lon, lon_shifted, lon[0]);
2900 Vector itw(4);
2901 interpweights(itw, gp_lat, gp_lon);
2902 for (Index iv = 0; iv < nf_in; iv++) {
2903 for (Index iz = 0; iz < nza; iz++) {
2904 for (Index is1 = 0; is1 < stokes_dim; is1++) {
2905 for (Index is2 = 0; is2 < stokes_dim; is2++) {
2906 r_f_za(iv, is1, is2, iz) =
2907 interp(itw,
2908 r_field.data(iv, is1, is2, iz, joker, joker),
2909 gp_lat,
2910 gp_lon);
2911 }
2912 }
2913 }
2914 }
2915 }
2916
2917 // Interpolate in incidence angle, cubic if possible
2918 Tensor3 r_f(nf_in, stokes_dim, stokes_dim);
2919 Index order = 3;
2920 if (nza < 4) {
2921 order = 1;
2922 }
2923 {
2925 "Incidence angle interpolation", r_field.get_numeric_grid(3), 180 - rtp_los[0]);
2926 const LagrangeInterpolation lag(0, 180 - rtp_los[0], r_field.get_numeric_grid(3), order);
2927 const auto itw = interpweights(lag);
2928 //
2929 for (Index i = 0; i < nf_in; i++) {
2930 for (Index is1 = 0; is1 < stokes_dim; is1++) {
2931 for (Index is2 = 0; is2 < stokes_dim; is2++) {
2932 r_f(i, is1, is2) = interp(r_f_za(i, is1, is2, joker), itw, lag);
2933 }
2934 }
2935 }
2936 }
2937
2938 // Extract or interpolate in frequency
2939 //
2940 if (nf_in == 1) {
2941 surface_reflectivity = r_f;
2942 } else {
2944 "Frequency interpolation", r_field.get_numeric_grid(0), f_grid);
2945 const Index nf_out = f_grid.nelem();
2946 surface_reflectivity.resize(nf_out, stokes_dim, stokes_dim);
2947 //
2948 ArrayOfGridPos gp(nf_out);
2949 Matrix itw(nf_out, 2);
2950 gridpos(gp, r_field.get_numeric_grid(0), f_grid);
2951 interpweights(itw, gp);
2952 for (Index is1 = 0; is1 < stokes_dim; is1++) {
2953 for (Index is2 = 0; is2 < stokes_dim; is2++) {
2954 interp(surface_reflectivity(joker, is1, is2),
2955 itw,
2956 r_f(joker, is1, is2),
2957 gp);
2958 }
2959 }
2960 }
2961}
2962
2963/* Workspace method: Doxygen documentation will be auto-generated */
2965 Vector& surface_scalar_reflectivity,
2966 const Index& stokes_dim,
2967 const Vector& f_grid,
2968 const Index& atmosphere_dim,
2969 const Vector& lat_grid,
2970 const Vector& lat_true,
2971 const Vector& lon_true,
2972 const Vector& rtp_pos,
2973 const Vector& rtp_los,
2974 const GriddedField4& r_field,
2975 const Verbosity&) {
2976 // Basic checks and sizes
2977 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2978 chk_if_in_range("stokes_dim", stokes_dim, 1, 1);
2979 chk_latlon_true(atmosphere_dim, lat_grid, lat_true, lon_true);
2980 chk_rte_pos(atmosphere_dim, rtp_pos);
2981 chk_rte_los(atmosphere_dim, rtp_los);
2982 r_field.checksize_strict();
2983 chk_griddedfield_gridname(r_field, 0, "Frequency");
2984 chk_griddedfield_gridname(r_field, 1, "Incidence angle");
2985 chk_griddedfield_gridname(r_field, 2, "Latitude");
2986 chk_griddedfield_gridname(r_field, 3, "Longitude");
2987 //
2988 const Index nf_in = r_field.data.nbooks();
2989 const Index nza = r_field.data.npages();
2990 const Index nlat = r_field.data.nrows();
2991 const Index nlon = r_field.data.ncols();
2992 //
2993 if (nlat < 2 || nlon < 2) {
2994 ostringstream os;
2995 os << "The data in *r_field* must span a geographical region. That is,\n"
2996 << "the latitude and longitude grids must have a length >= 2.";
2997 throw runtime_error(os.str());
2998 }
2999 //
3000 if (nza < 2) {
3001 ostringstream os;
3002 os << "The data in *r_field* must span a range of zenith angles. That\n"
3003 << "is the zenith angle grid must have a length >= 2.";
3004 throw runtime_error(os.str());
3005 }
3006
3007 // Determine true geographical position
3008 Vector lat(1), lon(1);
3010 lat[0], lon[0], atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
3011
3012 // Ensure correct coverage of lon grid
3013 Vector lon_shifted;
3014 lon_shiftgrid(lon_shifted, r_field.get_numeric_grid(3), lon[0]);
3015
3016 // Interpolate in lat and lon
3017 Matrix r_f_za(nf_in, nza);
3018 {
3020 "Latitude interpolation", r_field.get_numeric_grid(2), lat[0]);
3021 chk_interpolation_grids("Longitude interpolation", lon_shifted, lon[0]);
3022 GridPos gp_lat, gp_lon;
3023 gridpos(gp_lat, r_field.get_numeric_grid(2), lat[0]);
3024 gridpos(gp_lon, lon_shifted, lon[0]);
3025 Vector itw(4);
3026 interpweights(itw, gp_lat, gp_lon);
3027 for (Index iv = 0; iv < nf_in; iv++) {
3028 for (Index iz = 0; iz < nza; iz++) {
3029 r_f_za(iv, iz) =
3030 interp(itw, r_field.data(iv, iz, joker, joker), gp_lat, gp_lon);
3031 }
3032 }
3033 }
3034
3035 // Interpolate in incidence angle, cubic if possible
3036 Vector r_f(nf_in);
3037 Index order = 3;
3038 if (nza < 4) {
3039 order = 1;
3040 }
3041 {
3043 "Incidence angle interpolation", r_field.get_numeric_grid(1), 180 - rtp_los[0]);
3044 const LagrangeInterpolation lag(0, 180 - rtp_los[0], r_field.get_numeric_grid(1), order);
3045 const auto itw=interpweights(lag);
3046 for (Index i = 0; i < nf_in; i++) {
3047 r_f[i] = interp(r_f_za(i, joker), itw, lag);
3048 }
3049 }
3050
3051 // Extract or interpolate in frequency
3052 //
3053 if (nf_in == 1) {
3054 surface_scalar_reflectivity.resize(1);
3055 surface_scalar_reflectivity[0] = r_f[0];
3056 } else {
3058 "Frequency interpolation", r_field.get_numeric_grid(0), f_grid);
3059 const Index nf_out = f_grid.nelem();
3060 surface_scalar_reflectivity.resize(nf_out);
3061 //
3062 ArrayOfGridPos gp(nf_out);
3063 Matrix itw(nf_out, 2);
3064 gridpos(gp, r_field.get_numeric_grid(0), f_grid);
3065 interpweights(itw, gp);
3066 interp(surface_scalar_reflectivity, itw, r_f, gp);
3067 }
3068}
3069
3070/* Workspace method: Doxygen documentation will be auto-generated */
3072 Vector& surface_scalar_reflectivity,
3073 const Tensor4& surface_rmatrix,
3074 const Verbosity&) {
3075 const Index nf = surface_rmatrix.npages();
3076 const Index nlos = surface_rmatrix.nbooks();
3077
3078 surface_scalar_reflectivity.resize(nf);
3079 surface_scalar_reflectivity = 0;
3080
3081 for (Index i = 0; i < nf; i++) {
3082 for (Index l = 0; l < nlos; l++) {
3083 surface_scalar_reflectivity[i] += surface_rmatrix(l, i, 0, 0);
3084 }
3085 }
3086}
3087
3088/* Workspace method: Doxygen documentation will be auto-generated */
3089void InterpSurfaceTypeMask(Index& surface_type,
3090 const Index& atmosphere_dim,
3091 const Vector& lat_grid,
3092 const Vector& lat_true,
3093 const Vector& lon_true,
3094 const Vector& rtp_pos,
3095 const GriddedField2& surface_type_mask,
3096 const Verbosity&) {
3097 // Set expected order of grids
3098 Index gfield_latID = 0;
3099 Index gfield_lonID = 1;
3100
3101 // Basic checks and sizes
3102 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
3103 chk_latlon_true(atmosphere_dim, lat_grid, lat_true, lon_true);
3104 chk_rte_pos(atmosphere_dim, rtp_pos);
3105 surface_type_mask.checksize_strict();
3106 //
3107 chk_griddedfield_gridname(surface_type_mask, gfield_latID, "Latitude");
3108 chk_griddedfield_gridname(surface_type_mask, gfield_lonID, "Longitude");
3109 //
3110 const Index nlat = surface_type_mask.data.nrows();
3111 const Index nlon = surface_type_mask.data.ncols();
3112 //
3113 if (nlat < 2 || nlon < 2) {
3114 ostringstream os;
3115 os << "The data in *surface_type_mask* must span a geographical "
3116 << "region. Accordingly,\nthe latitude and longitude grids must "
3117 << "both have a length >= 2.";
3118 throw runtime_error(os.str());
3119 }
3120
3121 const Vector& GFlat = surface_type_mask.get_numeric_grid(gfield_latID);
3122 const Vector& GFlon = surface_type_mask.get_numeric_grid(gfield_lonID);
3123
3124 // Determine true geographical position
3125 Vector lat(1), lon(1);
3127 lat[0], lon[0], atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
3128
3129 // Ensure correct coverage of lon grid
3130 Vector lon_shifted;
3131 lon_shiftgrid(lon_shifted, GFlon, lon[0]);
3132
3133 // Check if lat/lon we need are actually covered
3134 chk_if_in_range("rtp_pos.lat", lat[0], GFlat[0], GFlat[nlat - 1]);
3135 chk_if_in_range("rtp_pos.lon", lon[0], lon_shifted[0], lon_shifted[nlon - 1]);
3136
3137 // Grid positions
3138 GridPos gp_lat, gp_lon;
3139 gridpos(gp_lat, GFlat, lat[0]);
3140 gridpos(gp_lon, lon_shifted, lon[0]);
3141
3142 // Extract closest point
3143 Index ilat, ilon;
3144 if (gp_lat.fd[0] < 0.5) {
3145 ilat = gp_lat.idx;
3146 } else {
3147 ilat = gp_lat.idx + 1;
3148 }
3149 if (gp_lon.fd[0] < 0.5) {
3150 ilon = gp_lon.idx;
3151 } else {
3152 ilon = gp_lon.idx + 1;
3153 }
3154 surface_type = (Index)round(surface_type_mask.data(ilat, ilon));
3155}
3156
3157/* Workspace method: Doxygen documentation will be auto-generated */
3159 Workspace& ws,
3160 Vector& surface_type_mix,
3161 Numeric& surface_skin_t,
3162 Matrix& surface_los,
3163 Tensor4& surface_rmatrix,
3164 Matrix& surface_emission,
3165 const Vector& f_grid,
3166 const Index& stokes_dim,
3167 const Index& atmosphere_dim,
3168 const Vector& lat_grid,
3169 const Vector& lon_grid,
3170 const Vector& lat_true,
3171 const Vector& lon_true,
3172 const Vector& rtp_pos,
3173 const Vector& rtp_los,
3174 const Vector& refellipsoid,
3175 const GriddedField2& surface_type_mask,
3176 const ArrayOfAgenda& surface_rtprop_agenda_array,
3177 const Numeric& z_sensor,
3178 const Matrix& dlos,
3179 const Vector& dlos_weight_vector,
3180 const Verbosity& verbosity)
3181{
3182 ARTS_USER_ERROR_IF(atmosphere_dim != 3, "This method only works for 3D.");
3183
3184 // Determine satellite sensor_pos and sensor_los
3185 Vector sat_pos = rtp_pos;
3186 Vector sat_los = rtp_los;
3188 sat_los,
3189 atmosphere_dim,
3190 lat_grid,
3191 lon_grid,
3192 refellipsoid,
3193 z_sensor,
3194 0,
3195 verbosity);
3196
3197 // Some sizes
3198 const Index nf = f_grid.nelem();
3199 const Index ntypes = surface_rtprop_agenda_array.nelem();
3200 const Index nlos = dlos.nrows();
3201
3202 // Find ground sample positions
3203 Matrix ground_pos, ground_los, sensor_pos(nlos, 3), sensor_los;
3204 //
3205 sensor_pos(joker, 0) = sat_pos[0];
3206 sensor_pos(joker, 1) = sat_pos[1];
3207 sensor_pos(joker, 2) = sat_pos[2];
3208 losAddLosAndDlos(sensor_los, sat_los, dlos, verbosity);
3209 //
3211 ground_los,
3212 sensor_pos,
3213 sensor_los,
3214 refellipsoid,
3215 lat_grid,
3216 lon_grid,
3217 rtp_pos[0],
3218 verbosity);
3219
3220 // Prepare output variables
3221 surface_type_mix.resize(ntypes);
3222 surface_type_mix = 0.;
3223 surface_skin_t = 0.;
3224 surface_los.resize(1,2);
3225 surface_los = 0.;
3226 surface_rmatrix.resize(1, nf, stokes_dim, stokes_dim);
3227 surface_rmatrix = 0.;
3228 surface_emission.resize(nf, stokes_dim);
3229 surface_emission = 0.;
3230
3231 // Help variables
3232 const Numeric weight_sum = sum(dlos_weight_vector);
3233 Numeric tmp_skin_t;
3234 Vector tmp_type_mix;
3235 Matrix tmp_emission, tmp_los;
3236 Tensor4 tmp_rmatrix;
3237
3238 // Loop los-es to sample
3239 for (Index i=0; i<nlos; ++i) {
3241 tmp_type_mix,
3242 tmp_skin_t,
3243 tmp_los,
3244 tmp_rmatrix,
3245 tmp_emission,
3246 f_grid,
3247 atmosphere_dim,
3248 lat_grid,
3249 lat_true,
3250 lon_true,
3251 Vector{ground_pos(i, joker)},
3252 Vector{ground_los(i, joker)},
3253 surface_type_mask,
3254 surface_rtprop_agenda_array,
3255 verbosity);
3256 ARTS_USER_ERROR_IF(tmp_los.nrows() != 1,
3257 "This method requires that all surface types "
3258 "returns a *surface_los* with one row.");
3259
3260 // Sum up
3261 const Numeric weight = dlos_weight_vector[i] / weight_sum;
3262 tmp_type_mix *= weight;
3263 surface_type_mix += tmp_type_mix;
3264 surface_skin_t += weight * tmp_skin_t;
3265 tmp_los *= weight;
3266 surface_los += tmp_los;
3267 tmp_rmatrix *= weight;
3268 surface_rmatrix += tmp_rmatrix;
3269 tmp_emission *= weight;
3270 surface_emission += tmp_emission;
3271 }
3272}
3273
3274/* Workspace method: Doxygen documentation will be auto-generated */
3276 Vector& surface_type_mix,
3277 Numeric& surface_skin_t,
3278 Matrix& surface_los,
3279 Tensor4& surface_rmatrix,
3280 Matrix& surface_emission,
3281 const Vector& f_grid,
3282 const Vector& rtp_pos,
3283 const Vector& rtp_los,
3284 const ArrayOfAgenda& surface_rtprop_agenda_array,
3285 const Index& surface_type,
3286 const Verbosity&)
3287{
3288 ARTS_USER_ERROR_IF(surface_type < 0 or
3289 surface_type >= surface_rtprop_agenda_array.nelem(),
3290 "Provided surface type index invalid (<0 or too high w.r.t. "
3291 "length of *surface_rtprop_agenda_array*).");
3292
3293 surface_type_mix.resize(surface_rtprop_agenda_array.nelem());
3294 surface_type_mix = 0.0;
3295 surface_type_mix[surface_type] = 1.0;
3296
3298 surface_skin_t,
3299 surface_emission,
3300 surface_los,
3301 surface_rmatrix,
3302 surface_type,
3303 f_grid,
3304 rtp_pos,
3305 rtp_los,
3306 surface_rtprop_agenda_array);
3307}
3308
3309/* Workspace method: Doxygen documentation will be auto-generated */
3311 Vector& surface_type_mix,
3312 Numeric& surface_skin_t,
3313 Matrix& surface_los,
3314 Tensor4& surface_rmatrix,
3315 Matrix& surface_emission,
3316 const Vector& f_grid,
3317 const Index& atmosphere_dim,
3318 const Vector& lat_grid,
3319 const Vector& lat_true,
3320 const Vector& lon_true,
3321 const Vector& rtp_pos,
3322 const Vector& rtp_los,
3323 const GriddedField2& surface_type_mask,
3324 const ArrayOfAgenda& surface_rtprop_agenda_array,
3325 const Verbosity& verbosity)
3326{
3327 Index surface_type;
3328 InterpSurfaceTypeMask(surface_type,
3329 atmosphere_dim,
3330 lat_grid,
3331 lat_true,
3332 lon_true,
3333 rtp_pos,
3334 surface_type_mask,
3335 verbosity);
3336 ARTS_USER_ERROR_IF(surface_type < 0 or
3337 surface_type >= surface_rtprop_agenda_array.nelem(),
3338 "Interpolation gave invalid surface type index (<0 or too "
3339 "high w.f.t. length of *surface_rtprop_agenda_array*).");
3340
3341 surface_type_mix.resize(surface_rtprop_agenda_array.nelem());
3342 surface_type_mix = 0.0;
3343 surface_type_mix[surface_type] = 1.0;
3344
3346 surface_skin_t,
3347 surface_emission,
3348 surface_los,
3349 surface_rmatrix,
3350 surface_type,
3351 f_grid,
3352 rtp_pos,
3353 rtp_los,
3354 surface_rtprop_agenda_array);
3355}
3356
3357/* Workspace method: Doxygen documentation will be auto-generated */
3358void surface_rtpropInterpFreq(Vector& f_grid,
3359 Tensor4& surface_rmatrix,
3360 Matrix& surface_emission,
3361 const Vector& f_new,
3362 const Verbosity&) {
3363 const Index nf = f_grid.nelem();
3364 const Index ns = surface_emission.ncols();
3365 const Index nlos = surface_rmatrix.nbooks();
3366 const Index nnew = f_new.nelem();
3367
3368 // Checks
3369 ARTS_USER_ERROR_IF(surface_emission.nrows() not_eq nf,
3370 "Different number of frequencies in *f_grid* and *surface_emission*.");
3371 ARTS_USER_ERROR_IF(surface_rmatrix.npages() not_eq nf,
3372 "Different number of frequencies in *f_grid* and *surface_rmatrix*.");
3373 ARTS_USER_ERROR_IF(surface_rmatrix.ncols() not_eq ns,
3374 "Different number of Stokes elements in *surface_emission* and *surface_rmatrix*.");
3375
3376 // Set up interpolation
3377 chk_interpolation_grids("Frequency interpolation", f_grid, f_new);
3378 ArrayOfGridPos gp(nnew);
3379 Matrix itw(nnew, 2);
3380 gridpos(gp, f_grid, f_new);
3381 interpweights(itw, gp);
3382
3383 // Interpolate
3384 Tensor4 rmatrix = surface_rmatrix;
3385 Matrix emission = surface_emission;
3386 surface_rmatrix.resize(nlos, nnew, ns, ns);
3387 surface_emission.resize(nnew, ns);
3388 //
3389 for (Index is = 0; is < ns; ++is) {
3390 interp(surface_emission(joker, is), itw, emission(joker, is), gp);
3391 for (Index il = 0; il < nlos; ++il) {
3392 for (Index is2 = 0; is2 < ns; ++is2) {
3393 interp(surface_rmatrix(il, joker, is, is2),
3394 itw,
3395 rmatrix(il, joker, is, is2),
3396 gp);
3397 }
3398 }
3399 }
3400 f_grid = f_new;
3401}
3402
3403/* Workspace method: Doxygen documentation will be auto-generated */
3404void SurfaceBlackbody(Matrix& surface_los,
3405 Tensor4& surface_rmatrix,
3406 ArrayOfTensor4& dsurface_rmatrix_dx,
3407 Matrix& surface_emission,
3408 ArrayOfMatrix& dsurface_emission_dx,
3409 const Index& stokes_dim,
3410 const Index& atmosphere_dim,
3411 const Vector& lat_grid,
3412 const Vector& lon_grid,
3413 const Vector& f_grid,
3414 const Vector& rtp_pos,
3415 const Vector& rtp_los,
3416 const Tensor3& surface_props_data,
3417 const ArrayOfString& surface_props_names,
3418 const ArrayOfString& dsurface_names,
3419 const Index& jacobian_do,
3420 const Verbosity& verbosity) {
3421 // Check surface_data
3422 surface_props_check(atmosphere_dim,
3423 lat_grid,
3424 lon_grid,
3425 surface_props_data,
3426 surface_props_names);
3427
3428 // Interplation grid positions and weights
3429 ArrayOfGridPos gp_lat(1), gp_lon(1);
3430 Matrix itw;
3432 gp_lat[0], gp_lon[0], atmosphere_dim, lat_grid, lon_grid, rtp_pos);
3433 interp_atmsurface_gp2itw(itw, atmosphere_dim, gp_lat, gp_lon);
3434
3435 // Skin temperature
3436 Vector skin_t(1);
3437 surface_props_interp(skin_t,
3438 "Skin temperature",
3439 atmosphere_dim,
3440 gp_lat,
3441 gp_lon,
3442 itw,
3443 surface_props_data,
3444 surface_props_names);
3445
3446 surfaceBlackbody(surface_los,
3447 surface_rmatrix,
3448 surface_emission,
3449 atmosphere_dim,
3450 f_grid,
3451 stokes_dim,
3452 rtp_pos,
3453 rtp_los,
3454 skin_t[0],
3455 verbosity);
3456
3457 surface_rmatrix.resize(1, f_grid.nelem(), stokes_dim, stokes_dim);
3458 surface_rmatrix = 0.0;
3459
3460 // Jacobian part
3461 if (jacobian_do) {
3462 dsurface_check(surface_props_names,
3463 dsurface_names,
3464 dsurface_rmatrix_dx,
3465 dsurface_emission_dx);
3466
3467 Index irq;
3468
3469 // Skin temperature
3470 irq = find_first(dsurface_names, String("Skin temperature"));
3471 if (irq >= 0) {
3472 Matrix dbdt(f_grid.nelem(), 1);
3473 dplanck_dt(dbdt(joker, 0), f_grid, skin_t[0]);
3474 dsurface_emission_dx[irq] = dbdt;
3475 dsurface_rmatrix_dx[irq].resize(surface_rmatrix.nbooks(),
3476 surface_rmatrix.npages(),
3477 surface_rmatrix.nrows(),
3478 surface_rmatrix.ncols());
3479 dsurface_rmatrix_dx[irq] = 0;
3480 }
3481 }
3482}
3483
3484/* Workspace method: Doxygen documentation will be auto-generated */
3485void SurfaceDummy(ArrayOfTensor4& dsurface_rmatrix_dx,
3486 ArrayOfMatrix& dsurface_emission_dx,
3487 const Index& atmosphere_dim,
3488 const Vector& lat_grid,
3489 const Vector& lon_grid,
3490 const Tensor3& surface_props_data,
3491 const ArrayOfString& surface_props_names,
3492 const ArrayOfString& dsurface_names,
3493 const Index& jacobian_do,
3494 const Verbosity&) {
3495 if (surface_props_names.nelem())
3496 throw runtime_error(
3497 "When calling this method, *surface_props_names* should be empty.");
3498
3499 surface_props_check(atmosphere_dim,
3500 lat_grid,
3501 lon_grid,
3502 surface_props_data,
3503 surface_props_names);
3504
3505 if (jacobian_do) {
3506 dsurface_check(surface_props_names,
3507 dsurface_names,
3508 dsurface_rmatrix_dx,
3509 dsurface_emission_dx);
3510 }
3511}
3512
3513/* Workspace method: Doxygen documentation will be auto-generated */
3514void SurfaceFastem(Matrix& surface_los,
3515 Tensor4& surface_rmatrix,
3516 ArrayOfTensor4& dsurface_rmatrix_dx,
3517 Matrix& surface_emission,
3518 ArrayOfMatrix& dsurface_emission_dx,
3519 const Index& stokes_dim,
3520 const Index& atmosphere_dim,
3521 const Vector& lat_grid,
3522 const Vector& lon_grid,
3523 const Vector& f_grid,
3524 const Vector& rtp_pos,
3525 const Vector& rtp_los,
3526 const Tensor3& surface_props_data,
3527 const ArrayOfString& surface_props_names,
3528 const ArrayOfString& dsurface_names,
3529 const Index& jacobian_do,
3530 const Vector& transmittance,
3531 const Index& fastem_version,
3532 const Verbosity& verbosity) {
3533 // Check surface_data
3534 surface_props_check(atmosphere_dim,
3535 lat_grid,
3536 lon_grid,
3537 surface_props_data,
3538 surface_props_names);
3539
3540 // Interplation grid positions and weights
3541 ArrayOfGridPos gp_lat(1), gp_lon(1);
3542 Matrix itw;
3544 gp_lat[0], gp_lon[0], atmosphere_dim, lat_grid, lon_grid, rtp_pos);
3545 interp_atmsurface_gp2itw(itw, atmosphere_dim, gp_lat, gp_lon);
3546
3547 // Skin temperature
3548 Vector skin_t(1);
3549 surface_props_interp(skin_t,
3550 "Water skin temperature",
3551 atmosphere_dim,
3552 gp_lat,
3553 gp_lon,
3554 itw,
3555 surface_props_data,
3556 surface_props_names);
3557
3558 // Wind speed
3559 Vector wind_speed(1);
3560 surface_props_interp(wind_speed,
3561 "Wind speed",
3562 atmosphere_dim,
3563 gp_lat,
3564 gp_lon,
3565 itw,
3566 surface_props_data,
3567 surface_props_names);
3568
3569 // Wind direction
3570 Vector wind_direction(1);
3571 surface_props_interp(wind_direction,
3572 "Wind direction",
3573 atmosphere_dim,
3574 gp_lat,
3575 gp_lon,
3576 itw,
3577 surface_props_data,
3578 surface_props_names);
3579
3580 // Salinity
3581 Vector salinity(1);
3582 surface_props_interp(salinity,
3583 "Salinity",
3584 atmosphere_dim,
3585 gp_lat,
3586 gp_lon,
3587 itw,
3588 surface_props_data,
3589 surface_props_names);
3590
3591 // Call FASTEM
3592 surfaceFastem(surface_los,
3593 surface_rmatrix,
3594 surface_emission,
3595 atmosphere_dim,
3596 stokes_dim,
3597 f_grid,
3598 rtp_pos,
3599 rtp_los,
3600 skin_t[0],
3601 salinity[0],
3602 wind_speed[0],
3603 wind_direction[0],
3604 transmittance,
3605 fastem_version,
3606 verbosity);
3607
3608 // Jacobian part
3609 if (jacobian_do) {
3610 dsurface_check(surface_props_names,
3611 dsurface_names,
3612 dsurface_rmatrix_dx,
3613 dsurface_emission_dx);
3614
3615 Index irq;
3616
3617 // Skin temperature
3618 irq = find_first(dsurface_names, String("Water skin temperature"));
3619 if (irq >= 0) {
3620 const Numeric dd = 0.1;
3621 Matrix surface_los2;
3622 surfaceFastem(surface_los2,
3623 dsurface_rmatrix_dx[irq],
3624 dsurface_emission_dx[irq],
3625 atmosphere_dim,
3626 stokes_dim,
3627 f_grid,
3628 rtp_pos,
3629 rtp_los,
3630 skin_t[0] + dd,
3631 salinity[0],
3632 wind_speed[0],
3633 wind_direction[0],
3634 transmittance,
3635 fastem_version,
3636 verbosity);
3637 //
3638 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
3639 dsurface_rmatrix_dx[irq] /= dd;
3640 dsurface_emission_dx[irq] -= surface_emission;
3641 dsurface_emission_dx[irq] /= dd;
3642 }
3643
3644 // Wind speed
3645 irq = find_first(dsurface_names, String("Wind speed"));
3646 if (irq >= 0) {
3647 const Numeric dd = 0.1;
3648 Matrix surface_los2;
3649 surfaceFastem(surface_los2,
3650 dsurface_rmatrix_dx[irq],
3651 dsurface_emission_dx[irq],
3652 atmosphere_dim,
3653 stokes_dim,
3654 f_grid,
3655 rtp_pos,
3656 rtp_los,
3657 skin_t[0],
3658 salinity[0],
3659 wind_speed[0] + dd,
3660 wind_direction[0],
3661 transmittance,
3662 fastem_version,
3663 verbosity);
3664 //
3665 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
3666 dsurface_rmatrix_dx[irq] /= dd;
3667 dsurface_emission_dx[irq] -= surface_emission;
3668 dsurface_emission_dx[irq] /= dd;
3669 }
3670
3671 // Wind direction
3672 irq = find_first(dsurface_names, String("Wind direction"));
3673 if (irq >= 0) {
3674 const Numeric dd = 1;
3675 Matrix surface_los2;
3676 surfaceFastem(surface_los2,
3677 dsurface_rmatrix_dx[irq],
3678 dsurface_emission_dx[irq],
3679 atmosphere_dim,
3680 stokes_dim,
3681 f_grid,
3682 rtp_pos,
3683 rtp_los,
3684 skin_t[0],
3685 salinity[0],
3686 wind_speed[0],
3687 wind_direction[0] + dd,
3688 transmittance,
3689 fastem_version,
3690 verbosity);
3691 //
3692 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
3693 dsurface_rmatrix_dx[irq] /= dd;
3694 dsurface_emission_dx[irq] -= surface_emission;
3695 dsurface_emission_dx[irq] /= dd;
3696 }
3697
3698 // Salinity
3699 irq = find_first(dsurface_names, String("Salinity"));
3700 if (irq >= 0) {
3701 const Numeric dd = 0.0005;
3702 Matrix surface_los2;
3703 surfaceFastem(surface_los2,
3704 dsurface_rmatrix_dx[irq],
3705 dsurface_emission_dx[irq],
3706 atmosphere_dim,
3707 stokes_dim,
3708 f_grid,
3709 rtp_pos,
3710 rtp_los,
3711 skin_t[0],
3712 salinity[0] + dd,
3713 wind_speed[0],
3714 wind_direction[0],
3715 transmittance,
3716 fastem_version,
3717 verbosity);
3718 //
3719 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
3720 dsurface_rmatrix_dx[irq] /= dd;
3721 dsurface_emission_dx[irq] -= surface_emission;
3722 dsurface_emission_dx[irq] /= dd;
3723 }
3724 }
3725}
3726
3727
3728/* Workspace method: Doxygen documentation will be auto-generated */
3729void SurfaceFlatScalarReflectivity(Matrix& surface_los,
3730 Tensor4& surface_rmatrix,
3731 ArrayOfTensor4& dsurface_rmatrix_dx,
3732 Matrix& surface_emission,
3733 ArrayOfMatrix& dsurface_emission_dx,
3734 const Index& stokes_dim,
3735 const Index& atmosphere_dim,
3736 const Vector& lat_grid,
3737 const Vector& lon_grid,
3738 const Vector& f_grid,
3739 const Vector& rtp_pos,
3740 const Vector& rtp_los,
3741 const Vector& specular_los,
3742 const Tensor3& surface_props_data,
3743 const ArrayOfString& surface_props_names,
3744 const ArrayOfString& dsurface_names,
3745 const Index& jacobian_do,
3746 const Vector& f_reflectivities,
3747 const Verbosity& verbosity) {
3748 // Check surface_data
3749 surface_props_check(atmosphere_dim,
3750 lat_grid,
3751 lon_grid,
3752 surface_props_data,
3753 surface_props_names);
3754
3755 // Check of GINs
3756 chk_if_increasing("GIN *f_reflectivities*", f_reflectivities);
3757
3758 // Sizes
3759 const Index nr = f_reflectivities.nelem();
3760 const Index nf = f_grid.nelem();
3761
3762 // Interplation grid positions and weights
3763 ArrayOfGridPos gp_lat(1), gp_lon(1);
3764 Matrix itw;
3766 gp_lat[0], gp_lon[0], atmosphere_dim, lat_grid, lon_grid, rtp_pos);
3767 interp_atmsurface_gp2itw(itw, atmosphere_dim, gp_lat, gp_lon);
3768
3769 // Skin temperature
3770 Vector surface_skin_t(1);
3771 surface_props_interp(surface_skin_t,
3772 "Skin temperature",
3773 atmosphere_dim,
3774 gp_lat,
3775 gp_lon,
3776 itw,
3777 surface_props_data,
3778 surface_props_names);
3779
3780 // Reflectivities
3781 Vector reflectivities(nr);
3782 //
3783 for (Index i=0; i<nr; i++) {
3784 Vector rv(1);
3785 ostringstream sstr;
3786 sstr << "Scalar reflectivity " << i;
3788 sstr.str(),
3789 atmosphere_dim,
3790 gp_lat,
3791 gp_lon,
3792 itw,
3793 surface_props_data,
3794 surface_props_names);
3795 reflectivities[i] = rv[0];
3796 }
3797
3798 // Create surface_scalar_reflectivity
3799 Vector surface_scalar_reflectivity(nf);
3800 ArrayOfGridPos gp(nf);
3801 if (nr==1) {
3802 gp4length1grid(gp);
3803 } else {
3804 gridpos(gp, f_reflectivities, f_grid, 1e99);
3806 }
3807 Matrix itw2(nf, 2);
3808 interpweights(itw2, gp);
3809 interp( surface_scalar_reflectivity, itw2, reflectivities, gp);
3810
3811 // Call non-Jacobian method
3813 surface_rmatrix,
3814 surface_emission,
3815 f_grid,
3816 stokes_dim,
3817 atmosphere_dim,
3818 rtp_pos,
3819 rtp_los,
3820 specular_los,
3821 surface_skin_t[0],
3822 surface_scalar_reflectivity,
3823 verbosity);
3824
3825 // Jacobian part
3826 if (jacobian_do) {
3827 dsurface_check(surface_props_names,
3828 dsurface_names,
3829 dsurface_rmatrix_dx,
3830 dsurface_emission_dx);
3831
3832 Index irq;
3833
3834 // Skin temperature
3835 irq = find_first(dsurface_names, String("Skin temperature"));
3836 if (irq >= 0) {
3837 const Numeric dd = 0.1;
3838 Matrix surface_los2;
3839 surfaceFlatScalarReflectivity(surface_los2,
3840 dsurface_rmatrix_dx[irq],
3841 dsurface_emission_dx[irq],
3842 f_grid,
3843 stokes_dim,
3844 atmosphere_dim,
3845 rtp_pos,
3846 rtp_los,
3847 specular_los,
3848 surface_skin_t[0]+dd,
3849 surface_scalar_reflectivity,
3850 verbosity);
3851 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
3852 dsurface_rmatrix_dx[irq] /= dd;
3853 dsurface_emission_dx[irq] -= surface_emission;
3854 dsurface_emission_dx[irq] /= dd;
3855 }
3856
3857 // Reflectivities
3858 for (Index i=0; i<nr; i++) {
3859 ostringstream sstr;
3860 sstr << "Scalar reflectivity " << i;
3861 irq = find_first(dsurface_names, String(sstr.str()));
3862 if (irq >= 0) {
3863 const Numeric dd = 1e-4;
3864 Vector rpert = reflectivities;
3865 rpert[i] += dd;
3866 interp( surface_scalar_reflectivity, itw2, rpert, gp);
3867 //
3868 Matrix surface_los2;
3869 surfaceFlatScalarReflectivity(surface_los2,
3870 dsurface_rmatrix_dx[irq],
3871 dsurface_emission_dx[irq],
3872 f_grid,
3873 stokes_dim,
3874 atmosphere_dim,
3875 rtp_pos,
3876 rtp_los,
3877 specular_los,
3878 surface_skin_t[0],
3879 surface_scalar_reflectivity,
3880 verbosity);
3881 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
3882 dsurface_rmatrix_dx[irq] /= dd;
3883 dsurface_emission_dx[irq] -= surface_emission;
3884 dsurface_emission_dx[irq] /= dd;
3885 }
3886 }
3887 }
3888}
3889
3890/* Workspace method: Doxygen documentation will be auto-generated */
3891void SurfaceTessem(Matrix& surface_los,
3892 Tensor4& surface_rmatrix,
3893 ArrayOfTensor4& dsurface_rmatrix_dx,
3894 Matrix& surface_emission,
3895 ArrayOfMatrix& dsurface_emission_dx,
3896 const Index& stokes_dim,
3897 const Index& atmosphere_dim,
3898 const Vector& lat_grid,
3899 const Vector& lon_grid,
3900 const Vector& f_grid,
3901 const Vector& rtp_pos,
3902 const Vector& rtp_los,
3903 const TessemNN& net_h,
3904 const TessemNN& net_v,
3905 const Tensor3& surface_props_data,
3906 const ArrayOfString& surface_props_names,
3907 const ArrayOfString& dsurface_names,
3908 const Index& jacobian_do,
3909 const Verbosity& verbosity) {
3910 // Check surface_data
3911 surface_props_check(atmosphere_dim,
3912 lat_grid,
3913 lon_grid,
3914 surface_props_data,
3915 surface_props_names);
3916
3917 // Interplation grid positions and weights
3918 ArrayOfGridPos gp_lat(1), gp_lon(1);
3919 Matrix itw;
3921 gp_lat[0], gp_lon[0], atmosphere_dim, lat_grid, lon_grid, rtp_pos);
3922 interp_atmsurface_gp2itw(itw, atmosphere_dim, gp_lat, gp_lon);
3923
3924 // Skin temperature
3925 Vector skin_t(1);
3926 surface_props_interp(skin_t,
3927 "Water skin temperature",
3928 atmosphere_dim,
3929 gp_lat,
3930 gp_lon,
3931 itw,
3932 surface_props_data,
3933 surface_props_names);
3934
3935 // Wind speed
3936 Vector wind_speed(1);
3937 surface_props_interp(wind_speed,
3938 "Wind speed",
3939 atmosphere_dim,
3940 gp_lat,
3941 gp_lon,
3942 itw,
3943 surface_props_data,
3944 surface_props_names);
3945
3946 // Salinity
3947 Vector salinity(1);
3948 surface_props_interp(salinity,
3949 "Salinity",
3950 atmosphere_dim,
3951 gp_lat,
3952 gp_lon,
3953 itw,
3954 surface_props_data,
3955 surface_props_names);
3956
3957 // Call TESSEM
3958 surfaceTessem(surface_los,
3959 surface_rmatrix,
3960 surface_emission,
3961 atmosphere_dim,
3962 stokes_dim,
3963 f_grid,
3964 rtp_pos,
3965 rtp_los,
3966 skin_t[0],
3967 net_h,
3968 net_v,
3969 salinity[0],
3970 wind_speed[0],
3971 verbosity);
3972
3973 // Jacobian part
3974 if (jacobian_do) {
3975 dsurface_check(surface_props_names,
3976 dsurface_names,
3977 dsurface_rmatrix_dx,
3978 dsurface_emission_dx);
3979
3980 Index irq;
3981
3982 // Skin temperature
3983 irq = find_first(dsurface_names, String("Water skin temperature"));
3984 if (irq >= 0) {
3985 const Numeric dd = 0.1;
3986 Matrix surface_los2;
3987 surfaceTessem(surface_los2,
3988 dsurface_rmatrix_dx[irq],
3989 dsurface_emission_dx[irq],
3990 atmosphere_dim,
3991 stokes_dim,
3992 f_grid,
3993 rtp_pos,
3994 rtp_los,
3995 skin_t[0] + dd,
3996 net_h,
3997 net_v,
3998 salinity[0],
3999 wind_speed[0],
4000 verbosity);
4001 //
4002 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
4003 dsurface_rmatrix_dx[irq] /= dd;
4004 dsurface_emission_dx[irq] -= surface_emission;
4005 dsurface_emission_dx[irq] /= dd;
4006 }
4007
4008 // Wind speed
4009 irq = find_first(dsurface_names, String("Wind speed"));
4010 if (irq >= 0) {
4011 const Numeric dd = 0.1;
4012 Matrix surface_los2;
4013 surfaceTessem(surface_los2,
4014 dsurface_rmatrix_dx[irq],
4015 dsurface_emission_dx[irq],
4016 atmosphere_dim,
4017 stokes_dim,
4018 f_grid,
4019 rtp_pos,
4020 rtp_los,
4021 skin_t[0],
4022 net_h,
4023 net_v,
4024 salinity[0],
4025 wind_speed[0] + dd,
4026 verbosity);
4027 //
4028 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
4029 dsurface_rmatrix_dx[irq] /= dd;
4030 dsurface_emission_dx[irq] -= surface_emission;
4031 dsurface_emission_dx[irq] /= dd;
4032 }
4033
4034 // Salinity
4035 irq = find_first(dsurface_names, String("Salinity"));
4036 if (irq >= 0) {
4037 const Numeric dd = 0.0005;
4038 Matrix surface_los2;
4039 surfaceTessem(surface_los2,
4040 dsurface_rmatrix_dx[irq],
4041 dsurface_emission_dx[irq],
4042 atmosphere_dim,
4043 stokes_dim,
4044 f_grid,
4045 rtp_pos,
4046 rtp_los,
4047 skin_t[0],
4048 net_h,
4049 net_v,
4050 salinity[0] + dd,
4051 wind_speed[0],
4052 verbosity);
4053 //
4054 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
4055 dsurface_rmatrix_dx[irq] /= dd;
4056 dsurface_emission_dx[irq] -= surface_emission;
4057 dsurface_emission_dx[irq] /= dd;
4058 }
4059 }
4060}
4061
4062/* Workspace method: Doxygen documentation will be auto-generated */
4063void transmittanceFromIy_aux(Vector& transmittance,
4064 const ArrayOfString& iy_aux_vars,
4065 const ArrayOfMatrix& iy_aux,
4066 const Verbosity&) {
4067 Index ihit = -1;
4068
4069 for (Index i = 0; i < iy_aux_vars.nelem(); i++) {
4070 if (iy_aux_vars[i] == "Optical depth") {
4071 ihit = i;
4072 break;
4073 }
4074 }
4075
4076 if (ihit < 0)
4077 throw runtime_error("No element in *iy_aux* holds optical depths.");
4078
4079 const Index n = iy_aux[ihit].nrows();
4080
4081 transmittance.resize(n);
4082
4083 for (Index i = 0; i < n; i++) {
4084 transmittance[i] = exp(-iy_aux[ihit](i, 0));
4085 }
4086}
This file contains the definition of Array.
Index find_first(const Array< base > &x, const base &w)
Find first occurance.
Definition: array.h:173
base max(const Array< base > &x)
Max function.
Definition: array.h:128
base min(const Array< base > &x)
Min function.
Definition: array.h:144
The global header file for ARTS.
Constants of physical expressions as constexpr.
Common ARTS conversions.
void iy_main_agendaExecute(Workspace &ws, Matrix &iy, ArrayOfMatrix &iy_aux, Ppath &ppath, ArrayOfTensor3 &diy_dx, Vector &geo_pos, const Index iy_agenda_call1, const Tensor3 &iy_transmittance, const ArrayOfString &iy_aux_vars, const Index iy_id, const String &iy_unit, const Index cloudbox_on, const Index jacobian_do, const Vector &f_grid, const EnergyLevelMap &nlte_field, const Vector &rte_pos, const Vector &rte_los, const Vector &rte_pos2, const Agenda &input_agenda)
Definition: auto_md.cc:25309
void surface_rtprop_agenda_arrayExecute(Workspace &ws, Numeric &surface_skin_t, Matrix &surface_emission, Matrix &surface_los, Tensor4 &surface_rmatrix, const Index agenda_array_index, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const ArrayOfAgenda &input_agenda_array)
Definition: auto_md.cc:26020
void surface_rtprop_agendaExecute(Workspace &ws, Numeric &surface_skin_t, Matrix &surface_emission, Matrix &surface_los, Tensor4 &surface_rmatrix, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Agenda &input_agenda)
Definition: auto_md.cc:25978
void chk_atm_surface(const String &x_name, const Matrix &x, const Index &dim, ConstVectorView lat_grid, ConstVectorView lon_grid)
chk_atm_surface
void chk_atm_grids(const Index &dim, ConstVectorView p_grid, ConstVectorView lat_grid, ConstVectorView lon_grid)
chk_atm_grids
void chk_interpolation_grids(const String &which_interpolation, ConstVectorView old_grid, ConstVectorView new_grid, const Index order, const Numeric &extpolfac, const bool islog)
Check interpolation grids.
Definition: check_input.cc:887
void chk_rte_pos(const Index &atmosphere_dim, ConstVectorView rte_pos, const bool &is_rte_pos2)
chk_rte_pos
void chk_rte_los(const Index &atmosphere_dim, ConstVectorView rte_los)
chk_rte_los
void chk_size(const String &x_name, ConstVectorView x, const Index &c)
Runtime check for size of Vector.
Definition: check_input.cc:402
void chk_griddedfield_gridname(const GriddedField &gf, const Index gridindex, const String &gridname)
Check name of grid in GriddedField.
void chk_latlon_true(const Index &atmosphere_dim, ConstVectorView lat_grid, ConstVectorView lat_true, ConstVectorView lon_true)
chk_latlon_true
void chk_if_in_range(const String &x_name, const Index &x, const Index &x_low, const Index &x_high)
chk_if_in_range
Definition: check_input.cc:67
void chk_if_in_range_exclude(const String &x_name, const Numeric &x, const Numeric &x_low, const Numeric &x_high)
chk_if_in_range_exclude
Definition: check_input.cc:211
void chk_if_increasing(const String &x_name, const ArrayOfIndex &x)
chk_if_increasing
Definition: check_input.cc:92
void chk_if_in_range_exclude_high(const String &x_name, const Numeric &x, const Numeric &x_low, const Numeric &x_high)
chk_if_in_range_exclude_high
Definition: check_input.cc:186
void chk_vector_length(const String &x_name, ConstVectorView x, const Index &l)
chk_vector_length
Definition: check_input.cc:238
void chk_not_negative(const String &x_name, const Numeric &x)
chk_not_negative
Definition: check_input.cc:115
The Agenda class.
Definition: agenda_class.h:52
String name() const
Agenda name.
Index nelem() const ARTS_NOEXCEPT
Definition: array.h:75
void checksize_strict() const final
Strict consistency check.
void resize(const GriddedField3 &gf)
Make this GriddedField3 the same size as the given one.
void checksize_strict() const final
Strict consistency check.
void checksize_strict() const final
Strict consistency check.
void checksize_strict() const final
Strict consistency check.
void set_grid_name(Index i, const String &s)
Set grid name.
void set_grid(Index i, const Vector &g)
Set a numeric grid.
const Vector & get_numeric_grid(Index i) const
Get a numeric grid.
A telsem atlas.
Definition: telsem.h:42
Index calc_cellnum(Numeric lat, Numeric lon) const
Definition: telsem.cc:124
bool contains(Index cellnumber) const
Definition: telsem.h:68
Index get_class2(Index cellnumber) const
Definition: telsem.h:103
Index calc_cellnum_nearest_neighbor(Numeric lat, Numeric lon) const
Definition: telsem.cc:152
Vector get_emis_h(Index cellnum) const
Definition: telsem.h:142
std::pair< Numeric, Numeric > get_coordinates(Index cellnum) const
Definition: telsem.cc:207
Index get_class1(Index cellnumber) const
Definition: telsem.h:85
std::pair< Numeric, Numeric > emis_interp(Numeric theta, Numeric freq, Index class1, Index class2, const ConstVectorView &ev, const ConstVectorView &eh) const
Definition: telsem.cc:269
Vector get_emis_v(Index i) const
Definition: telsem.h:120
Workspace class.
Definition: workspace_ng.h:36
void mult(MatrixView C, ConstMatrixView A, const Block &B)
#define ARTS_ASSERT(condition,...)
Definition: debug.h:84
#define ARTS_USER_ERROR_IF(condition,...)
Definition: debug.h:135
void fastem(Vector &emissivity, Vector &reflectivity, const Numeric frequency, const Numeric za, const Numeric temperature, const Numeric salinity, const Numeric wind_speed, const Numeric transmittance, const Numeric rel_azimuth, const Index fastem_version)
Calculate the surface emissivity using FASTEM.
Definition: fastem.cc:89
This file contains functions that are adapted from FASTEM code which is used to calculate surface emi...
Header file for functions related to gas scattering.
Numeric sphdist(const Numeric &lat1, const Numeric &lon1, const Numeric &lat2, const Numeric &lon2)
sphdist
Definition: geodetic.cc:1318
void lon_shiftgrid(Vector &longrid_out, ConstVectorView longrid_in, const Numeric lon)
lon_shiftgrid
Definition: geodetic.cc:1391
Numeric refell2d(ConstVectorView refellipsoid, ConstVectorView lat_grid, const GridPos gp)
refell2d
Definition: geodetic.cc:1287
void gridpos(ArrayOfGridPos &gp, ConstVectorView old_grid, ConstVectorView new_grid, const Numeric &extpolfac)
Set up a grid position Array.
void interpweights(VectorView itw, const GridPos &tc)
Red 1D interpolation weights.
void gp4length1grid(ArrayOfGridPos &gp)
Grid position matching a grid of length 1.
Numeric interp(ConstVectorView itw, ConstVectorView a, const GridPos &tc)
Red 1D Interpolate.
Header file for interpolation.cc.
void diy_from_pos_to_rgrids(Tensor3View diy_dx, const RetrievalQuantity &jacobian_quantity, ConstMatrixView diy_dpos, const Index &atmosphere_dim, ConstVectorView rtp_pos)
diy_from_pos_to_rgrids
Definition: jacobian.cc:458
void jacobian_type_extrapol(ArrayOfGridPos &gp)
Adopts grid positions to extrapolation used for jacobians.
Definition: jacobian.cc:819
void AngularGridsSetFluxCalc(Vector &za_grid, Vector &aa_grid, Vector &za_grid_weights, const Index &N_za_grid, const Index &N_aa_grid, const String &za_grid_type, const Verbosity &)
WORKSPACE METHOD: AngularGridsSetFluxCalc.
Definition: m_fluxes.cc:42
void IntersectionGeometricalWithAltitude(Matrix &pos, Matrix &los, const Matrix &sensor_pos, const Matrix &sensor_los, const Vector &refellipsoid, const Vector &lat_grid, const Vector &lon_grid, const Numeric &altitude, const Verbosity &)
WORKSPACE METHOD: IntersectionGeometricalWithAltitude.
Definition: m_geodetic.cc:297
void losAddLosAndDlos(Matrix &new_los, const Vector &ref_los, const Matrix &dlos, const Verbosity &)
WORKSPACE METHOD: losAddLosAndDlos.
Definition: m_ppath.cc:291
void rte_pos_losBackwardToAltitude(Vector &rte_pos, Vector &rte_los, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Vector &refellipsoid, const Numeric &altitude, const Index &los_is_reversed, const Verbosity &verbosity)
WORKSPACE METHOD: rte_pos_losBackwardToAltitude.
Definition: m_ppath.cc:2121
void surfaceTelsem(Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Index &atmosphere_dim, const Index &stokes_dim, const Vector &f_grid, const Vector &lat_grid, const Vector &lat_true, const Vector &lon_true, const Vector &rtp_pos, const Vector &rtp_los, const Vector &specular_los, const Numeric &surface_skin_t, const TelsemAtlas &atlas, const Numeric &r_min, const Numeric &r_max, const Numeric &d_max, const Verbosity &verbosity)
WORKSPACE METHOD: surfaceTelsem.
Definition: m_surface.cc:2146
void surface_scalar_reflectivityFromGriddedField4(Vector &surface_scalar_reflectivity, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lat_true, const Vector &lon_true, const Vector &rtp_pos, const Vector &rtp_los, const GriddedField4 &r_field, const Verbosity &)
WORKSPACE METHOD: surface_scalar_reflectivityFromGriddedField4.
Definition: m_surface.cc:2964
void surface_complex_refr_indexFromGriddedField5(GriddedField3 &surface_complex_refr_index, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lat_true, const Vector &lon_true, const Vector &rtp_pos, const GriddedField5 &complex_n_field, const Verbosity &)
WORKSPACE METHOD: surface_complex_refr_indexFromGriddedField5.
Definition: m_surface.cc:2734
void InterpSurfaceFieldToPosition(Numeric &outvalue, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Vector &rtp_pos, const Matrix &z_surface, const Matrix &field, const Verbosity &verbosity)
WORKSPACE METHOD: InterpSurfaceFieldToPosition.
Definition: m_surface.cc:187
void surfaceTessem(Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Index &atmosphere_dim, const Index &stokes_dim, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Numeric &surface_skin_t, const TessemNN &net_h, const TessemNN &net_v, const Numeric &salinity, const Numeric &wind_speed, const Verbosity &verbosity)
WORKSPACE METHOD: surfaceTessem.
Definition: m_surface.cc:2253
constexpr Numeric EARTH_RADIUS
Definition: m_surface.cc:45
void specular_losCalc(Vector &specular_los, Vector &surface_normal, const Vector &rtp_pos, const Vector &rtp_los, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Vector &refellipsoid, const Matrix &z_surface, const Index &ignore_surface_slope, const Verbosity &verbosity)
WORKSPACE METHOD: specular_losCalc.
Definition: m_surface.cc:1828
constexpr Numeric DEG2RAD
Definition: m_surface.cc:46
void surfaceFastem(Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Index &atmosphere_dim, const Index &stokes_dim, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Numeric &surface_skin_t, const Numeric &salinity, const Numeric &wind_speed, const Numeric &wind_direction, const Vector &transmittance, const Index &fastem_version, const Verbosity &verbosity)
WORKSPACE METHOD: surfaceFastem.
Definition: m_surface.cc:1976
void surfaceMapToLinearPolarisation(Matrix &surface_emission, Tensor4 &surface_rmatrix, const Index &stokes_dim, const Numeric &pol_angle, const Verbosity &)
WORKSPACE METHOD: surfaceMapToLinearPolarisation.
Definition: m_surface.cc:2085
void surfaceFlatRefractiveIndex(Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Vector &f_grid, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &rtp_pos, const Vector &rtp_los, const Vector &specular_los, const Numeric &surface_skin_t, const GriddedField3 &surface_complex_refr_index, const Verbosity &verbosity)
WORKSPACE METHOD: surfaceFlatRefractiveIndex.
Definition: m_surface.cc:2332
void surface_rtpropFromTypesNearest(Workspace &ws, Vector &surface_type_mix, Numeric &surface_skin_t, Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Vector &f_grid, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lat_true, const Vector &lon_true, const Vector &rtp_pos, const Vector &rtp_los, const GriddedField2 &surface_type_mask, const ArrayOfAgenda &surface_rtprop_agenda_array, const Verbosity &verbosity)
WORKSPACE METHOD: surface_rtpropFromTypesNearest.
Definition: m_surface.cc:3310
void SurfaceBlackbody(Matrix &surface_los, Tensor4 &surface_rmatrix, ArrayOfTensor4 &dsurface_rmatrix_dx, Matrix &surface_emission, ArrayOfMatrix &dsurface_emission_dx, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Tensor3 &surface_props_data, const ArrayOfString &surface_props_names, const ArrayOfString &dsurface_names, const Index &jacobian_do, const Verbosity &verbosity)
WORKSPACE METHOD: SurfaceBlackbody.
Definition: m_surface.cc:3404
void surface_rtpropInterpFreq(Vector &f_grid, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Vector &f_new, const Verbosity &)
WORKSPACE METHOD: surface_rtpropInterpFreq.
Definition: m_surface.cc:3358
void specular_losCalcNoTopography(Vector &specular_los, Vector &surface_normal, const Vector &rtp_pos, const Vector &rtp_los, const Index &atmosphere_dim, const Verbosity &)
WORKSPACE METHOD: specular_losCalcNoTopography.
Definition: m_surface.cc:1791
void iySurfaceFlatRefractiveIndex(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, ArrayOfTensor4 &dsurface_rmatrix_dx, ArrayOfMatrix &dsurface_emission_dx, const Tensor3 &iy_transmittance, const Index &iy_id, const Index &jacobian_do, const Index &suns_do, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, const Vector &f_grid, const Vector &lat_grid, const Vector &lon_grid, const Matrix &z_surface, const Vector &refellipsoid, const Vector &rtp_pos, const Vector &rtp_los, const Vector &rte_pos2, const String &iy_unit, const GriddedField3 &surface_complex_refr_index, const Tensor3 &surface_props_data, const ArrayOfString &surface_props_names, const ArrayOfString &dsurface_names, const ArrayOfRetrievalQuantity &jacobian_quantities, const Agenda &iy_main_agenda, const Verbosity &verbosity)
WORKSPACE METHOD: iySurfaceFlatRefractiveIndex.
Definition: m_surface.cc:666
void iySurfaceLambertian(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, const Tensor3 &iy_transmittance, const Index &iy_id, const Index &jacobian_do, const Index &suns_do, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, const Vector &f_grid, const Vector &lat_grid, const Vector &lon_grid, const Matrix &z_surface, const Vector &refellipsoid, const Vector &rtp_pos, const Vector &rtp_los, const Vector &rte_pos2, const String &iy_unit, const Vector &surface_scalar_reflectivity, const Tensor3 &surface_props_data, const ArrayOfString &surface_props_names, const ArrayOfString &dsurface_names, const ArrayOfRetrievalQuantity &jacobian_quantities, const Agenda &iy_main_agenda, const Index &N_za, const Index &N_aa, const Verbosity &verbosity)
WORKSPACE METHOD: iySurfaceLambertian.
Definition: m_surface.cc:987
void surfaceBlackbody(Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Index &atmosphere_dim, const Vector &f_grid, const Index &stokes_dim, const Vector &rtp_pos, const Vector &rtp_los, const Numeric &surface_skin_t, const Verbosity &verbosity)
WORKSPACE METHOD: surfaceBlackbody.
Definition: m_surface.cc:1937
void InterpSurfaceTypeMask(Index &surface_type, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lat_true, const Vector &lon_true, const Vector &rtp_pos, const GriddedField2 &surface_type_mask, const Verbosity &)
WORKSPACE METHOD: InterpSurfaceTypeMask.
Definition: m_surface.cc:3089
void SurfaceFastem(Matrix &surface_los, Tensor4 &surface_rmatrix, ArrayOfTensor4 &dsurface_rmatrix_dx, Matrix &surface_emission, ArrayOfMatrix &dsurface_emission_dx, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Tensor3 &surface_props_data, const ArrayOfString &surface_props_names, const ArrayOfString &dsurface_names, const Index &jacobian_do, const Vector &transmittance, const Index &fastem_version, const Verbosity &verbosity)
WORKSPACE METHOD: SurfaceFastem.
Definition: m_surface.cc:3514
void iySurfaceRtpropAgenda(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, Numeric &surface_skin_t, Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Tensor3 &iy_transmittance, const Index &iy_id, const Index &jacobian_do, const Index &suns_do, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Vector &rte_pos2, const String &iy_unit, const Agenda &iy_main_agenda, const Agenda &surface_rtprop_agenda, const Verbosity &)
WORKSPACE METHOD: iySurfaceRtpropAgenda.
Definition: m_surface.cc:1474
void surfaceFlatScalarReflectivity(Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Vector &f_grid, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &rtp_pos, const Vector &rtp_los, const Vector &specular_los, const Numeric &surface_skin_t, const Vector &surface_scalar_reflectivity, const Verbosity &)
WORKSPACE METHOD: surfaceFlatScalarReflectivity.
Definition: m_surface.cc:2566
void surface_reflectivityFromGriddedField6(Tensor3 &surface_reflectivity, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lat_true, const Vector &lon_true, const Vector &rtp_pos, const Vector &rtp_los, const GriddedField6 &r_field, const Verbosity &)
WORKSPACE METHOD: surface_reflectivityFromGriddedField6.
Definition: m_surface.cc:2830
void surfaceFlatReflectivity(Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Vector &f_grid, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &rtp_pos, const Vector &rtp_los, const Vector &specular_los, const Numeric &surface_skin_t, const Tensor3 &surface_reflectivity, const Verbosity &verbosity)
WORKSPACE METHOD: surfaceFlatReflectivity.
Definition: m_surface.cc:2402
void surfaceFlatRvRh(Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Vector &f_grid, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &rtp_pos, const Vector &rtp_los, const Vector &specular_los, const Numeric &surface_skin_t, const Matrix &surface_rv_rh, const Verbosity &)
WORKSPACE METHOD: surfaceFlatRvRh.
Definition: m_surface.cc:2486
void surface_rtpropFromTypesAverage(Workspace &ws, Vector &surface_type_mix, Numeric &surface_skin_t, Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Vector &f_grid, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Vector &lat_true, const Vector &lon_true, const Vector &rtp_pos, const Vector &rtp_los, const Vector &refellipsoid, const GriddedField2 &surface_type_mask, const ArrayOfAgenda &surface_rtprop_agenda_array, const Numeric &z_sensor, const Matrix &dlos, const Vector &dlos_weight_vector, const Verbosity &verbosity)
WORKSPACE METHOD: surface_rtpropFromTypesAverage.
Definition: m_surface.cc:3158
void iySurfaceLambertianDirect(Workspace &ws, Matrix &iy, const Vector &rtp_pos, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, 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 Matrix &z_surface, const Vector &surface_scalar_reflectivity, const Vector &refellipsoid, const Tensor4 &pnd_field, const ArrayOfTensor4 &dpnd_field_dx, const ArrayOfString &scat_species, const ArrayOfArrayOfSingleScatteringData &scat_data, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Index &suns_do, const Index &gas_scattering_do, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfSun &suns, const Numeric &rte_alonglos_v, const String &iy_unit, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &gas_scattering_agenda, const Agenda &ppath_step_agenda, const Verbosity &verbosity)
WORKSPACE METHOD: iySurfaceLambertianDirect.
Definition: m_surface.cc:1289
void transmittanceFromIy_aux(Vector &transmittance, const ArrayOfString &iy_aux_vars, const ArrayOfMatrix &iy_aux, const Verbosity &)
WORKSPACE METHOD: transmittanceFromIy_aux.
Definition: m_surface.cc:4063
void iySurfaceFlatReflectivity(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, ArrayOfTensor4 &dsurface_rmatrix_dx, ArrayOfMatrix &dsurface_emission_dx, const Tensor3 &iy_transmittance, const Index &iy_id, const Index &jacobian_do, const Index &suns_do, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, const Vector &f_grid, const Vector &lat_grid, const Vector &lon_grid, const Matrix &z_surface, const Vector &refellipsoid, const Vector &rtp_pos, const Vector &rtp_los, const Vector &rte_pos2, const String &iy_unit, const Tensor3 &surface_reflectivity, const Tensor3 &surface_props_data, const ArrayOfString &surface_props_names, const ArrayOfString &dsurface_names, const ArrayOfRetrievalQuantity &jacobian_quantities, const Agenda &iy_main_agenda, const Verbosity &verbosity)
WORKSPACE METHOD: iySurfaceFlatReflectivity.
Definition: m_surface.cc:353
void surfaceLambertianSimple(Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Vector &f_grid, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &rtp_pos, const Vector &rtp_los, const Vector &surface_normal, const Numeric &surface_skin_t, const Vector &surface_scalar_reflectivity, const Index &lambertian_nza, const Numeric &za_pos, const Verbosity &)
WORKSPACE METHOD: surfaceLambertianSimple.
Definition: m_surface.cc:2632
void iySurfaceFastem(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, const Tensor3 &iy_transmittance, const Index &iy_id, const Index &jacobian_do, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Vector &rte_pos2, const String &iy_unit, const Agenda &iy_main_agenda, const Numeric &surface_skin_t, const Numeric &salinity, const Numeric &wind_speed, const Numeric &wind_direction, const Index &fastem_version, const Verbosity &verbosity)
WORKSPACE METHOD: iySurfaceFastem.
Definition: m_surface.cc:234
void surface_scalar_reflectivityFromSurface_rmatrix(Vector &surface_scalar_reflectivity, const Tensor4 &surface_rmatrix, const Verbosity &)
WORKSPACE METHOD: surface_scalar_reflectivityFromSurface_rmatrix.
Definition: m_surface.cc:3071
void iySurfaceRtpropCalc(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, const Matrix &surface_los, const Tensor4 &surface_rmatrix, const Matrix &surface_emission, const ArrayOfString &dsurface_names, const ArrayOfTensor4 &dsurface_rmatrix_dx, const ArrayOfMatrix &dsurface_emission_dx, const Tensor3 &iy_transmittance, const Index &iy_id, const Index &jacobian_do, const Index &suns_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Vector &rte_pos2, const String &iy_unit, const Agenda &iy_main_agenda, const Verbosity &)
WORKSPACE METHOD: iySurfaceRtpropCalc.
Definition: m_surface.cc:1617
void SurfaceTessem(Matrix &surface_los, Tensor4 &surface_rmatrix, ArrayOfTensor4 &dsurface_rmatrix_dx, Matrix &surface_emission, ArrayOfMatrix &dsurface_emission_dx, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const TessemNN &net_h, const TessemNN &net_v, const Tensor3 &surface_props_data, const ArrayOfString &surface_props_names, const ArrayOfString &dsurface_names, const Index &jacobian_do, const Verbosity &verbosity)
WORKSPACE METHOD: SurfaceTessem.
Definition: m_surface.cc:3891
void SurfaceDummy(ArrayOfTensor4 &dsurface_rmatrix_dx, ArrayOfMatrix &dsurface_emission_dx, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &surface_props_data, const ArrayOfString &surface_props_names, const ArrayOfString &dsurface_names, const Index &jacobian_do, const Verbosity &)
WORKSPACE METHOD: SurfaceDummy.
Definition: m_surface.cc:3485
void surface_rtpropFromTypesManual(Workspace &ws, Vector &surface_type_mix, Numeric &surface_skin_t, Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const ArrayOfAgenda &surface_rtprop_agenda_array, const Index &surface_type, const Verbosity &)
WORKSPACE METHOD: surface_rtpropFromTypesManual.
Definition: m_surface.cc:3275
void iySurfaceInit(Matrix &iy, const Vector &f_grid, const Index &stokes_dim, const Verbosity &)
WORKSPACE METHOD: iySurfaceInit.
Definition: m_surface.cc:978
void iySurfaceFlatReflectivityDirect(Workspace &ws, Matrix &iy, const Vector &rtp_pos, const Vector &rtp_los, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, 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 Matrix &z_surface, const Tensor3 &surface_reflectivity, const Vector &refellipsoid, const Tensor4 &pnd_field, const ArrayOfTensor4 &dpnd_field_dx, const ArrayOfString &scat_species, const ArrayOfArrayOfSingleScatteringData &scat_data, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Index &ppath_inside_cloudbox_do, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Index &suns_do, const Index &gas_scattering_do, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfSun &suns, const Numeric &rte_alonglos_v, const String &iy_unit, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &gas_scattering_agenda, const Agenda &ppath_step_agenda, const Verbosity &verbosity)
WORKSPACE METHOD: iySurfaceFlatReflectivityDirect.
Definition: m_surface.cc:527
void InterpGriddedField2ToPosition(Numeric &outvalue, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lat_true, const Vector &lon_true, const Vector &rtp_pos, const GriddedField2 &gfield2, const Verbosity &)
WORKSPACE METHOD: InterpGriddedField2ToPosition.
Definition: m_surface.cc:128
void iySurfaceFlatRefractiveIndexDirect(Workspace &ws, Matrix &iy, const Vector &rtp_pos, const Vector &rtp_los, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, 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 Matrix &z_surface, const GriddedField3 &surface_complex_refr_index, const Vector &refellipsoid, const Tensor4 &pnd_field, const ArrayOfTensor4 &dpnd_field_dx, const ArrayOfString &scat_species, const ArrayOfArrayOfSingleScatteringData &scat_data, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Index &ppath_inside_cloudbox_do, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Index &suns_do, const Index &gas_scattering_do, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfSun &suns, const Numeric &rte_alonglos_v, const String &iy_unit, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &gas_scattering_agenda, const Agenda &ppath_step_agenda, const Verbosity &verbosity)
WORKSPACE METHOD: iySurfaceFlatRefractiveIndexDirect.
Definition: m_surface.cc:839
void FastemStandAlone(Matrix &emissivity, Matrix &reflectivity, const Vector &f_grid, const Numeric &surface_skin_t, const Numeric &za, const Numeric &salinity, const Numeric &wind_speed, const Numeric &rel_aa, const Vector &transmittance, const Index &fastem_version, const Verbosity &)
WORKSPACE METHOD: FastemStandAlone.
Definition: m_surface.cc:53
void SurfaceFlatScalarReflectivity(Matrix &surface_los, Tensor4 &surface_rmatrix, ArrayOfTensor4 &dsurface_rmatrix_dx, Matrix &surface_emission, ArrayOfMatrix &dsurface_emission_dx, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Vector &specular_los, const Tensor3 &surface_props_data, const ArrayOfString &surface_props_names, const ArrayOfString &dsurface_names, const Index &jacobian_do, const Vector &f_reflectivities, const Verbosity &verbosity)
WORKSPACE METHOD: SurfaceFlatScalarReflectivity.
Definition: m_surface.cc:3729
Numeric fac(const Index n)
fac
Definition: math_funcs.cc:46
Numeric sign(const Numeric &x)
sign
Definition: math_funcs.cc:445
void calculate_int_weights_arbitrary_grid(Vector &w, const Vector &x)
Calculates trapezoidal integration weights for arbitray grid.
Definition: math_funcs.cc:833
Declarations having to do with the four output streams.
#define CREATE_OUT3
Definition: messages.h:189
#define CREATE_OUT2
Definition: messages.h:188
my_basic_string< char > String
The String type for ARTS.
Definition: mystring.h:199
constexpr Numeric pi
The following mathematical constants are generated in python Decimal package by the code:
constexpr Numeric earth_radius
Global constant, the radius of the Earth [m].
constexpr auto deg2rad(auto x) noexcept
Converts degrees to radians.
Numeric planck(const Numeric &f, const Numeric &t)
planck
Numeric dplanck_dt(const Numeric &f, const Numeric &t)
dplanck_dt
void fresnel(Complex &Rv, Complex &Rh, const Complex &n1, const Complex &n2, const Numeric &theta)
fresnel
This file contains declerations of functions of physical character.
void zaaa2cart(Numeric &dx, Numeric &dy, Numeric &dz, const Numeric &za, const Numeric &aa)
Converts zenith and azimuth angles to a cartesian unit vector.
Definition: ppath.cc:296
void cart2zaaa(Numeric &za, Numeric &aa, const Numeric &dx, const Numeric &dy, const Numeric &dz)
Converts a cartesian directional vector to zenith and azimuth.
Definition: ppath.cc:283
void plevel_slope_3d(Numeric &c1, Numeric &c2, const Numeric &lat1, const Numeric &lat3, const Numeric &lon5, const Numeric &lon6, const Numeric &r15, const Numeric &r35, const Numeric &r36, const Numeric &r16, const Numeric &lat, const Numeric &lon, const Numeric &aa)
Definition: ppath.cc:1051
void resolve_lon(Numeric &lon, const Numeric &lon5, const Numeric &lon6)
Resolves which longitude angle that shall be used.
Definition: ppath.cc:484
Numeric plevel_angletilt(const Numeric &r, const Numeric &c1)
Calculates the angular tilt of the surface or a pressure level.
Definition: ppath.cc:600
void plevel_slope_2d(Numeric &c1, ConstVectorView lat_grid, ConstVectorView refellipsoid, ConstVectorView z_surf, const GridPos &gp, const Numeric &za)
Calculates the radial slope of the surface or a pressure level for 2D.
Definition: ppath.cc:563
Propagation path structure and functions.
void mueller_stokes2modif(Matrix &Cm, const Index &stokes_dim)
mueller_stokes2modif
Definition: rte.cc:1893
void pos2true_latlon(Numeric &lat, Numeric &lon, const Index &atmosphere_dim, const ConstVectorView &lat_grid, const ConstVectorView &lat_true, const ConstVectorView &lon_true, const ConstVectorView &pos)
Determines the true alt and lon for an "ARTS position".
Definition: rte.cc:1914
void iy_transmittance_mult(Tensor3 &iy_trans_total, const ConstTensor3View &iy_trans_old, const ConstTensor3View &iy_trans_new)
Multiplicates iy_transmittance with transmissions.
Definition: rte.cc:1755
void mueller_modif2stokes(Matrix &Cs, const Index &stokes_dim)
mueller_modif2stokes
Definition: rte.cc:1845
void mueller_rotation(Matrix &L, const Index &stokes_dim, const Numeric &rotangle)
mueller_rotation
Definition: rte.cc:1866
void mirror_los(Vector &los_mirrored, const ConstVectorView &los, const Index &atmosphere_dim)
Determines the backward direction for a given line-of-sight.
Definition: rte.cc:1792
Declaration of functions in rte.cc.
void complex_n_interp(MatrixView n_real, MatrixView n_imag, const GriddedField3 &complex_n, const String &varname, ConstVectorView f_grid, ConstVectorView t_grid)
General function for interpolating data of complex n type.
void interp_atmsurface_gp2itw(Matrix &itw, const Index &atmosphere_dim, const ArrayOfGridPos &gp_lat, const ArrayOfGridPos &gp_lon)
Converts atmospheric grid positions to weights for interpolation of a surface-type variable.
void rte_pos2gridpos(GridPos &gp_p, GridPos &gp_lat, GridPos &gp_lon, const Index &atmosphere_dim, ConstVectorView p_grid, ConstVectorView lat_grid, ConstVectorView lon_grid, ConstTensor3View z_field, ConstVectorView rte_pos)
Converts a geographical position (rte_pos) to grid positions for p, lat and lon.
void interp_atmsurface_by_gp(VectorView x, const Index &atmosphere_dim, ConstMatrixView x_surface, const ArrayOfGridPos &gp_lat, const ArrayOfGridPos &gp_lon)
Interpolates a surface-type variable given the grid positions.
Header file for special_interp.cc.
Structure to store a grid position.
Definition: interpolation.h:56
std::array< Numeric, 2 > fd
Definition: interpolation.h:58
Index idx
Definition: interpolation.h:57
The structure to describe a propagation path and releated quantities.
Definition: ppath_struct.h:17
void get_direct_radiation(Workspace &ws, ArrayOfMatrix &direct_radiation, ArrayOfArrayOfTensor3 &ddirect_radiation_dx, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_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 Index &gas_scattering_do, const Index &irradiance_flag, const ArrayOfPpath &sun_ppaths, const ArrayOfSun &suns, const ArrayOfIndex &suns_visible, const Vector &refellipsoid, const Tensor4 &pnd_field, const ArrayOfTensor4 &dpnd_field_dx, const ArrayOfString &scat_species, const ArrayOfArrayOfSingleScatteringData &scat_data, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &gas_scattering_agenda, const Numeric &rte_alonglos_v, const Verbosity &verbosity)
Calculates the transmitted sun radiation at the end position of the ppath.
Definition: sun.cc:195
void get_sun_ppaths(Workspace &ws, ArrayOfPpath &sun_ppaths, ArrayOfIndex &suns_visible, ArrayOfVector &sun_rte_los, const Vector &rte_pos, const ArrayOfSun &suns, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const Matrix &z_surface, const Vector &refellipsoid, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Agenda &ppath_step_agenda, const Verbosity &verbosity)
Calculates the ppath towards the suns from a given position and indicates if sun is visible or not.
Definition: sun.cc:371
void surface_get_incoming_direct(Workspace &ws, Matrix &iy_incoming, Index &stars_visible, Vector &specular_los, const Vector &rtp_pos, const Vector &rtp_los, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, 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 Matrix &z_surface, const Vector &refellipsoid, const Tensor4 &pnd_field, const ArrayOfTensor4 &dpnd_field_dx, const ArrayOfString &scat_species, const ArrayOfArrayOfSingleScatteringData &scat_data, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Index &ppath_inside_cloudbox_do, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Index &gas_scattering_do, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfSun &suns, const Numeric &rte_alonglos_v, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &gas_scattering_agenda, const Agenda &ppath_step_agenda, const Verbosity &verbosity)
Calculate the incoming direct radiation at the surface for a given line of sight.
Definition: surface.cc:213
void dsurface_check(const ArrayOfString &surface_props_names, const ArrayOfString &dsurface_names, const ArrayOfTensor4 dsurface_rmatrix_dx, const ArrayOfMatrix &dsurface_emission_dx)
Peforms basic checks of the dsurface variables.
Definition: surface.cc:185
void surface_calc(Matrix &iy, ConstTensor3View I, ConstMatrixView surface_los, ConstTensor4View surface_rmatrix, ConstMatrixView surface_emission)
Weights together downwelling radiation and surface emission.
Definition: surface.cc:46
void surface_props_check(const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &surface_props_data, const ArrayOfString &surface_props_names)
Peforms basic checks of surface_props_data and surface_props_names
Definition: surface.cc:122
Numeric calc_incang(ConstVectorView rte_los, ConstVectorView specular_los)
Calculates the incidence angle for a flat surface, based on rte_los and specular_los.
Definition: surface.cc:33
void surface_props_interp(Vector &v, const String &vname, const Index &atmosphere_dim, const ArrayOfGridPos &gp_lat, const ArrayOfGridPos &gp_lon, const Matrix &itw, const Tensor3 &surface_props_data, const ArrayOfString &surface_props_names)
Peforms an interpolation of surface_props_data
Definition: surface.cc:156
void surface_specular_R_and_b(MatrixView surface_rmatrix, VectorView surface_emission, const Complex &Rv, const Complex &Rh, const Numeric &f, const Index &stokes_dim, const Numeric &surface_skin_t)
Sets up the surface reflection matrix and emission vector for the case of specular reflection.
Definition: surface.cc:71
void tessem_prop_nn(VectorView ny, const TessemNN &net, ConstVectorView nx)
Definition: tessem.cc:68
This file contains functions that are adapted from TESSEM code which is used to calculate surface emi...
#define d
#define w
#define b