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