ARTS 2.5.0 (git: 9ee3ac6c)
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 "auto_md.h"
43#include "check_input.h"
44#include "matpack_complex.h"
45#include "fastem.h"
46#include "geodetic.h"
47#include "interpolation.h"
49#include "math_funcs.h"
50#include "messages.h"
51#include "physics_funcs.h"
52#include "ppath.h"
53#include "rte.h"
54#include "special_interp.h"
55#include "surface.h"
56#include "tessem.h"
57
59extern const Numeric DEG2RAD;
60
61/*===========================================================================
62 === The functions (in alphabetical order)
63 ===========================================================================*/
64
65/* Workspace method: Doxygen documentation will be auto-generated */
66void FastemStandAlone(Matrix& emissivity,
67 Matrix& reflectivity,
68 const Vector& f_grid,
69 const Numeric& surface_skin_t,
70 const Numeric& za,
71 const Numeric& salinity,
72 const Numeric& wind_speed,
73 const Numeric& rel_aa,
74 const Vector& transmittance,
75 const Index& fastem_version,
76 const Verbosity&) {
77 const Index nf = f_grid.nelem();
78
79 chk_if_in_range("zenith angle", za, 90, 180);
80 chk_if_in_range_exclude("surface skin temperature", surface_skin_t, 260, 373);
81 chk_if_in_range_exclude_high("salinity", salinity, 0, 1);
82 chk_if_in_range_exclude_high("wind speed", wind_speed, 0, 100);
83 chk_if_in_range("azimuth angle", rel_aa, -180, 180);
84 chk_vector_length("transmittance", "f_grid", transmittance, f_grid);
85 if (fastem_version < 3 || fastem_version > 6)
86 throw std::runtime_error(
87 "Invalid fastem version: 3 <= fastem_version <= 6");
88
89 emissivity.resize(nf, 4);
90 reflectivity.resize(nf, 4);
91
92 const Numeric t = max(surface_skin_t, Numeric(270));
93
94 for (Index i = 0; i < nf; i++) {
95 if (f_grid[i] > 250e9)
96 throw std::runtime_error("Only frequency <= 250 GHz are allowed");
97 chk_if_in_range("transmittance", transmittance[i], 0, 1);
98
99 Vector e, r;
100
101 fastem(e,
102 r,
103 f_grid[i],
104 za,
105 t,
106 salinity,
107 wind_speed,
108 transmittance[i],
109 rel_aa,
110 fastem_version);
111
112 emissivity(i, joker) = e;
113 reflectivity(i, joker) = r;
114 }
115
116 // FASTEM does not work close to the horizon (at least v6). Make sure values
117 // are inside [0,1]. Then seems best to make sure that e+r=1.
118 for (Index i = 0; i < nf; i++) {
119 for (Index s = 0; s < 2; s++) {
120 if (emissivity(i, s) > 1) {
121 emissivity(i, s) = 1;
122 reflectivity(i, s) = 0;
123 }
124 if (emissivity(i, s) < 0) {
125 emissivity(i, s) = 0;
126 reflectivity(i, s) = 1;
127 }
128 if (reflectivity(i, s) > 1) {
129 emissivity(i, s) = 0;
130 reflectivity(i, s) = 1;
131 }
132 if (reflectivity(i, s) < 0) {
133 emissivity(i, s) = 1;
134 reflectivity(i, s) = 0;
135 }
136 }
137 }
138}
139
140/* Workspace method: Doxygen documentation will be auto-generated */
142 const Index& atmosphere_dim,
143 const Vector& lat_grid,
144 const Vector& lat_true,
145 const Vector& lon_true,
146 const Vector& rtp_pos,
147 const GriddedField2& gfield2,
148 const Verbosity&) {
149 // Set expected order of grids
150 Index gfield_latID = 0;
151 Index gfield_lonID = 1;
152
153 // Basic checks and sizes
154 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
155 chk_latlon_true(atmosphere_dim, lat_grid, lat_true, lon_true);
156 chk_rte_pos(atmosphere_dim, rtp_pos);
157 gfield2.checksize_strict();
158 //
159 chk_griddedfield_gridname(gfield2, gfield_latID, "Latitude");
160 chk_griddedfield_gridname(gfield2, gfield_lonID, "Longitude");
161 //
162 const Index nlat = gfield2.data.nrows();
163 const Index nlon = gfield2.data.ncols();
164 //
165 if (nlat < 2 || nlon < 2) {
166 ostringstream os;
167 os << "The data in *gfield2* must span a geographical region. That is,\n"
168 << "the latitude and longitude grids must have a length >= 2.";
169 throw runtime_error(os.str());
170 }
171
172 const Vector& GFlat = gfield2.get_numeric_grid(gfield_latID);
173 const Vector& GFlon = gfield2.get_numeric_grid(gfield_lonID);
174
175 // Determine true geographical position
176 Vector lat(1), lon(1);
178 lat[0], lon[0], atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
179
180 // Ensure correct coverage of lon grid
181 Vector lon_shifted;
182 lon_shiftgrid(lon_shifted, GFlon, lon[0]);
183
184 // Check if lat/lon we need are actually covered
185 chk_if_in_range("rtp_pos.lat", lat[0], GFlat[0], GFlat[nlat - 1]);
186 chk_if_in_range("rtp_pos.lon", lon[0], lon_shifted[0], lon_shifted[nlon - 1]);
187
188 // Interpolate in lat and lon
189 //
190 GridPos gp_lat, gp_lon;
191 gridpos(gp_lat, GFlat, lat[0]);
192 gridpos(gp_lon, lon_shifted, lon[0]);
193 Vector itw(4);
194 interpweights(itw, gp_lat, gp_lon);
195 outvalue = interp(itw, gfield2.data, gp_lat, gp_lon);
196}
197
198
199/* Workspace method: Doxygen documentation will be auto-generated */
201 const Index& atmosphere_dim,
202 const Vector& lat_grid,
203 const Vector& lon_grid,
204 const Vector& rtp_pos,
205 const Matrix& z_surface,
206 const Matrix& field,
207 const Verbosity& verbosity) {
208 // Input checks (dummy p_grid)
209 chk_atm_grids(atmosphere_dim, Vector(2, 2, -1), lat_grid, lon_grid);
211 "input argument *field*", field, atmosphere_dim, lat_grid, lon_grid);
212 chk_rte_pos(atmosphere_dim, rtp_pos);
213 //
214 const Numeric zmax = max(z_surface);
215 const Numeric zmin = min(z_surface);
216 const Numeric dzok = 1;
217 if (rtp_pos[0] < zmin - dzok || rtp_pos[0] > zmax + dzok) {
218 ostringstream os;
219 os << "The given position does not match *z_surface*.\nThe altitude in "
220 << "*rtp_pos* is " << rtp_pos[0] / 1e3 << " km.\n"
221 << "The altitude range covered by *z_surface* is [" << zmin / 1e3 << ","
222 << zmax / 1e3 << "] km.\n"
223 << "One possible mistake is to mix up *rtp_pos* and *rte_los*.";
224 throw runtime_error(os.str());
225 }
226
227 if (atmosphere_dim == 1) {
228 outvalue = field(0, 0);
229 } else {
230 chk_interpolation_grids("Latitude interpolation", lat_grid, rtp_pos[1]);
231 GridPos gp_lat, gp_lon;
232 gridpos(gp_lat, lat_grid, rtp_pos[1]);
233 if (atmosphere_dim == 3) {
234 chk_interpolation_grids("Longitude interpolation", lon_grid, rtp_pos[2]);
235 gridpos(gp_lon, lon_grid, rtp_pos[2]);
236 }
237 //
238 outvalue = interp_atmsurface_by_gp(atmosphere_dim, field, gp_lat, gp_lon);
239 }
240
241 // Interpolate
243 out3 << " Result = " << outvalue << "\n";
244}
245
246/* Workspace method: Doxygen documentation will be auto-generated */
248 Matrix& iy,
249 ArrayOfTensor3& diy_dx,
250 const String& iy_unit,
251 const Tensor3& iy_transmittance,
252 const Index& iy_id,
253 const Index& cloudbox_on,
254 const Index& jacobian_do,
255 const Vector& f_grid,
256 const Agenda& iy_main_agenda,
257 const Vector& rtp_pos,
258 const Vector& rtp_los,
259 const Vector& rte_pos2,
260 const ArrayOfAgenda& iy_surface_agenda_array,
261 const ArrayOfIndex& surface_types,
262 const Vector& surface_types_aux,
263 const Vector& surface_types_weights,
264 const Verbosity&) {
265
266 const Index ntypes = surface_types.nelem();
267 if (ntypes < 1)
268 throw runtime_error("*surface_types* is empty!");
269
270 // Loop surface types and sum up
271 //
272 Numeric wtot = 0;
273 //
274 for (Index t=0; t<ntypes; t++ ) {
275
276 if (surface_types[t] < 0)
277 throw runtime_error(
278 "No element in *surface_types* is allowed to be negative.");
279 if (surface_types[t] >= iy_surface_agenda_array.nelem()) {
280 ostringstream os;
281 os << "*iy_surface_agenda_array* has only "
282 << iy_surface_agenda_array.nelem()
283 << " elements,\n while you have selected element " << surface_types[t];
284 throw runtime_error(os.str());
285 }
286
287 Matrix iy1;
288 ArrayOfTensor3 diy_dx1;
289
291 iy1,
292 diy_dx1,
293 surface_types[t],
294 iy_unit,
295 iy_transmittance,
296 iy_id,
297 cloudbox_on,
298 jacobian_do,
299 iy_main_agenda,
300 f_grid,
301 rtp_pos,
302 rtp_los,
303 rte_pos2,
304 surface_types_aux[t],
305 iy_surface_agenda_array);
306
307 iy1 *= surface_types_weights[t];
308 for (Index i=0; diy_dx1.nelem(); i++ )
309 diy_dx1[i] *= surface_types_weights[t];
310 wtot += surface_types_weights[t];;
311
312 if (t==0) {
313 iy = iy1;
314 diy_dx = diy_dx1;
315 } else {
316 iy += iy1;
317 for (Index i=0; diy_dx1.nelem(); i++ )
318 diy_dx[i] += diy_dx1[i];
319 }
320 }
321 if (abs(wtot-1)>1e-4)
322 throw runtime_error("Sum of *surface_types_weights* deviates from 1.");
323}
324
325/* Workspace method: Doxygen documentation will be auto-generated */
327 Matrix& iy,
328 ArrayOfTensor3& diy_dx,
329 const Tensor3& iy_transmittance,
330 const Index& iy_id,
331 const Index& jacobian_do,
332 const Index& atmosphere_dim,
333 const EnergyLevelMap& nlte_field,
334 const Index& cloudbox_on,
335 const Index& stokes_dim,
336 const Vector& f_grid,
337 const Vector& rtp_pos,
338 const Vector& rtp_los,
339 const Vector& rte_pos2,
340 const String& iy_unit,
341 const Agenda& iy_main_agenda,
342 const Numeric& surface_skin_t,
343 const Numeric& salinity,
344 const Numeric& wind_speed,
345 const Numeric& wind_direction,
346 const Index& fastem_version,
347 const Verbosity& verbosity) {
348 // Most obvious input checks are performed in specular_losCalc and surfaceFastem
349
350 // Obtain radiance and transmission for specular direction
351
352 // Determine specular direction
353 Vector specular_los, surface_normal;
354 specular_losCalcNoTopography(specular_los,
355 surface_normal,
356 rtp_pos,
357 rtp_los,
358 atmosphere_dim,
359 verbosity);
360
361 // Use iy_aux to get optical depth for downwelling radiation.
362 ArrayOfString iy_aux_vars(1);
363 iy_aux_vars[0] = "Optical depth";
364
365 // Calculate iy for downwelling radiation
366 // Note that iy_transmittance used here lacks surface R. Fixed below.
367 //
368 const Index nf = f_grid.nelem();
369 Vector transmittance(nf);
370 ArrayOfMatrix iy_aux;
371 Ppath ppath;
372 //
374 iy,
375 iy_aux,
376 ppath,
377 diy_dx,
378 0,
379 iy_transmittance,
380 iy_aux_vars,
381 iy_id,
382 iy_unit,
383 cloudbox_on,
384 jacobian_do,
385 f_grid,
386 nlte_field,
387 rtp_pos,
388 specular_los,
389 rte_pos2,
390 iy_main_agenda);
391
392 // Convert tau to transmissions
393 for (Index i = 0; i < nf; i++) {
394 transmittance[i] = exp(-iy_aux[0](i, 0));
395 }
396
397 // Call Fastem by surface_RTprop version
398 //
399 Matrix surface_los;
400 Tensor4 surface_rmatrix;
401 Matrix surface_emission;
402 //
403 surfaceFastem(surface_los,
404 surface_rmatrix,
405 surface_emission,
406 atmosphere_dim,
407 stokes_dim,
408 f_grid,
409 rtp_pos,
410 rtp_los,
411 surface_skin_t,
412 salinity,
413 wind_speed,
414 wind_direction,
415 transmittance,
416 fastem_version,
417 verbosity);
418
419 // Add up
420 //
421 Tensor3 I(1, nf, stokes_dim);
422 I(0, joker, joker) = iy;
423 Matrix sensor_los_dummy(1, 1, 0);
424 //
425 surface_calc(iy, I, sensor_los_dummy, surface_rmatrix, surface_emission);
426
427 // Adjust diy_dx, if necessary.
428 // For vector cases this is a slight approximation, as the order of the
429 // different transmission and reflectivities matters.
430 if (iy_transmittance.npages()) {
431 for (Index q = 0; q < diy_dx.nelem(); q++) {
432 for (Index p = 0; p < diy_dx[q].npages(); p++) {
433 for (Index i = 0; i < nf; i++) {
434 Vector x = diy_dx[q](p, i, joker);
435 mult(diy_dx[q](p, i, joker), surface_rmatrix(0, i, joker, joker), x);
436 }
437 }
438 }
439 }
440}
441
442/* Workspace method: Doxygen documentation will be auto-generated */
444 Matrix& iy,
445 ArrayOfTensor3& diy_dx,
446 Numeric& surface_skin_t,
447 Matrix& surface_los,
448 Tensor4& surface_rmatrix,
449 Matrix& surface_emission,
450 const Tensor3& iy_transmittance,
451 const Index& iy_id,
452 const Index& jacobian_do,
453 const Index& atmosphere_dim,
454 const EnergyLevelMap& nlte_field,
455 const Index& cloudbox_on,
456 const Index& stokes_dim,
457 const Vector& f_grid,
458 const Vector& rtp_pos,
459 const Vector& rtp_los,
460 const Vector& rte_pos2,
461 const String& iy_unit,
462 const Agenda& iy_main_agenda,
463 const Agenda& surface_rtprop_agenda,
464 const Verbosity&) {
465 // Input checks
466 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
467 chk_rte_pos(atmosphere_dim, rtp_pos);
468 chk_rte_los(atmosphere_dim, rtp_los);
469
470 // Call *surface_rtprop_agenda*
472 surface_skin_t,
473 surface_emission,
474 surface_los,
475 surface_rmatrix,
476 f_grid,
477 rtp_pos,
478 rtp_los,
479 surface_rtprop_agenda);
480
481 // Check output of *surface_rtprop_agenda*
482 const Index nlos = surface_los.nrows();
483 const Index nf = f_grid.nelem();
484 //
485 if (nlos) // if 0, blackbody ground and not all checks are needed
486 {
487 if (surface_los.ncols() != rtp_los.nelem())
488 throw runtime_error("Number of columns in *surface_los* is not correct.");
489 if (nlos != surface_rmatrix.nbooks())
490 throw runtime_error(
491 "Mismatch in size of *surface_los* and *surface_rmatrix*.");
492 if (surface_rmatrix.npages() != nf)
493 throw runtime_error(
494 "Mismatch in size of *surface_rmatrix* and *f_grid*.");
495 if (surface_rmatrix.nrows() != stokes_dim ||
496 surface_rmatrix.ncols() != stokes_dim)
497 throw runtime_error(
498 "Mismatch between size of *surface_rmatrix* and *stokes_dim*.");
499 }
500 if (surface_emission.ncols() != stokes_dim)
501 throw runtime_error(
502 "Mismatch between size of *surface_emission* and *stokes_dim*.");
503 if (surface_emission.nrows() != nf)
504 throw runtime_error("Mismatch in size of *surface_emission* and f_grid*.");
505
506 // Variable to hold down-welling radiation
507 Tensor3 I(nlos, nf, stokes_dim);
508
509 // Loop *surface_los*-es. If no such LOS, we are ready.
510 if (nlos > 0) {
511 for (Index ilos = 0; ilos < nlos; ilos++) {
512 Vector los = surface_los(ilos, joker);
513
514 // Include surface reflection matrix in *iy_transmittance*
515 // If iy_transmittance is empty, this is interpreted as the
516 // variable is not needed.
517 //
518 Tensor3 iy_trans_new;
519 //
520 if (iy_transmittance.npages()) {
521 iy_transmittance_mult(iy_trans_new,
522 iy_transmittance,
523 surface_rmatrix(ilos, joker, joker, joker));
524 }
525
526 // Calculate downwelling radiation for LOS ilos
527 //
528 {
529 ArrayOfMatrix iy_aux;
530 Ppath ppath;
531 Index iy_id_new = iy_id + ilos + 1;
533 iy,
534 iy_aux,
535 ppath,
536 diy_dx,
537 0,
538 iy_trans_new,
539 ArrayOfString(0),
540 iy_id_new,
541 iy_unit,
542 cloudbox_on,
543 jacobian_do,
544 f_grid,
545 nlte_field,
546 rtp_pos,
547 los,
548 rte_pos2,
549 iy_main_agenda);
550 }
551
552 if (iy.ncols() != stokes_dim || iy.nrows() != nf) {
553 ostringstream os;
554 os << "The size of *iy* returned from *" << iy_main_agenda.name()
555 << "* is\n"
556 << "not correct:\n"
557 << " expected size = [" << nf << "," << stokes_dim << "]\n"
558 << " size of iy = [" << iy.nrows() << "," << iy.ncols() << "]\n";
559 throw runtime_error(os.str());
560 }
561
562 I(ilos, joker, joker) = iy;
563 }
564 }
565
566 // Add up
567 surface_calc(iy, I, surface_los, surface_rmatrix, surface_emission);
568}
569
570/* Workspace method: Doxygen documentation will be auto-generated */
572 Matrix& iy,
573 ArrayOfTensor3& diy_dx,
574 const Matrix& surface_los,
575 const Tensor4& surface_rmatrix,
576 const Matrix& surface_emission,
577 const ArrayOfString& dsurface_names,
578 const ArrayOfTensor4& dsurface_rmatrix_dx,
579 const ArrayOfMatrix& dsurface_emission_dx,
580 const Tensor3& iy_transmittance,
581 const Index& iy_id,
582 const Index& jacobian_do,
583 const ArrayOfRetrievalQuantity& jacobian_quantities,
584 const Index& atmosphere_dim,
585 const EnergyLevelMap& nlte_field,
586 const Index& cloudbox_on,
587 const Index& stokes_dim,
588 const Vector& f_grid,
589 const Vector& rtp_pos,
590 const Vector& rtp_los,
591 const Vector& rte_pos2,
592 const String& iy_unit,
593 const Agenda& iy_main_agenda,
594 const Verbosity&) {
595 // Input checks
596 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
597 chk_rte_pos(atmosphere_dim, rtp_pos);
598 chk_rte_los(atmosphere_dim, rtp_los);
599
600 // Check provided surface rtprop variables
601 const Index nlos = surface_los.nrows();
602 const Index nf = f_grid.nelem();
603 //
604 if (nlos) // if 0, blackbody ground and not all checks are needed
605 {
606 if (surface_los.ncols() != rtp_los.nelem())
607 throw runtime_error("Number of columns in *surface_los* is not correct.");
608 if (nlos != surface_rmatrix.nbooks())
609 throw runtime_error(
610 "Mismatch in size of *surface_los* and *surface_rmatrix*.");
611 if (surface_rmatrix.npages() != nf)
612 throw runtime_error(
613 "Mismatch in size of *surface_rmatrix* and *f_grid*.");
614 if (surface_rmatrix.nrows() != stokes_dim ||
615 surface_rmatrix.ncols() != stokes_dim)
616 throw runtime_error(
617 "Mismatch between size of *surface_rmatrix* and *stokes_dim*.");
618 }
619 if (surface_emission.ncols() != stokes_dim)
620 throw runtime_error(
621 "Mismatch between size of *surface_emission* and *stokes_dim*.");
622 if (surface_emission.nrows() != nf)
623 throw runtime_error("Mismatch in size of *surface_emission* and f_grid*.");
624
625 // Variable to hold down-welling radiation
626 Tensor3 I(nlos, nf, stokes_dim);
627
628 // Loop *surface_los*-es.
629 if (nlos > 0) {
630 for (Index ilos = 0; ilos < nlos; ilos++) {
631 Vector los = surface_los(ilos, joker);
632
633 // Include surface reflection matrix in *iy_transmittance*
634 // If iy_transmittance is empty, this is interpreted as the
635 // variable is not needed.
636 //
637 Tensor3 iy_trans_new;
638 //
639 if (iy_transmittance.npages()) {
640 iy_transmittance_mult(iy_trans_new,
641 iy_transmittance,
642 surface_rmatrix(ilos, joker, joker, joker));
643 }
644
645 // Calculate downwelling radiation for LOS ilos
646 //
647 {
648 ArrayOfMatrix iy_aux;
649 Ppath ppath;
651 iy,
652 iy_aux,
653 ppath,
654 diy_dx,
655 0,
656 iy_trans_new,
657 ArrayOfString(0),
658 iy_id,
659 iy_unit,
660 cloudbox_on,
661 jacobian_do,
662 f_grid,
663 nlte_field,
664 rtp_pos,
665 los,
666 rte_pos2,
667 iy_main_agenda);
668 }
669
670 if (iy.ncols() != stokes_dim || iy.nrows() != nf) {
671 ostringstream os;
672 os << "The size of *iy* returned from *" << iy_main_agenda.name()
673 << "* is\n"
674 << "not correct:\n"
675 << " expected size = [" << nf << "," << stokes_dim << "]\n"
676 << " size of iy = [" << iy.nrows() << "," << iy.ncols() << "]\n";
677 throw runtime_error(os.str());
678 }
679
680 I(ilos, joker, joker) = iy;
681 }
682 }
683
684 // Add up
685 surface_calc(iy, I, surface_los, surface_rmatrix, surface_emission);
686
687 // Surface Jacobians
688 if (jacobian_do && dsurface_names.nelem()) {
689 // Loop dsurface_names
690 for (Index i = 0; i < dsurface_names.nelem(); i++) {
691 // Error if derivatives not calculated
692 // Or should we accept this?
693 if (dsurface_emission_dx[i].empty() || dsurface_rmatrix_dx[i].empty()) {
694 ostringstream os;
695 os << "The derivatives for surface quantity: " << dsurface_names[i]
696 << "\nwere not calculated by *iy_surface_agenda*.\n"
697 << "That is, *dsurface_emission_dx* and/or *dsurface_rmatrix_dx*\n"
698 << "are empty.";
699 throw runtime_error(os.str());
700 } else {
701 // Find index among jacobian quantities
702 Index ihit = -1;
703 for (Index j = 0; j < jacobian_quantities.nelem() && ihit < 0; j++) {
704 if (dsurface_names[i] == jacobian_quantities[j].Subtag()) {
705 ihit = j;
706 }
707 }
708 ARTS_ASSERT(ihit >= 0);
709 // Derivative, as observed at the surface
710 Matrix diy_dpos0, diy_dpos;
711 surface_calc(diy_dpos0,
712 I,
713 surface_los,
714 dsurface_rmatrix_dx[i],
715 dsurface_emission_dx[i]);
716 // Weight with transmission to sensor
717 iy_transmittance_mult(diy_dpos, iy_transmittance, diy_dpos0);
718 // Put into diy_dx
719 diy_from_pos_to_rgrids(diy_dx[ihit],
720 jacobian_quantities[ihit],
721 diy_dpos,
722 atmosphere_dim,
723 rtp_pos);
724 }
725 }
726 }
727}
728
729/* Workspace method: Doxygen documentation will be auto-generated */
731 Vector& surface_normal,
732 const Vector& rtp_pos,
733 const Vector& rtp_los,
734 const Index& atmosphere_dim,
735 const Verbosity&) {
736 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
737 chk_rte_pos(atmosphere_dim, rtp_pos);
738 chk_rte_los(atmosphere_dim, rtp_los);
739
740 surface_normal.resize(max(Index(1), atmosphere_dim - 1));
741 specular_los.resize(max(Index(1), atmosphere_dim - 1));
742
743 if (atmosphere_dim == 1) {
744 surface_normal[0] = 0;
745 if (rtp_los[0] < 90) {
746 throw runtime_error(
747 "Invalid zenith angle. The zenith angle corresponds "
748 "to observe the surface from below.");
749 }
750 specular_los[0] = 180 - rtp_los[0];
751 }
752
753 else if (atmosphere_dim == 2) {
754 specular_los[0] = sign(rtp_los[0]) * 180 - rtp_los[0];
755 surface_normal[0] = 0;
756 }
757
758 else if (atmosphere_dim == 3) {
759 specular_los[0] = 180 - rtp_los[0];
760 specular_los[1] = rtp_los[1];
761 surface_normal[0] = 0;
762 surface_normal[1] = 0;
763 }
764}
765
766/* Workspace method: Doxygen documentation will be auto-generated */
767void specular_losCalc(Vector& specular_los,
768 Vector& surface_normal,
769 const Vector& rtp_pos,
770 const Vector& rtp_los,
771 const Index& atmosphere_dim,
772 const Vector& lat_grid,
773 const Vector& lon_grid,
774 const Vector& refellipsoid,
775 const Matrix& z_surface,
776 const Index& ignore_surface_slope,
777 const Verbosity& verbosity) {
778 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
779 chk_rte_pos(atmosphere_dim, rtp_pos);
780 chk_rte_los(atmosphere_dim, rtp_los);
781 chk_if_in_range("ignore_surface_slope", ignore_surface_slope, 0, 1);
782
783 // Use special function if there is no slope, or it is ignored
784 if (atmosphere_dim == 1 || ignore_surface_slope) {
785 specular_losCalcNoTopography(specular_los,
786 surface_normal,
787 rtp_pos,
788 rtp_los,
789 atmosphere_dim,
790 verbosity);
791 return;
792 }
793
794 surface_normal.resize(max(Index(1), atmosphere_dim - 1));
795 specular_los.resize(max(Index(1), atmosphere_dim - 1));
796
797 if (atmosphere_dim == 2) {
798 chk_interpolation_grids("Latitude interpolation", lat_grid, rtp_pos[1]);
799 GridPos gp_lat;
800 gridpos(gp_lat, lat_grid, rtp_pos[1]);
801 Numeric c1; // Radial slope of the surface
803 c1, lat_grid, refellipsoid, z_surface(joker, 0), gp_lat, rtp_los[0]);
804 Vector itw(2);
805 interpweights(itw, gp_lat);
806 const Numeric rv_surface = refell2d(refellipsoid, lat_grid, gp_lat) +
807 interp(itw, z_surface(joker, 0), gp_lat);
808 surface_normal[0] = -plevel_angletilt(rv_surface, c1);
809 if (abs(rtp_los[0] - surface_normal[0]) < 90) {
810 throw runtime_error(
811 "Invalid zenith angle. The zenith angle corresponds "
812 "to observe the surface from below.");
813 }
814 specular_los[0] =
815 sign(rtp_los[0]) * 180 - rtp_los[0] + 2 * surface_normal[0];
816 }
817
818 else {
819 // Calculate surface normal in South-North direction
820 chk_interpolation_grids("Latitude interpolation", lat_grid, rtp_pos[1]);
821 chk_interpolation_grids("Longitude interpolation", lon_grid, rtp_pos[2]);
822 GridPos gp_lat, gp_lon;
823 gridpos(gp_lat, lat_grid, rtp_pos[1]);
824 gridpos(gp_lon, lon_grid, rtp_pos[2]);
825 Numeric c1, c2;
827 c1, c2, lat_grid, lon_grid, refellipsoid, z_surface, gp_lat, gp_lon, 0);
828 Vector itw(4);
829 interpweights(itw, gp_lat, gp_lon);
830 const Numeric rv_surface = refell2d(refellipsoid, lat_grid, gp_lat) +
831 interp(itw, z_surface, gp_lat, gp_lon);
832 const Numeric zaSN = 90 - plevel_angletilt(rv_surface, c1);
833 // The same for East-West
835 c2,
836 lat_grid,
837 lon_grid,
838 refellipsoid,
839 z_surface,
840 gp_lat,
841 gp_lon,
842 90);
843 const Numeric zaEW = 90 - plevel_angletilt(rv_surface, c1);
844 // Convert to Cartesian, and determine normal by cross-product
845 Vector tangentSN(3), tangentEW(3), normal(3);
846 zaaa2cart(tangentSN[0], tangentSN[1], tangentSN[2], zaSN, 0);
847 zaaa2cart(tangentEW[0], tangentEW[1], tangentEW[2], zaEW, 90);
848 cross3(normal, tangentSN, tangentEW);
849 // Convert rtp_los to cartesian and flip direction
850 Vector di(3);
851 zaaa2cart(di[0], di[1], di[2], rtp_los[0], rtp_los[1]);
852 di *= -1;
853 // Set LOS normal vector
854 cart2zaaa(
855 surface_normal[0], surface_normal[1], normal[0], normal[1], normal[2]);
856 if (abs(rtp_los[0] - surface_normal[0]) < 90) {
857 throw runtime_error(
858 "Invalid zenith angle. The zenith angle corresponds "
859 "to observe the surface from below.");
860 }
861 // Specular direction is 2(dn*di)dn-di, where dn is the normal vector
862 Vector speccart(3);
863 const Numeric fac = 2 * (normal * di);
864 for (Index i = 0; i < 3; i++) {
865 speccart[i] = fac * normal[i] - di[i];
866 }
867 cart2zaaa(specular_los[0],
868 specular_los[1],
869 speccart[0],
870 speccart[1],
871 speccart[2]);
872 }
873}
874
875/* Workspace method: Doxygen documentation will be auto-generated */
876void surfaceBlackbody(Matrix& surface_los,
877 Tensor4& surface_rmatrix,
878 Matrix& surface_emission,
879 const Index& atmosphere_dim,
880 const Vector& f_grid,
881 const Index& stokes_dim,
882 const Vector& rtp_pos,
883 const Vector& rtp_los,
884 const Numeric& surface_skin_t,
885 const Verbosity& verbosity) {
886 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
887 chk_rte_pos(atmosphere_dim, rtp_pos);
888 chk_rte_los(atmosphere_dim, rtp_los);
889 chk_not_negative("surface_skin_t", surface_skin_t);
890
892 out2 << " Sets variables to model a blackbody surface with a temperature "
893 << " of " << surface_skin_t << " K.\n";
894
895 surface_los.resize(0, 0);
896 surface_rmatrix.resize(0, 0, 0, 0);
897
898 const Index nf = f_grid.nelem();
899
900 Vector b(nf);
901 planck(b, f_grid, surface_skin_t);
902
903 surface_emission.resize(nf, stokes_dim);
904 surface_emission = 0.0;
905
906 for (Index iv = 0; iv < nf; iv++) {
907 surface_emission(iv, 0) = b[iv];
908 for (Index is = 1; is < stokes_dim; is++) {
909 surface_emission(iv, is) = 0;
910 }
911 }
912}
913
914/* Workspace method: Doxygen documentation will be auto-generated */
915void surfaceFastem(Matrix& surface_los,
916 Tensor4& surface_rmatrix,
917 Matrix& surface_emission,
918 const Index& atmosphere_dim,
919 const Index& stokes_dim,
920 const Vector& f_grid,
921 const Vector& rtp_pos,
922 const Vector& rtp_los,
923 const Numeric& surface_skin_t,
924 const Numeric& salinity,
925 const Numeric& wind_speed,
926 const Numeric& wind_direction,
927 const Vector& transmittance,
928 const Index& fastem_version,
929 const Verbosity& verbosity) {
930 // Input checks
931 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
932 chk_rte_pos(atmosphere_dim, rtp_pos);
933 chk_rte_los(atmosphere_dim, rtp_los);
934 chk_if_in_range("wind_direction", wind_direction, -180, 180);
935
936 const Index nf = f_grid.nelem();
937
938 // Determine specular direction
939 Vector specular_los, surface_normal;
940 specular_losCalcNoTopography(specular_los,
941 surface_normal,
942 rtp_pos,
943 rtp_los,
944 atmosphere_dim,
945 verbosity);
946
947 // Determine relative azimuth
948 //
949 // According to email from James Hocking, UkMet:
950 // All azimuth angles are defined as being measured clockwise from north
951 // (i.e. if the projection onto the Earth's surface of the view path lies due
952 // north the azimuth angle is 0 and if it lies due east the azimuth angle is
953 // 90 degrees). The relative azimuth is the wind direction (azimuth) angle
954 // minus the satellite azimuth angle.
955 Numeric rel_azimuth = wind_direction; // Always true for 1D
956 if (atmosphere_dim == 2 && rtp_los[0] < 0) {
957 rel_azimuth -= 180;
958 resolve_lon(rel_azimuth, -180, 180);
959 } else if (atmosphere_dim == 3) {
960 rel_azimuth -= rtp_los[1];
961 resolve_lon(rel_azimuth, -180, 180);
962 }
963
964 // Call FASTEM
965 Matrix emissivity, reflectivity;
966 FastemStandAlone(emissivity,
967 reflectivity,
968 f_grid,
969 surface_skin_t,
970 abs(rtp_los[0]),
971 salinity,
972 wind_speed,
973 rel_azimuth,
974 transmittance,
975 fastem_version,
976 verbosity);
977
978 // Set surface_los
979 surface_los.resize(1, specular_los.nelem());
980 surface_los(0, joker) = specular_los;
981
982 // Surface emission
983 //
984 Vector b(nf);
985 planck(b, f_grid, surface_skin_t);
986 //
987 surface_emission.resize(nf, stokes_dim);
988 for (Index i = 0; i < nf; i++) {
989 // I
990 surface_emission(i, 0) = b[i] * 0.5 * (emissivity(i, 0) + emissivity(i, 1));
991 // Q
992 if (stokes_dim >= 2) {
993 surface_emission(i, 1) =
994 b[i] * 0.5 * (emissivity(i, 0) - emissivity(i, 1));
995 }
996 // U and V
997 for (Index j = 2; j < stokes_dim; j++) {
998 surface_emission(i, j) = b[i] * emissivity(i, j);
999 }
1000 }
1001
1002 // Surface reflectivity matrix
1003 //
1004 surface_rmatrix.resize(1, nf, stokes_dim, stokes_dim);
1005 surface_rmatrix = 0.0;
1006 for (Index iv = 0; iv < nf; iv++) {
1007 surface_rmatrix(0, iv, 0, 0) =
1008 0.5 * (reflectivity(iv, 0) + reflectivity(iv, 1));
1009 if (stokes_dim >= 2) {
1010 surface_rmatrix(0, iv, 0, 1) =
1011 0.5 * (reflectivity(iv, 0) - reflectivity(iv, 1));
1012 ;
1013 surface_rmatrix(0, iv, 1, 0) = surface_rmatrix(0, iv, 0, 1);
1014 surface_rmatrix(0, iv, 1, 1) = surface_rmatrix(0, iv, 0, 0);
1015
1016 for (Index i = 2; i < stokes_dim; i++) {
1017 surface_rmatrix(0, iv, i, i) = surface_rmatrix(0, iv, 0, 0);
1018 }
1019 }
1020 }
1021}
1022
1023/* Workspace method: Doxygen documentation will be auto-generated */
1024void surfaceTelsem(Matrix& surface_los,
1025 Tensor4& surface_rmatrix,
1026 Matrix& surface_emission,
1027 const Index& atmosphere_dim,
1028 const Index& stokes_dim,
1029 const Vector& f_grid,
1030 const Vector& lat_grid,
1031 const Vector& lat_true,
1032 const Vector& lon_true,
1033 const Vector& rtp_pos,
1034 const Vector& rtp_los,
1035 const Vector& specular_los,
1036 const Numeric& surface_skin_t,
1037 const TelsemAtlas& atlas,
1038 const Numeric& r_min,
1039 const Numeric& r_max,
1040 const Numeric& d_max,
1041 const Verbosity& verbosity) {
1042 // Input checks
1043 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1044 chk_rte_pos(atmosphere_dim, rtp_pos);
1045 chk_rte_los(atmosphere_dim, rtp_los);
1046 chk_rte_los(atmosphere_dim, specular_los);
1048 "surface skin temperature", surface_skin_t, 190.0, 373.0);
1049
1050 const Index nf = f_grid.nelem();
1051 Matrix surface_rv_rh(nf, 2);
1052
1053 // Lookup cell in atlas.
1054
1055 Numeric lat, lon;
1057 lat, lon, atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
1058 chk_if_in_range("Latitude input to TELSEM2", lat, -90.0, 90.0);
1059
1060 lon = fmod(lon, 360.0);
1061 if (lon < 0) {
1062 lon += 360.0;
1063 }
1064 chk_if_in_range("Longitude input to TELSEM2", lon, 0.0, 360.0);
1065
1066 Index cellnumber = atlas.calc_cellnum(lat, lon);
1067 // Check if cell is in atlas.
1068 if (!atlas.contains(cellnumber)) {
1069 if (d_max <= 0.0) {
1070 throw std::runtime_error(
1071 "Given coordinates are not contained in "
1072 " TELSEM atlas. To enable nearest neighbour"
1073 "interpolation set *d_max* to a positive "
1074 "value.");
1075 } else {
1076 cellnumber = atlas.calc_cellnum_nearest_neighbor(lat, lon);
1077 Numeric lat_nn, lon_nn;
1078 std::tie(lat_nn, lon_nn) = atlas.get_coordinates(cellnumber);
1079 Numeric d = sphdist(lat, lon, lat_nn, lon_nn);
1080 if (d > d_max) {
1081 std::ostringstream out{};
1082 out << "Distance of nearest neighbour exceeds provided limit (";
1083 out << d << " > " << d_max << ").";
1084 throw std::runtime_error(out.str());
1085 }
1086 }
1087 }
1088
1089 Index class1 = atlas.get_class1(cellnumber);
1090 Index class2 = atlas.get_class2(cellnumber);
1091
1092 Vector emis_h = atlas.get_emis_h(cellnumber);
1093 Vector emis_v = atlas.get_emis_v(cellnumber);
1094 Numeric e_v, e_h;
1095
1096 // Incidence angle
1097 const Numeric theta = calc_incang(rtp_los, specular_los);
1098 ARTS_ASSERT(theta <= 90);
1099
1100 for (Index i = 0; i < nf; ++i) {
1101 if (f_grid[i] < 5e9)
1102 throw std::runtime_error("Only frequency >= 5 GHz are allowed");
1103 if (f_grid[i] > 900e9)
1104 throw std::runtime_error("Only frequency <= 900 GHz are allowed");
1105
1106 // For frequencies 700 <= f <= 900 GHz use 700 GHz to avoid extrapolation
1107 // outside of TELSEM specifications.
1108 Numeric f = std::min(f_grid[i], 700e9) * 1e-9;
1109 std::tie(e_v, e_h) =
1110 atlas.emis_interp(theta, f, class1, class2, emis_v, emis_h);
1111
1112 surface_rv_rh(i, 0) = min(max(1.0 - e_v, r_min), r_max);
1113 surface_rv_rh(i, 1) = min(max(1.0 - e_h, r_min), r_max);
1114 }
1115
1116 surfaceFlatRvRh(surface_los,
1117 surface_rmatrix,
1118 surface_emission,
1119 f_grid,
1120 stokes_dim,
1121 atmosphere_dim,
1122 rtp_pos,
1123 rtp_los,
1124 specular_los,
1125 surface_skin_t,
1126 surface_rv_rh,
1127 verbosity);
1128}
1129
1130/* Workspace method: Doxygen documentation will be auto-generated */
1131void surfaceTessem(Matrix& surface_los,
1132 Tensor4& surface_rmatrix,
1133 Matrix& surface_emission,
1134 const Index& atmosphere_dim,
1135 const Index& stokes_dim,
1136 const Vector& f_grid,
1137 const Vector& rtp_pos,
1138 const Vector& rtp_los,
1139 const Numeric& surface_skin_t,
1140 const TessemNN& net_h,
1141 const TessemNN& net_v,
1142 const Numeric& salinity,
1143 const Numeric& wind_speed,
1144 const Verbosity& verbosity) {
1145 // Input checks
1146 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1147 chk_rte_pos(atmosphere_dim, rtp_pos);
1148 chk_rte_los(atmosphere_dim, rtp_los);
1150 "surface skin temperature", surface_skin_t, 260.0, 373.0);
1151 chk_if_in_range_exclude_high("salinity", salinity, 0, 1);
1152 chk_if_in_range_exclude_high("wind speed", wind_speed, 0, 100);
1153
1154 // Determine specular direction
1155 Vector specular_los, surface_normal;
1156 specular_losCalcNoTopography(specular_los,
1157 surface_normal,
1158 rtp_pos,
1159 rtp_los,
1160 atmosphere_dim,
1161 verbosity);
1162
1163 // TESSEM in and out
1164 //
1165 Vector out(2);
1166 VectorView e_h = out[Range(0, 1)];
1167 VectorView e_v = out[Range(1, 1)];
1168 //
1169 Vector in(5);
1170 in[1] = 180.0 - abs(rtp_los[0]);
1171 in[2] = wind_speed;
1172 in[3] = surface_skin_t;
1173 in[4] = salinity;
1174
1175 // Get Rv and Rh
1176 //
1177 const Index nf = f_grid.nelem();
1178 Matrix surface_rv_rh(nf, 2);
1179 //
1180 for (Index i = 0; i < nf; ++i) {
1181 if (f_grid[i] < 5e9)
1182 throw std::runtime_error("Only frequency >= 5 GHz are allowed");
1183 if (f_grid[i] > 900e9)
1184 throw std::runtime_error("Only frequency <= 900 GHz are allowed");
1185
1186 in[0] = f_grid[i];
1187
1188 tessem_prop_nn(e_h, net_h, in);
1189 tessem_prop_nn(e_v, net_v, in);
1190
1191 surface_rv_rh(i, 0) = min(max(1 - e_v[0], (Numeric)0), (Numeric)1);
1192 surface_rv_rh(i, 1) = min(max(1 - e_h[0], (Numeric)0), (Numeric)1);
1193 }
1194
1195 surfaceFlatRvRh(surface_los,
1196 surface_rmatrix,
1197 surface_emission,
1198 f_grid,
1199 stokes_dim,
1200 atmosphere_dim,
1201 rtp_pos,
1202 rtp_los,
1203 specular_los,
1204 surface_skin_t,
1205 surface_rv_rh,
1206 verbosity);
1207}
1208
1209/* Workspace method: Doxygen documentation will be auto-generated */
1211 Tensor4& surface_rmatrix,
1212 Matrix& surface_emission,
1213 const Vector& f_grid,
1214 const Index& stokes_dim,
1215 const Index& atmosphere_dim,
1216 const Vector& rtp_pos,
1217 const Vector& rtp_los,
1218 const Vector& specular_los,
1219 const Numeric& surface_skin_t,
1220 const GriddedField3& surface_complex_refr_index,
1221 const Verbosity& verbosity) {
1224
1225 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1226 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
1227 chk_rte_pos(atmosphere_dim, rtp_pos);
1228 chk_rte_los(atmosphere_dim, rtp_los);
1229 chk_rte_los(atmosphere_dim, specular_los);
1230 chk_not_negative("surface_skin_t", surface_skin_t);
1231
1232 // Interpolate *surface_complex_refr_index*
1233 //
1234 const Index nf = f_grid.nelem();
1235 //
1236 Matrix n_real(nf, 1), n_imag(nf, 1);
1237 //
1238 complex_n_interp(n_real,
1239 n_imag,
1240 surface_complex_refr_index,
1241 "surface_complex_refr_index",
1242 f_grid,
1243 Vector(1, surface_skin_t));
1244
1245 out2 << " Sets variables to model a flat surface\n";
1246 out3 << " surface temperature: " << surface_skin_t << " K.\n";
1247
1248 surface_los.resize(1, specular_los.nelem());
1249 surface_los(0, joker) = specular_los;
1250
1251 surface_emission.resize(nf, stokes_dim);
1252 surface_rmatrix.resize(1, nf, stokes_dim, stokes_dim);
1253
1254 // Incidence angle
1255 const Numeric incang = calc_incang(rtp_los, specular_los);
1256 ARTS_ASSERT(incang <= 90);
1257
1258 // Complex (amplitude) reflection coefficients
1259 Complex Rv, Rh;
1260
1261 for (Index iv = 0; iv < nf; iv++) {
1262 // Set n2 (refractive index of surface medium)
1263 Complex n2(n_real(iv, 0), n_imag(iv, 0));
1264
1265 // Amplitude reflection coefficients
1266 fresnel(Rv, Rh, Numeric(1.0), n2, incang);
1267
1268 // Fill reflection matrix and emission vector
1269 surface_specular_R_and_b(surface_rmatrix(0, iv, joker, joker),
1270 surface_emission(iv, joker),
1271 Rv,
1272 Rh,
1273 f_grid[iv],
1274 stokes_dim,
1275 surface_skin_t);
1276 }
1277}
1278
1279/* Workspace method: Doxygen documentation will be auto-generated */
1281 Tensor4& surface_rmatrix,
1282 Matrix& surface_emission,
1283 const Vector& f_grid,
1284 const Index& stokes_dim,
1285 const Index& atmosphere_dim,
1286 const Vector& rtp_pos,
1287 const Vector& rtp_los,
1288 const Vector& specular_los,
1289 const Numeric& surface_skin_t,
1290 const Tensor3& surface_reflectivity,
1291 const Verbosity& verbosity) {
1293
1294 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1295 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
1296 chk_rte_pos(atmosphere_dim, rtp_pos);
1297 chk_rte_los(atmosphere_dim, rtp_los);
1298 chk_rte_los(atmosphere_dim, specular_los);
1299 chk_not_negative("surface_skin_t", surface_skin_t);
1300
1301 const Index nf = f_grid.nelem();
1302
1303 if (surface_reflectivity.nrows() != stokes_dim &&
1304 surface_reflectivity.ncols() != stokes_dim) {
1305 ostringstream os;
1306 os << "The number of rows and columnss in *surface_reflectivity* must\n"
1307 << "match *stokes_dim*."
1308 << "\n stokes_dim : " << stokes_dim
1309 << "\n number of rows in *surface_reflectivity* : "
1310 << surface_reflectivity.nrows()
1311 << "\n number of columns in *surface_reflectivity* : "
1312 << surface_reflectivity.ncols() << "\n";
1313 throw runtime_error(os.str());
1314 }
1315
1316 if (surface_reflectivity.npages() != nf &&
1317 surface_reflectivity.npages() != 1) {
1318 ostringstream os;
1319 os << "The number of pages in *surface_reflectivity* should\n"
1320 << "match length of *f_grid* or be 1."
1321 << "\n length of *f_grid* : " << nf
1322 << "\n dimension of *surface_reflectivity* : "
1323 << surface_reflectivity.npages() << "\n";
1324 throw runtime_error(os.str());
1325 }
1326
1327 out2 << " Sets variables to model a flat surface\n";
1328
1329 surface_los.resize(1, specular_los.nelem());
1330 surface_los(0, joker) = specular_los;
1331
1332 surface_emission.resize(nf, stokes_dim);
1333 surface_rmatrix.resize(1, nf, stokes_dim, stokes_dim);
1334
1335 Matrix R, IR(stokes_dim, stokes_dim);
1336
1337 Vector b(nf);
1338 planck(b, f_grid, surface_skin_t);
1339
1340 Vector B(stokes_dim, 0);
1341
1342 for (Index iv = 0; iv < nf; iv++) {
1343 if (iv == 0 || surface_reflectivity.npages() > 1) {
1344 R = surface_reflectivity(iv, joker, joker);
1345 for (Index i = 0; i < stokes_dim; i++) {
1346 for (Index j = 0; j < stokes_dim; j++) {
1347 if (i == j) {
1348 IR(i, j) = 1 - R(i, j);
1349 } else {
1350 IR(i, j) = -R(i, j);
1351 }
1352 }
1353 }
1354 }
1355
1356 surface_rmatrix(0, iv, joker, joker) = R;
1357
1358 B[0] = b[iv];
1359 mult(surface_emission(iv, joker), IR, B);
1360 }
1361}
1362
1363/* Workspace method: Doxygen documentation will be auto-generated */
1364void surfaceFlatRvRh(Matrix& surface_los,
1365 Tensor4& surface_rmatrix,
1366 Matrix& surface_emission,
1367 const Vector& f_grid,
1368 const Index& stokes_dim,
1369 const Index& atmosphere_dim,
1370 const Vector& rtp_pos,
1371 const Vector& rtp_los,
1372 const Vector& specular_los,
1373 const Numeric& surface_skin_t,
1374 const Matrix& surface_rv_rh,
1375 const Verbosity&) {
1376 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1377 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
1378 chk_rte_pos(atmosphere_dim, rtp_pos);
1379 chk_rte_los(atmosphere_dim, rtp_los);
1380 chk_rte_los(atmosphere_dim, specular_los);
1381 chk_not_negative("surface_skin_t", surface_skin_t);
1382
1383 const Index nf = f_grid.nelem();
1384
1385 if (surface_rv_rh.ncols() != 2) {
1386 ostringstream os;
1387 os << "The number of columns in *surface_rv_rh* must be two,\n"
1388 << "but the actual number of columns is " << surface_rv_rh.ncols()
1389 << "\n";
1390 throw runtime_error(os.str());
1391 }
1392
1393 if (surface_rv_rh.nrows() != nf && surface_rv_rh.nrows() != 1) {
1394 ostringstream os;
1395 os << "The number of rows in *surface_rv_rh* should\n"
1396 << "match length of *f_grid* or be 1."
1397 << "\n length of *f_grid* : " << nf
1398 << "\n rows in *surface_rv_rh* : " << surface_rv_rh.nrows() << "\n";
1399 throw runtime_error(os.str());
1400 }
1401
1402 if (min(surface_rv_rh) < 0 || max(surface_rv_rh) > 1) {
1403 throw runtime_error("All values in *surface_rv_rh* must be inside [0,1].");
1404 }
1405
1406 surface_los.resize(1, specular_los.nelem());
1407 surface_los(0, joker) = specular_los;
1408
1409 surface_emission.resize(nf, stokes_dim);
1410 surface_rmatrix.resize(1, nf, stokes_dim, stokes_dim);
1411
1412 surface_emission = 0;
1413 surface_rmatrix = 0;
1414
1415 Vector b(nf);
1416 planck(b, f_grid, surface_skin_t);
1417
1418 Numeric rmean = 0.0, rdiff = 0.0;
1419
1420 for (Index iv = 0; iv < nf; iv++) {
1421 if (iv == 0 || surface_rv_rh.nrows() > 1) {
1422 rmean = 0.5 * (surface_rv_rh(iv, 0) + surface_rv_rh(iv, 1));
1423 rdiff = 0.5 * (surface_rv_rh(iv, 0) - surface_rv_rh(iv, 1));
1424 }
1425
1426 surface_emission(iv, 0) = (1.0 - rmean) * b[iv];
1427 surface_rmatrix(0, iv, 0, 0) = rmean;
1428
1429 if (stokes_dim > 1) {
1430 surface_emission(iv, 1) = -rdiff * b[iv];
1431
1432 surface_rmatrix(0, iv, 0, 1) = rdiff;
1433 surface_rmatrix(0, iv, 1, 0) = rdiff;
1434 surface_rmatrix(0, iv, 1, 1) = rmean;
1435
1436 for (Index i = 2; i < stokes_dim; i++) {
1437 surface_rmatrix(0, iv, i, i) = rmean;
1438 }
1439 }
1440 }
1441}
1442
1443/* Workspace method: Doxygen documentation will be auto-generated */
1445 Tensor4& surface_rmatrix,
1446 Matrix& surface_emission,
1447 const Vector& f_grid,
1448 const Index& stokes_dim,
1449 const Index& atmosphere_dim,
1450 const Vector& rtp_pos,
1451 const Vector& rtp_los,
1452 const Vector& specular_los,
1453 const Numeric& surface_skin_t,
1454 const Vector& surface_scalar_reflectivity,
1455 const Verbosity&) {
1456 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1457 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
1458 chk_rte_pos(atmosphere_dim, rtp_pos);
1459 chk_rte_los(atmosphere_dim, rtp_los);
1460 chk_rte_los(atmosphere_dim, specular_los);
1461 chk_not_negative("surface_skin_t", surface_skin_t);
1462
1463 const Index nf = f_grid.nelem();
1464
1465 if (surface_scalar_reflectivity.nelem() != nf &&
1466 surface_scalar_reflectivity.nelem() != 1) {
1467 ostringstream os;
1468 os << "The number of elements in *surface_scalar_reflectivity* should\n"
1469 << "match length of *f_grid* or be 1."
1470 << "\n length of *f_grid* : " << nf
1471 << "\n length of *surface_scalar_reflectivity* : "
1472 << surface_scalar_reflectivity.nelem() << "\n";
1473 throw runtime_error(os.str());
1474 }
1475
1476 if (min(surface_scalar_reflectivity) < 0 ||
1477 max(surface_scalar_reflectivity) > 1) {
1478 throw runtime_error(
1479 "All values in *surface_scalar_reflectivity* must be inside [0,1].");
1480 }
1481
1482 surface_los.resize(1, specular_los.nelem());
1483 surface_los(0, joker) = specular_los;
1484
1485 surface_emission.resize(nf, stokes_dim);
1486 surface_rmatrix.resize(1, nf, stokes_dim, stokes_dim);
1487
1488 surface_emission = 0;
1489 surface_rmatrix = 0;
1490
1491 Vector b(nf);
1492 planck(b, f_grid, surface_skin_t);
1493
1494 Numeric r = 0.0;
1495
1496 for (Index iv = 0; iv < nf; iv++) {
1497 if (iv == 0 || surface_scalar_reflectivity.nelem() > 1) {
1498 r = surface_scalar_reflectivity[iv];
1499 }
1500
1501 surface_emission(iv, 0) = (1.0 - r) * b[iv];
1502 surface_rmatrix(0, iv, 0, 0) = r;
1503 for (Index i = 1; i < stokes_dim; i++) {
1504 surface_rmatrix(0, iv, i, i) = r;
1505 }
1506 }
1507}
1508
1509/* Workspace method: Doxygen documentation will be auto-generated */
1511 Tensor4& surface_rmatrix,
1512 Matrix& surface_emission,
1513 const Vector& f_grid,
1514 const Index& stokes_dim,
1515 const Index& atmosphere_dim,
1516 const Vector& rtp_pos,
1517 const Vector& rtp_los,
1518 const Vector& surface_normal,
1519 const Numeric& surface_skin_t,
1520 const Vector& surface_scalar_reflectivity,
1521 const Index& lambertian_nza,
1522 const Numeric& za_pos,
1523 const Verbosity&) {
1524 const Index nf = f_grid.nelem();
1525
1526 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1527 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
1528 chk_rte_pos(atmosphere_dim, rtp_pos);
1529 chk_rte_los(atmosphere_dim, rtp_los);
1530 chk_not_negative("surface_skin_t", surface_skin_t);
1531 chk_if_in_range("za_pos", za_pos, 0, 1);
1532
1533 if (surface_scalar_reflectivity.nelem() != nf &&
1534 surface_scalar_reflectivity.nelem() != 1) {
1535 ostringstream os;
1536 os << "The number of elements in *surface_scalar_reflectivity* should\n"
1537 << "match length of *f_grid* or be 1."
1538 << "\n length of *f_grid* : " << nf
1539 << "\n length of *surface_scalar_reflectivity* : "
1540 << surface_scalar_reflectivity.nelem() << "\n";
1541 throw runtime_error(os.str());
1542 }
1543
1544 if (min(surface_scalar_reflectivity) < 0 ||
1545 max(surface_scalar_reflectivity) > 1) {
1546 throw runtime_error(
1547 "All values in *surface_scalar_reflectivity* must be inside [0,1].");
1548 }
1549
1550 // Allocate and init everything to zero
1551 //
1552 surface_los.resize(lambertian_nza, rtp_los.nelem());
1553 surface_rmatrix.resize(lambertian_nza, nf, stokes_dim, stokes_dim);
1554 surface_emission.resize(nf, stokes_dim);
1555 //
1556 surface_los = 0.0;
1557 surface_rmatrix = 0.0;
1558 surface_emission = 0.0;
1559
1560 // Help variables
1561 //
1562 const Numeric dza = (90.0 - abs(surface_normal[0])) / (Numeric)lambertian_nza;
1563 const Vector za_lims(0.0, lambertian_nza + 1, dza);
1564
1565 // surface_los
1566 for (Index ip = 0; ip < lambertian_nza; ip++) {
1567 surface_los(ip, 0) = za_lims[ip] + za_pos * dza;
1568 if (atmosphere_dim == 2) {
1569 if (rtp_los[0] < 0) {
1570 surface_los(ip, 0) *= -1.0;
1571 }
1572
1573 } else if (atmosphere_dim == 3) {
1574 surface_los(ip, 1) = rtp_los[1];
1575 }
1576 }
1577
1578 Vector b(nf);
1579 planck(b, f_grid, surface_skin_t);
1580
1581 // Loop frequencies and set remaining values
1582 //
1583 Numeric r = 0.0;
1584 //
1585 for (Index iv = 0; iv < nf; iv++) {
1586 // Get reflectivity
1587 if (iv == 0 || surface_scalar_reflectivity.nelem() > 1) {
1588 r = surface_scalar_reflectivity[iv];
1589 }
1590
1591 // surface_rmatrix:
1592 // Only element (0,0) is set to be non-zero. This follows VDISORT
1593 // that refers to: K. L. Coulson, Polarization and Intensity of Light in
1594 // the Atmosphere (1989), page 229
1595 // (Thanks to Michael Kahnert for providing this information!)
1596 // Update: Is the above for a later edition? We have not found a copy of
1597 // that edition. In a 1988 version of the book, the relevant page seems
1598 // to be 232.
1599 for (Index ip = 0; ip < lambertian_nza; ip++) {
1600 const Numeric w =
1601 r * 0.5 *
1602 (cos(2 * DEG2RAD * za_lims[ip]) - cos(2 * DEG2RAD * za_lims[ip + 1]));
1603 surface_rmatrix(ip, iv, 0, 0) = w;
1604 }
1605
1606 // surface_emission
1607 surface_emission(iv, 0) = (1 - r) * b[iv];
1608 }
1609}
1610
1611/* Workspace method: Doxygen documentation will be auto-generated */
1613 Numeric& surface_skin_t,
1614 Matrix& surface_los,
1615 Tensor4& surface_rmatrix,
1616 Matrix& surface_emission,
1617 const Index& atmosphere_dim,
1618 const Vector& f_grid,
1619 const Vector& rtp_pos,
1620 const Vector& rtp_los,
1621 const Agenda& surface_rtprop_sub_agenda,
1622 const Numeric& specular_factor,
1623 const Numeric& dza,
1624 const Verbosity&) {
1625 chk_rte_pos(atmosphere_dim, rtp_pos);
1626 chk_rte_los(atmosphere_dim, rtp_los);
1627
1628 // Checks of GIN variables
1629 if (specular_factor > 1 || specular_factor < 1.0 / 3.0)
1630 throw runtime_error("The valid range for *specular_factor* is [1/3,1].");
1631 if (dza > 45 || dza <= 0)
1632 throw runtime_error("The valid range for *dza* is ]0,45].");
1633
1634 // Obtain data for specular direction
1635 //
1636 Matrix los1, emission1;
1637 Tensor4 rmatrix1;
1638 //
1640 surface_skin_t,
1641 emission1,
1642 los1,
1643 rmatrix1,
1644 f_grid,
1645 rtp_pos,
1646 rtp_los,
1647 surface_rtprop_sub_agenda);
1648 if (los1.nrows() != 1)
1649 throw runtime_error(
1650 "*surface_rtprop_sub_agenda* must return data "
1651 "describing a specular surface.");
1652
1653 // Standard number of beams. Set to 2 if try/catch below fails
1654 Index nbeams = 3;
1655
1656 // Test if some lower zenith angle works.
1657 // It will fail if a higher za results in looking at the surface from below.
1658 //
1659 Matrix los2, emission2;
1660 Tensor4 rmatrix2;
1661 //
1662 Numeric skin_t_dummy;
1663 Numeric dza_try = dza;
1664 bool failed = true;
1665 while (failed && dza_try > 0) {
1666 try {
1667 Vector los_new = rtp_los;
1668 los_new[0] -=
1669 sign(rtp_los[0]) * dza_try; // Sign to also handle 2D negative za
1670 adjust_los(los_new, atmosphere_dim);
1672 skin_t_dummy,
1673 emission2,
1674 los2,
1675 rmatrix2,
1676 f_grid,
1677 rtp_pos,
1678 los_new,
1679 surface_rtprop_sub_agenda);
1680 failed = false;
1681 } catch (const runtime_error& e) {
1682 dza_try -= 1.0;
1683 }
1684 }
1685 if (failed) {
1686 nbeams = 2;
1687 }
1688
1689 // Allocate output WSVs
1690 //
1691 surface_emission.resize(emission1.nrows(), emission1.ncols());
1692 surface_emission = 0;
1693 surface_los.resize(nbeams, los1.ncols());
1694 surface_rmatrix.resize(
1695 nbeams, rmatrix1.npages(), rmatrix1.nrows(), rmatrix1.ncols());
1696
1697 // Put in specular direction at index 1
1698 //
1699 Numeric w;
1700 if (nbeams == 3) {
1701 w = specular_factor;
1702 } else {
1703 w = specular_factor + (1.0 - specular_factor) / 2.0;
1704 }
1705 //
1706 surface_los(1, joker) = los1(0, joker);
1707 for (Index p = 0; p < rmatrix1.npages(); p++) {
1708 for (Index r = 0; r < rmatrix1.nrows(); r++) {
1709 surface_emission(p, r) += w * emission1(p, r);
1710 for (Index c = 0; c < rmatrix1.ncols(); c++) {
1711 surface_rmatrix(1, p, r, c) = w * rmatrix1(0, p, r, c);
1712 }
1713 }
1714 }
1715
1716 // Put in lower za as index 2, if worked
1717 //
1718 w = (1.0 - specular_factor) / 2.0;
1719 //
1720 if (nbeams == 3) {
1721 surface_los(2, joker) = los2(0, joker);
1722 for (Index p = 0; p < rmatrix2.npages(); p++) {
1723 for (Index r = 0; r < rmatrix2.nrows(); r++) {
1724 surface_emission(p, r) += w * emission2(p, r);
1725 for (Index c = 0; c < rmatrix1.ncols(); c++) {
1726 surface_rmatrix(2, p, r, c) = w * rmatrix2(0, p, r, c);
1727 }
1728 }
1729 }
1730 }
1731
1732 // Do higher za and put in as index 0 (reusing variables for beam 2)
1733 //
1734 Vector los_new = rtp_los;
1735 los_new[0] += sign(rtp_los[0]) * dza; // Sign to also handle 2D negative za
1736 adjust_los(los_new, atmosphere_dim);
1738 skin_t_dummy,
1739 emission2,
1740 los2,
1741 rmatrix2,
1742 f_grid,
1743 rtp_pos,
1744 los_new,
1745 surface_rtprop_sub_agenda);
1746 //
1747 surface_los(0, joker) = los2(0, joker);
1748 for (Index p = 0; p < rmatrix2.npages(); p++) {
1749 for (Index r = 0; r < rmatrix2.nrows(); r++) {
1750 surface_emission(p, r) += w * emission2(p, r);
1751 for (Index c = 0; c < rmatrix1.ncols(); c++) {
1752 surface_rmatrix(0, p, r, c) = w * rmatrix2(0, p, r, c);
1753 }
1754 }
1755 }
1756}
1757
1758/* Workspace method: Doxygen documentation will be auto-generated */
1760 Tensor4& surface_rmatrix,
1761 const Index& atmosphere_dim,
1762 const Vector& rtp_pos,
1763 const Vector& rtp_los,
1764 const Numeric& specular_factor,
1765 const Numeric& dza,
1766 const Verbosity&) {
1767 chk_rte_pos(atmosphere_dim, rtp_pos);
1768 chk_rte_los(atmosphere_dim, rtp_los);
1769
1770 // Check that input surface data are of specular type
1771 if (surface_los.nrows() != 1)
1772 throw runtime_error(
1773 "Input surface data must be of specular type. That is, "
1774 "*surface_los* must contain a single direction.");
1775 if (surface_rmatrix.nbooks() != 1)
1776 throw runtime_error(
1777 "*surface_rmatrix* describes a different number of "
1778 "directions than *surface_los*.");
1779
1780 // Checks of GIN variables
1781 if (specular_factor > 1 || specular_factor < 1.0 / 3.0)
1782 throw runtime_error("The valid range for *specular_factor* is [1/3,1].");
1783 if (dza > 45 || dza <= 0)
1784 throw runtime_error("The valid range for *dza* is ]0,45].");
1785
1786 // Make copies of input data
1787 const Matrix los1 = surface_los;
1788 const Tensor4 rmatrix1 = surface_rmatrix;
1789
1790 // Use abs(za) in all expressions below, to also handle 2D
1791
1792 // Calculate highest possible za for downwelling radiation, with 1 degree
1793 // margin to the surface. This can be derived from za in rtp_los and surface_los.
1794 // (The directions in surface_los are not allowed to point into the surface)
1795 const Numeric za_max = 89 + (180 - abs(los1(0, 0)) - abs(rtp_los[0])) / 2.0;
1796
1797 // Number of downwelling beams
1798 Index nbeams = 3;
1799 if (abs(los1(0, 0)) > za_max) {
1800 nbeams = 2;
1801 }
1802
1803 // New los-s
1804 //
1805 surface_los.resize(nbeams, los1.ncols());
1806 //
1807 for (Index r = 0; r < nbeams; r++) {
1808 surface_los(r, 0) = ((Numeric)r - 1.0) * dza + abs(los1(0, 0));
1809 if (r == 2 && surface_los(r, 0) > za_max) {
1810 surface_los(r, 0) = za_max;
1811 }
1812 for (Index c = 1; c < los1.ncols(); c++) {
1813 surface_los(r, c) = los1(0, c);
1814 }
1815 }
1816
1817 // New rmatrix
1818 //
1819 surface_rmatrix.resize(
1820 nbeams, rmatrix1.npages(), rmatrix1.nrows(), rmatrix1.ncols());
1821 //
1822 for (Index b = 0; b < nbeams; b++) {
1823 Numeric w;
1824 if (b == 1 && nbeams == 3) // Specular direction with nbeams==3
1825 {
1826 w = specular_factor;
1827 } else if (b == 1) // Specular direction with nbeams==2
1828 {
1829 w = specular_factor + (1.0 - specular_factor) / 2.0;
1830 } else // Side directions
1831 {
1832 w = (1.0 - specular_factor) / 2.0;
1833 }
1834
1835 for (Index p = 0; p < rmatrix1.npages(); p++) {
1836 for (Index r = 0; r < rmatrix1.nrows(); r++) {
1837 for (Index c = 0; c < rmatrix1.ncols(); c++) {
1838 surface_rmatrix(b, p, r, c) = w * rmatrix1(0, p, r, c);
1839 }
1840 }
1841 }
1842 }
1843
1844 // Handle sign of za
1845 if (atmosphere_dim == 1) {
1846 // We only need to make sure that first direction has positive za
1847 surface_los(0, 0) = abs(surface_los(0, 0));
1848 } else if (atmosphere_dim == 2) {
1849 // Change sign if specular direction has za < 0
1850 if (los1(0, 0) < 0) {
1851 for (Index r = 0; r < rmatrix1.nrows(); r++) {
1852 surface_los(r, 0) = -surface_los(r, 0);
1853 }
1854 }
1855 } else if (atmosphere_dim == 1) {
1856 // We only need to make sure that first direction has positive za
1857 if (surface_los(0, 0) < 0) {
1858 surface_los(0, 0) = -surface_los(0, 0);
1859 surface_los(0, 1) += 180;
1860 if (surface_los(0, 1) > 180) {
1861 surface_los(0, 1) -= 360;
1862 }
1863 }
1864 }
1865}
1866
1867/* Workspace method: Doxygen documentation will be auto-generated */
1869 GriddedField3& surface_complex_refr_index,
1870 const Index& atmosphere_dim,
1871 const Vector& lat_grid,
1872 const Vector& lat_true,
1873 const Vector& lon_true,
1874 const Vector& rtp_pos,
1875 const GriddedField5& complex_n_field,
1876 const Verbosity&) {
1877 // Set expected order of grids
1878 Index gfield_fID = 0;
1879 Index gfield_tID = 1;
1880 Index gfield_compID = 2;
1881 Index gfield_latID = 3;
1882 Index gfield_lonID = 4;
1883
1884 // Basic checks and sizes
1885 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1886 chk_latlon_true(atmosphere_dim, lat_grid, lat_true, lon_true);
1887 chk_rte_pos(atmosphere_dim, rtp_pos);
1888 complex_n_field.checksize_strict();
1889 //
1890 chk_griddedfield_gridname(complex_n_field, gfield_fID, "Frequency");
1891 chk_griddedfield_gridname(complex_n_field, gfield_tID, "Temperature");
1892 chk_griddedfield_gridname(complex_n_field, gfield_compID, "Complex");
1893 chk_griddedfield_gridname(complex_n_field, gfield_latID, "Latitude");
1894 chk_griddedfield_gridname(complex_n_field, gfield_lonID, "Longitude");
1895 //
1896 const Index nf = complex_n_field.data.nshelves();
1897 const Index nt = complex_n_field.data.nbooks();
1898 const Index nn = complex_n_field.data.npages();
1899 const Index nlat = complex_n_field.data.nrows();
1900 const Index nlon = complex_n_field.data.ncols();
1901 //
1902 if (nlat < 2 || nlon < 2) {
1903 ostringstream os;
1904 os << "The data in *complex_refr_index_field* must span a geographical "
1905 << "region. That is,\nthe latitude and longitude grids must have a "
1906 << "length >= 2.";
1907 throw runtime_error(os.str());
1908 }
1909 //
1910 if (nn != 2) {
1911 ostringstream os;
1912 os << "The data in *complex_refr_index_field* must have exactly two "
1913 << "pages. One page each\nfor the real and imaginary part of the "
1914 << "complex refractive index.";
1915 throw runtime_error(os.str());
1916 }
1917
1918 const Vector& GFlat = complex_n_field.get_numeric_grid(gfield_latID);
1919 const Vector& GFlon = complex_n_field.get_numeric_grid(gfield_lonID);
1920
1921 // Determine true geographical position
1922 Vector lat(1), lon(1);
1924 lat[0], lon[0], atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
1925
1926 // Ensure correct coverage of lon grid
1927 Vector lon_shifted;
1928 lon_shiftgrid(lon_shifted, GFlon, lon[0]);
1929
1930 // Check if lat/lon we need are actually covered
1931 chk_if_in_range("rtp_pos.lat", lat[0], GFlat[0], GFlat[nlat - 1]);
1932 chk_if_in_range("rtp_pos.lon", lon[0], lon_shifted[0], lon_shifted[nlon - 1]);
1933
1934 // Size and fills grids of *surface_complex_refr_index*
1935 surface_complex_refr_index.resize(nf, nt, 2);
1936 surface_complex_refr_index.set_grid_name(0, "Frequency");
1937 surface_complex_refr_index.set_grid(
1938 0, complex_n_field.get_numeric_grid(gfield_fID));
1939 surface_complex_refr_index.set_grid_name(1, "Temperature");
1940 surface_complex_refr_index.set_grid(
1941 1, complex_n_field.get_numeric_grid(gfield_tID));
1942 surface_complex_refr_index.set_grid_name(2, "Complex");
1943 surface_complex_refr_index.set_grid(2, {"real", "imaginary"});
1944
1945 // Interpolate in lat and lon
1946 //
1947 GridPos gp_lat, gp_lon;
1948 gridpos(gp_lat, GFlat, lat[0]);
1949 gridpos(gp_lon, lon_shifted, lon[0]);
1950 Vector itw(4);
1951 interpweights(itw, gp_lat, gp_lon);
1952 //
1953 for (Index iv = 0; iv < nf; iv++) {
1954 for (Index it = 0; it < nt; it++) {
1955 surface_complex_refr_index.data(iv, it, 0) = interp(
1956 itw, complex_n_field.data(iv, it, 0, joker, joker), gp_lat, gp_lon);
1957 surface_complex_refr_index.data(iv, it, 1) = interp(
1958 itw, complex_n_field.data(iv, it, 1, joker, joker), gp_lat, gp_lon);
1959 }
1960 }
1961}
1962
1963/* Workspace method: Doxygen documentation will be auto-generated */
1965 const Index& stokes_dim,
1966 const Vector& f_grid,
1967 const Index& atmosphere_dim,
1968 const Vector& lat_grid,
1969 const Vector& lat_true,
1970 const Vector& lon_true,
1971 const Vector& rtp_pos,
1972 const Vector& rtp_los,
1973 const GriddedField6& r_field,
1974 const Verbosity&) {
1975 // Basic checks and sizes
1976 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1977 chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
1978 chk_latlon_true(atmosphere_dim, lat_grid, lat_true, lon_true);
1979 chk_rte_pos(atmosphere_dim, rtp_pos);
1980 chk_rte_los(atmosphere_dim, rtp_los);
1981 r_field.checksize_strict();
1982 chk_griddedfield_gridname(r_field, 0, "Frequency");
1983 chk_griddedfield_gridname(r_field, 1, "Stokes element");
1984 chk_griddedfield_gridname(r_field, 2, "Stokes element");
1985 chk_griddedfield_gridname(r_field, 3, "Incidence angle");
1986 chk_griddedfield_gridname(r_field, 4, "Latitude");
1987 chk_griddedfield_gridname(r_field, 5, "Longitude");
1988 //
1989 const Index nf_in = r_field.data.nvitrines();
1990 const Index ns2 = r_field.data.nshelves();
1991 const Index ns1 = r_field.data.nbooks();
1992 const Index nza = r_field.data.npages();
1993 const Index nlat = r_field.data.nrows();
1994 const Index nlon = r_field.data.ncols();
1995 //
1996 if (nlat < 2 || nlon < 2) {
1997 ostringstream os;
1998 os << "The data in *r_field* must span a geographical region. That is,\n"
1999 << "the latitude and longitude grids must have a length >= 2.";
2000 throw runtime_error(os.str());
2001 }
2002 //
2003 if (nza < 2) {
2004 ostringstream os;
2005 os << "The data in *r_field* must span a range of zenith angles. That\n"
2006 << "is the zenith angle grid must have a length >= 2.";
2007 throw runtime_error(os.str());
2008 }
2009 if (ns1 < stokes_dim || ns2 < stokes_dim || ns1 > 4 || ns2 > 4) {
2010 ostringstream os;
2011 os << "The \"Stokes dimensions\" must have a size that is >= "
2012 << "*stokes_dim* (but not exceeding 4).";
2013 throw runtime_error(os.str());
2014 }
2015
2016 // Determine true geographical position
2017 Vector lat(1), lon(1);
2019 lat[0], lon[0], atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
2020
2021 // Ensure correct coverage of lon grid
2022 Vector lon_shifted;
2023 lon_shiftgrid(lon_shifted, r_field.get_numeric_grid(5), lon[0]);
2024
2025 // Interpolate in lat and lon
2026 Tensor4 r_f_za(nf_in, stokes_dim, stokes_dim, nza);
2027 {
2029 "Latitude interpolation", r_field.get_numeric_grid(4), lat[0]);
2030 chk_interpolation_grids("Longitude interpolation", lon_shifted, lon[0]);
2031 GridPos gp_lat, gp_lon;
2032 gridpos(gp_lat, r_field.get_numeric_grid(4), lat[0]);
2033 gridpos(gp_lon, lon_shifted, lon[0]);
2034 Vector itw(4);
2035 interpweights(itw, gp_lat, gp_lon);
2036 for (Index iv = 0; iv < nf_in; iv++) {
2037 for (Index iz = 0; iz < nza; iz++) {
2038 for (Index is1 = 0; is1 < stokes_dim; is1++) {
2039 for (Index is2 = 0; is2 < stokes_dim; is2++) {
2040 r_f_za(iv, is1, is2, iz) =
2041 interp(itw,
2042 r_field.data(iv, is1, is2, iz, joker, joker),
2043 gp_lat,
2044 gp_lon);
2045 }
2046 }
2047 }
2048 }
2049 }
2050
2051 // Interpolate in incidence angle, cubic if possible
2052 Tensor3 r_f(nf_in, stokes_dim, stokes_dim);
2053 Index order = 3;
2054 if (nza < 4) {
2055 order = 1;
2056 }
2057 {
2059 "Incidence angle interpolation", r_field.get_numeric_grid(3), 180 - rtp_los[0]);
2060 const LagrangeInterpolation lag(0, 180 - rtp_los[0], r_field.get_numeric_grid(3), order);
2061 const auto itw = interpweights(lag);
2062 //
2063 for (Index i = 0; i < nf_in; i++) {
2064 for (Index is1 = 0; is1 < stokes_dim; is1++) {
2065 for (Index is2 = 0; is2 < stokes_dim; is2++) {
2066 r_f(i, is1, is2) = interp(r_f_za(i, is1, is2, joker), itw, lag);
2067 }
2068 }
2069 }
2070 }
2071
2072 // Extract or interpolate in frequency
2073 //
2074 if (nf_in == 1) {
2075 surface_reflectivity = r_f;
2076 } else {
2078 "Frequency interpolation", r_field.get_numeric_grid(0), f_grid);
2079 const Index nf_out = f_grid.nelem();
2080 surface_reflectivity.resize(nf_out, stokes_dim, stokes_dim);
2081 //
2082 ArrayOfGridPos gp(nf_out);
2083 Matrix itw(nf_out, 2);
2084 gridpos(gp, r_field.get_numeric_grid(0), f_grid);
2085 interpweights(itw, gp);
2086 for (Index is1 = 0; is1 < stokes_dim; is1++) {
2087 for (Index is2 = 0; is2 < stokes_dim; is2++) {
2088 interp(surface_reflectivity(joker, is1, is2),
2089 itw,
2090 r_f(joker, is1, is2),
2091 gp);
2092 }
2093 }
2094 }
2095}
2096
2097/* Workspace method: Doxygen documentation will be auto-generated */
2099 Vector& surface_scalar_reflectivity,
2100 const Index& stokes_dim,
2101 const Vector& f_grid,
2102 const Index& atmosphere_dim,
2103 const Vector& lat_grid,
2104 const Vector& lat_true,
2105 const Vector& lon_true,
2106 const Vector& rtp_pos,
2107 const Vector& rtp_los,
2108 const GriddedField4& r_field,
2109 const Verbosity&) {
2110 // Basic checks and sizes
2111 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2112 chk_if_in_range("stokes_dim", stokes_dim, 1, 1);
2113 chk_latlon_true(atmosphere_dim, lat_grid, lat_true, lon_true);
2114 chk_rte_pos(atmosphere_dim, rtp_pos);
2115 chk_rte_los(atmosphere_dim, rtp_los);
2116 r_field.checksize_strict();
2117 chk_griddedfield_gridname(r_field, 0, "Frequency");
2118 chk_griddedfield_gridname(r_field, 1, "Incidence angle");
2119 chk_griddedfield_gridname(r_field, 2, "Latitude");
2120 chk_griddedfield_gridname(r_field, 3, "Longitude");
2121 //
2122 const Index nf_in = r_field.data.nbooks();
2123 const Index nza = r_field.data.npages();
2124 const Index nlat = r_field.data.nrows();
2125 const Index nlon = r_field.data.ncols();
2126 //
2127 if (nlat < 2 || nlon < 2) {
2128 ostringstream os;
2129 os << "The data in *r_field* must span a geographical region. That is,\n"
2130 << "the latitude and longitude grids must have a length >= 2.";
2131 throw runtime_error(os.str());
2132 }
2133 //
2134 if (nza < 2) {
2135 ostringstream os;
2136 os << "The data in *r_field* must span a range of zenith angles. That\n"
2137 << "is the zenith angle grid must have a length >= 2.";
2138 throw runtime_error(os.str());
2139 }
2140
2141 // Determine true geographical position
2142 Vector lat(1), lon(1);
2144 lat[0], lon[0], atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
2145
2146 // Ensure correct coverage of lon grid
2147 Vector lon_shifted;
2148 lon_shiftgrid(lon_shifted, r_field.get_numeric_grid(3), lon[0]);
2149
2150 // Interpolate in lat and lon
2151 Matrix r_f_za(nf_in, nza);
2152 {
2154 "Latitude interpolation", r_field.get_numeric_grid(2), lat[0]);
2155 chk_interpolation_grids("Longitude interpolation", lon_shifted, lon[0]);
2156 GridPos gp_lat, gp_lon;
2157 gridpos(gp_lat, r_field.get_numeric_grid(2), lat[0]);
2158 gridpos(gp_lon, lon_shifted, lon[0]);
2159 Vector itw(4);
2160 interpweights(itw, gp_lat, gp_lon);
2161 for (Index iv = 0; iv < nf_in; iv++) {
2162 for (Index iz = 0; iz < nza; iz++) {
2163 r_f_za(iv, iz) =
2164 interp(itw, r_field.data(iv, iz, joker, joker), gp_lat, gp_lon);
2165 }
2166 }
2167 }
2168
2169 // Interpolate in incidence angle, cubic if possible
2170 Vector r_f(nf_in);
2171 Index order = 3;
2172 if (nza < 4) {
2173 order = 1;
2174 }
2175 {
2177 "Incidence angle interpolation", r_field.get_numeric_grid(1), 180 - rtp_los[0]);
2178 const LagrangeInterpolation lag(0, 180 - rtp_los[0], r_field.get_numeric_grid(1), order);
2179 const auto itw=interpweights(lag);
2180 for (Index i = 0; i < nf_in; i++) {
2181 r_f[i] = interp(r_f_za(i, joker), itw, lag);
2182 }
2183 }
2184
2185 // Extract or interpolate in frequency
2186 //
2187 if (nf_in == 1) {
2188 surface_scalar_reflectivity.resize(1);
2189 surface_scalar_reflectivity[0] = r_f[0];
2190 } else {
2192 "Frequency interpolation", r_field.get_numeric_grid(0), f_grid);
2193 const Index nf_out = f_grid.nelem();
2194 surface_scalar_reflectivity.resize(nf_out);
2195 //
2196 ArrayOfGridPos gp(nf_out);
2197 Matrix itw(nf_out, 2);
2198 gridpos(gp, r_field.get_numeric_grid(0), f_grid);
2199 interpweights(itw, gp);
2200 interp(surface_scalar_reflectivity, itw, r_f, gp);
2201 }
2202}
2203
2204/* Workspace method: Doxygen documentation will be auto-generated */
2206 Vector& surface_scalar_reflectivity,
2207 const Tensor4& surface_rmatrix,
2208 const Verbosity&) {
2209 const Index nf = surface_rmatrix.npages();
2210 const Index nlos = surface_rmatrix.nbooks();
2211
2212 surface_scalar_reflectivity.resize(nf);
2213 surface_scalar_reflectivity = 0;
2214
2215 for (Index i = 0; i < nf; i++) {
2216 for (Index l = 0; l < nlos; l++) {
2217 surface_scalar_reflectivity[i] += surface_rmatrix(l, i, 0, 0);
2218 }
2219 }
2220}
2221
2222/* Workspace method: Doxygen documentation will be auto-generated */
2224 Vector& surface_types_aux,
2225 Vector& surface_types_weights,
2226 const Index& atmosphere_dim,
2227 const Vector& lat_grid,
2228 const Vector& lat_true,
2229 const Vector& lon_true,
2230 const Vector& rtp_pos,
2231 const GriddedField2& surface_type_mask,
2232 const String& method,
2233 const Verbosity&) {
2234 // Set expected order of grids
2235 Index gfield_latID = 0;
2236 Index gfield_lonID = 1;
2237
2238 // Basic checks and sizes
2239 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
2240 chk_latlon_true(atmosphere_dim, lat_grid, lat_true, lon_true);
2241 chk_rte_pos(atmosphere_dim, rtp_pos);
2242 surface_type_mask.checksize_strict();
2243 //
2244 chk_griddedfield_gridname(surface_type_mask, gfield_latID, "Latitude");
2245 chk_griddedfield_gridname(surface_type_mask, gfield_lonID, "Longitude");
2246 //
2247 const Index nlat = surface_type_mask.data.nrows();
2248 const Index nlon = surface_type_mask.data.ncols();
2249 //
2250 if (nlat < 2 || nlon < 2) {
2251 ostringstream os;
2252 os << "The data in *surface_type_mask* must span a geographical "
2253 << "region. Accordingly,\nthe latitude and longitude grids must "
2254 << "both have a length >= 2.";
2255 throw runtime_error(os.str());
2256 }
2257
2258 const Vector& GFlat = surface_type_mask.get_numeric_grid(gfield_latID);
2259 const Vector& GFlon = surface_type_mask.get_numeric_grid(gfield_lonID);
2260
2261 // Determine true geographical position
2262 Vector lat(1), lon(1);
2264 lat[0], lon[0], atmosphere_dim, lat_grid, lat_true, lon_true, rtp_pos);
2265
2266 // Ensure correct coverage of lon grid
2267 Vector lon_shifted;
2268 lon_shiftgrid(lon_shifted, GFlon, lon[0]);
2269
2270 // Check if lat/lon we need are actually covered
2271 chk_if_in_range("rtp_pos.lat", lat[0], GFlat[0], GFlat[nlat - 1]);
2272 chk_if_in_range("rtp_pos.lon", lon[0], lon_shifted[0], lon_shifted[nlon - 1]);
2273
2274 // Grid positions
2275 GridPos gp_lat, gp_lon;
2276 gridpos(gp_lat, GFlat, lat[0]);
2277 gridpos(gp_lon, lon_shifted, lon[0]);
2278
2279 if (method == "nearest" ) {
2280 // Extract closest point
2281 Index ilat, ilon;
2282 if (gp_lat.fd[0] < 0.5) {
2283 ilat = gp_lat.idx;
2284 } else {
2285 ilat = gp_lat.idx + 1;
2286 }
2287 if (gp_lon.fd[0] < 0.5) {
2288 ilon = gp_lon.idx;
2289 } else {
2290 ilon = gp_lon.idx + 1;
2291 }
2292 //
2293 surface_types.resize(1);
2294 surface_types_aux.resize(1);
2295 surface_types_weights.resize(1);
2296 surface_types[0] = (Index)floor(surface_type_mask.data(ilat, ilon));
2297 surface_types_aux[0] = surface_type_mask.data(ilat, ilon) -
2298 Numeric(surface_types[0]);
2299 surface_types_weights[0] = 1.0;
2300 }
2301
2302 else if (method == "linear" ) {
2303 // Determine types involved
2304 ArrayOfIndex types0(4), types(4);
2305 Index ntypes = 1;
2306 //
2307 types0[0] = (Index)floor(surface_type_mask.data(gp_lat.idx,gp_lon.idx));
2308 types0[1] = (Index)floor(surface_type_mask.data(gp_lat.idx,gp_lon.idx+1));
2309 types0[2] = (Index)floor(surface_type_mask.data(gp_lat.idx+1,gp_lon.idx));
2310 types0[3] = (Index)floor(surface_type_mask.data(gp_lat.idx+1,gp_lon.idx+1));
2311 //
2312 types[0] = types0[0];
2313 for (Index t=1; t<4; t++) {
2314 bool unique = true;
2315 for (Index n=0; n<ntypes && unique; n++)
2316 if (types0[t] == types[n] ) { unique = false; }
2317 if (unique) {
2318 types[ntypes] = types0[t];
2319 ntypes += 1;
2320 }
2321 }
2322 //
2323 surface_types.resize(ntypes);
2324 surface_types_aux.resize(ntypes);
2325 surface_types_weights.resize(ntypes);
2326
2327 // Interpolation weights
2328 Vector itw(4);
2329 interpweights(itw, gp_lat, gp_lon);
2330
2331 // Determine weight for each type, and make an interpolation of aux value
2332 for (Index t=0; t<ntypes; t++ ) {
2333 Numeric wtot = 0, auxtot = 0;
2334 Numeric ntype = (Numeric)types[t];
2335 if (types[t] == types0[0]) { wtot += itw[0];
2336 auxtot += itw[0]*(surface_type_mask.data(gp_lat.idx,gp_lon.idx)-ntype);
2337 };
2338 if (types[t] == types0[1]) { wtot += itw[1];
2339 auxtot += itw[1]*(surface_type_mask.data(gp_lat.idx,gp_lon.idx+1)-ntype);
2340 };
2341 if (types[t] == types0[2]) { wtot += itw[2];
2342 auxtot += itw[2]*(surface_type_mask.data(gp_lat.idx+1,gp_lon.idx)-ntype);
2343 };
2344 if (types[t] == types0[3]) { wtot += itw[3];
2345 auxtot += itw[3]*(surface_type_mask.data(gp_lat.idx+1,gp_lon.idx+1)-ntype);
2346 };
2347 //
2348 surface_types[t] = types[t];
2349 surface_types_aux[t] = auxtot / wtot;
2350 surface_types_weights[t] = wtot;
2351 }
2352 }
2353
2354 else {
2355 ostringstream os;
2356 os << "The allowed options for *method are: \"nearest\" and \"linear\",\n"
2357 << "but you have selected: " << method;
2358 throw runtime_error(os.str());
2359 }
2360}
2361
2362/* Workspace method: Doxygen documentation will be auto-generated */
2364 Numeric& surface_skin_t,
2365 Matrix& surface_los,
2366 Tensor4& surface_rmatrix,
2367 Matrix& surface_emission,
2368 const Vector& f_grid,
2369 const Vector& rtp_pos,
2370 const Vector& rtp_los,
2371 const ArrayOfAgenda& surface_rtprop_agenda_array,
2372 const ArrayOfIndex& surface_types,
2373 const Vector& surface_types_aux,
2374 const Vector& surface_types_weights,
2375 const Verbosity&) {
2376 if (surface_types.nelem() != 1) {
2377 ostringstream os;
2378 os << "This method requires that *surface_types* have length 1.\n"
2379 << "If you are trying to apply a mixture model, you have to use\n"
2380 << "a setup based on *iySurfaceCallAgendaX* instead.";
2381 throw runtime_error(os.str());
2382 }
2383 if (surface_types[0] < 0)
2384 throw runtime_error("*surface_type* is not allowed to be negative.");
2385 if (surface_types[0] >= surface_rtprop_agenda_array.nelem()) {
2386 ostringstream os;
2387 os << "*surface_rtprop_agenda_array* has only "
2388 << surface_rtprop_agenda_array.nelem()
2389 << " elements,\n while you have selected element " << surface_types[0];
2390 throw runtime_error(os.str());
2391 }
2392 if (abs(surface_types_weights[0]-1)>1e-4)
2393 throw runtime_error("Sum of *surface_types_weights* deviates from 1.");
2394
2396 surface_skin_t,
2397 surface_emission,
2398 surface_los,
2399 surface_rmatrix,
2400 surface_types[0],
2401 f_grid,
2402 rtp_pos,
2403 rtp_los,
2404 surface_types_aux[0],
2405 surface_rtprop_agenda_array);
2406}
2407
2408/* Workspace method: Doxygen documentation will be auto-generated */
2409void SurfaceBlackbody(Matrix& surface_los,
2410 Tensor4& surface_rmatrix,
2411 ArrayOfTensor4& dsurface_rmatrix_dx,
2412 Matrix& surface_emission,
2413 ArrayOfMatrix& dsurface_emission_dx,
2414 const Index& stokes_dim,
2415 const Index& atmosphere_dim,
2416 const Vector& lat_grid,
2417 const Vector& lon_grid,
2418 const Vector& f_grid,
2419 const Vector& rtp_pos,
2420 const Vector& rtp_los,
2421 const Tensor3& surface_props_data,
2422 const ArrayOfString& surface_props_names,
2423 const ArrayOfString& dsurface_names,
2424 const Index& jacobian_do,
2425 const Verbosity& verbosity) {
2426 // Check surface_data
2427 surface_props_check(atmosphere_dim,
2428 lat_grid,
2429 lon_grid,
2430 surface_props_data,
2431 surface_props_names);
2432
2433 // Interplation grid positions and weights
2434 ArrayOfGridPos gp_lat(1), gp_lon(1);
2435 Matrix itw;
2437 gp_lat[0], gp_lon[0], atmosphere_dim, lat_grid, lon_grid, rtp_pos);
2438 interp_atmsurface_gp2itw(itw, atmosphere_dim, gp_lat, gp_lon);
2439
2440 // Skin temperature
2441 Vector skin_t(1);
2442 surface_props_interp(skin_t,
2443 "Skin temperature",
2444 atmosphere_dim,
2445 gp_lat,
2446 gp_lon,
2447 itw,
2448 surface_props_data,
2449 surface_props_names);
2450
2451 surfaceBlackbody(surface_los,
2452 surface_rmatrix,
2453 surface_emission,
2454 atmosphere_dim,
2455 f_grid,
2456 stokes_dim,
2457 rtp_pos,
2458 rtp_los,
2459 skin_t[0],
2460 verbosity);
2461
2462 surface_rmatrix.resize(1, f_grid.nelem(), stokes_dim, stokes_dim);
2463 surface_rmatrix = 0.0;
2464
2465 // Jacobian part
2466 if (jacobian_do) {
2467 dsurface_check(surface_props_names,
2468 dsurface_names,
2469 dsurface_rmatrix_dx,
2470 dsurface_emission_dx);
2471
2472 Index irq;
2473
2474 // Skin temperature
2475 irq = find_first(dsurface_names, String("Skin temperature"));
2476 if (irq >= 0) {
2477 Matrix dbdt(f_grid.nelem(), 1);
2478 dplanck_dt(dbdt(joker, 0), f_grid, skin_t[0]);
2479 dsurface_emission_dx[irq] = dbdt;
2480 dsurface_rmatrix_dx[irq].resize(surface_rmatrix.nbooks(),
2481 surface_rmatrix.npages(),
2482 surface_rmatrix.nrows(),
2483 surface_rmatrix.ncols());
2484 dsurface_rmatrix_dx[irq] = 0;
2485 }
2486 }
2487}
2488
2489/* Workspace method: Doxygen documentation will be auto-generated */
2490void SurfaceDummy(ArrayOfTensor4& dsurface_rmatrix_dx,
2491 ArrayOfMatrix& dsurface_emission_dx,
2492 const Index& atmosphere_dim,
2493 const Vector& lat_grid,
2494 const Vector& lon_grid,
2495 const Tensor3& surface_props_data,
2496 const ArrayOfString& surface_props_names,
2497 const ArrayOfString& dsurface_names,
2498 const Index& jacobian_do,
2499 const Verbosity&) {
2500 if (surface_props_names.nelem())
2501 throw runtime_error(
2502 "When calling this method, *surface_props_names* should be empty.");
2503
2504 surface_props_check(atmosphere_dim,
2505 lat_grid,
2506 lon_grid,
2507 surface_props_data,
2508 surface_props_names);
2509
2510 if (jacobian_do) {
2511 dsurface_check(surface_props_names,
2512 dsurface_names,
2513 dsurface_rmatrix_dx,
2514 dsurface_emission_dx);
2515 }
2516}
2517
2518/* Workspace method: Doxygen documentation will be auto-generated */
2519void SurfaceFastem(Matrix& surface_los,
2520 Tensor4& surface_rmatrix,
2521 ArrayOfTensor4& dsurface_rmatrix_dx,
2522 Matrix& surface_emission,
2523 ArrayOfMatrix& dsurface_emission_dx,
2524 const Index& stokes_dim,
2525 const Index& atmosphere_dim,
2526 const Vector& lat_grid,
2527 const Vector& lon_grid,
2528 const Vector& f_grid,
2529 const Vector& rtp_pos,
2530 const Vector& rtp_los,
2531 const Tensor3& surface_props_data,
2532 const ArrayOfString& surface_props_names,
2533 const ArrayOfString& dsurface_names,
2534 const Index& jacobian_do,
2535 const Vector& transmittance,
2536 const Index& fastem_version,
2537 const Verbosity& verbosity) {
2538 // Check surface_data
2539 surface_props_check(atmosphere_dim,
2540 lat_grid,
2541 lon_grid,
2542 surface_props_data,
2543 surface_props_names);
2544
2545 // Interplation grid positions and weights
2546 ArrayOfGridPos gp_lat(1), gp_lon(1);
2547 Matrix itw;
2549 gp_lat[0], gp_lon[0], atmosphere_dim, lat_grid, lon_grid, rtp_pos);
2550 interp_atmsurface_gp2itw(itw, atmosphere_dim, gp_lat, gp_lon);
2551
2552 // Skin temperature
2553 Vector skin_t(1);
2554 surface_props_interp(skin_t,
2555 "Water skin temperature",
2556 atmosphere_dim,
2557 gp_lat,
2558 gp_lon,
2559 itw,
2560 surface_props_data,
2561 surface_props_names);
2562
2563 // Wind speed
2564 Vector wind_speed(1);
2565 surface_props_interp(wind_speed,
2566 "Wind speed",
2567 atmosphere_dim,
2568 gp_lat,
2569 gp_lon,
2570 itw,
2571 surface_props_data,
2572 surface_props_names);
2573
2574 // Wind direction
2575 Vector wind_direction(1);
2576 surface_props_interp(wind_direction,
2577 "Wind direction",
2578 atmosphere_dim,
2579 gp_lat,
2580 gp_lon,
2581 itw,
2582 surface_props_data,
2583 surface_props_names);
2584
2585 // Salinity
2586 Vector salinity(1);
2587 surface_props_interp(salinity,
2588 "Salinity",
2589 atmosphere_dim,
2590 gp_lat,
2591 gp_lon,
2592 itw,
2593 surface_props_data,
2594 surface_props_names);
2595
2596 // Call FASTEM
2597 surfaceFastem(surface_los,
2598 surface_rmatrix,
2599 surface_emission,
2600 atmosphere_dim,
2601 stokes_dim,
2602 f_grid,
2603 rtp_pos,
2604 rtp_los,
2605 skin_t[0],
2606 salinity[0],
2607 wind_speed[0],
2608 wind_direction[0],
2609 transmittance,
2610 fastem_version,
2611 verbosity);
2612
2613 // Jacobian part
2614 if (jacobian_do) {
2615 dsurface_check(surface_props_names,
2616 dsurface_names,
2617 dsurface_rmatrix_dx,
2618 dsurface_emission_dx);
2619
2620 Index irq;
2621
2622 // Skin temperature
2623 irq = find_first(dsurface_names, String("Water skin temperature"));
2624 if (irq >= 0) {
2625 const Numeric dd = 0.1;
2626 Matrix surface_los2;
2627 surfaceFastem(surface_los2,
2628 dsurface_rmatrix_dx[irq],
2629 dsurface_emission_dx[irq],
2630 atmosphere_dim,
2631 stokes_dim,
2632 f_grid,
2633 rtp_pos,
2634 rtp_los,
2635 skin_t[0] + dd,
2636 salinity[0],
2637 wind_speed[0],
2638 wind_direction[0],
2639 transmittance,
2640 fastem_version,
2641 verbosity);
2642 //
2643 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
2644 dsurface_rmatrix_dx[irq] /= dd;
2645 dsurface_emission_dx[irq] -= surface_emission;
2646 dsurface_emission_dx[irq] /= dd;
2647 }
2648
2649 // Wind speed
2650 irq = find_first(dsurface_names, String("Wind speed"));
2651 if (irq >= 0) {
2652 const Numeric dd = 0.1;
2653 Matrix surface_los2;
2654 surfaceFastem(surface_los2,
2655 dsurface_rmatrix_dx[irq],
2656 dsurface_emission_dx[irq],
2657 atmosphere_dim,
2658 stokes_dim,
2659 f_grid,
2660 rtp_pos,
2661 rtp_los,
2662 skin_t[0],
2663 salinity[0],
2664 wind_speed[0] + dd,
2665 wind_direction[0],
2666 transmittance,
2667 fastem_version,
2668 verbosity);
2669 //
2670 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
2671 dsurface_rmatrix_dx[irq] /= dd;
2672 dsurface_emission_dx[irq] -= surface_emission;
2673 dsurface_emission_dx[irq] /= dd;
2674 }
2675
2676 // Wind direction
2677 irq = find_first(dsurface_names, String("Wind direction"));
2678 if (irq >= 0) {
2679 const Numeric dd = 1;
2680 Matrix surface_los2;
2681 surfaceFastem(surface_los2,
2682 dsurface_rmatrix_dx[irq],
2683 dsurface_emission_dx[irq],
2684 atmosphere_dim,
2685 stokes_dim,
2686 f_grid,
2687 rtp_pos,
2688 rtp_los,
2689 skin_t[0],
2690 salinity[0],
2691 wind_speed[0],
2692 wind_direction[0] + dd,
2693 transmittance,
2694 fastem_version,
2695 verbosity);
2696 //
2697 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
2698 dsurface_rmatrix_dx[irq] /= dd;
2699 dsurface_emission_dx[irq] -= surface_emission;
2700 dsurface_emission_dx[irq] /= dd;
2701 }
2702
2703 // Salinity
2704 irq = find_first(dsurface_names, String("Salinity"));
2705 if (irq >= 0) {
2706 const Numeric dd = 0.0005;
2707 Matrix surface_los2;
2708 surfaceFastem(surface_los2,
2709 dsurface_rmatrix_dx[irq],
2710 dsurface_emission_dx[irq],
2711 atmosphere_dim,
2712 stokes_dim,
2713 f_grid,
2714 rtp_pos,
2715 rtp_los,
2716 skin_t[0],
2717 salinity[0] + dd,
2718 wind_speed[0],
2719 wind_direction[0],
2720 transmittance,
2721 fastem_version,
2722 verbosity);
2723 //
2724 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
2725 dsurface_rmatrix_dx[irq] /= dd;
2726 dsurface_emission_dx[irq] -= surface_emission;
2727 dsurface_emission_dx[irq] /= dd;
2728 }
2729 }
2730}
2731
2732
2733/* Workspace method: Doxygen documentation will be auto-generated */
2735 Tensor4& surface_rmatrix,
2736 ArrayOfTensor4& dsurface_rmatrix_dx,
2737 Matrix& surface_emission,
2738 ArrayOfMatrix& dsurface_emission_dx,
2739 const Index& stokes_dim,
2740 const Index& atmosphere_dim,
2741 const Vector& lat_grid,
2742 const Vector& lon_grid,
2743 const Vector& f_grid,
2744 const Vector& rtp_pos,
2745 const Vector& rtp_los,
2746 const Vector& specular_los,
2747 const Tensor3& surface_props_data,
2748 const ArrayOfString& surface_props_names,
2749 const ArrayOfString& dsurface_names,
2750 const Index& jacobian_do,
2751 const Vector& f_reflectivities,
2752 const Verbosity& verbosity) {
2753 // Check surface_data
2754 surface_props_check(atmosphere_dim,
2755 lat_grid,
2756 lon_grid,
2757 surface_props_data,
2758 surface_props_names);
2759
2760 // Check of GINs
2761 chk_if_increasing("GIN *f_reflectivities*", f_reflectivities);
2762
2763 // Sizes
2764 const Index nr = f_reflectivities.nelem();
2765 const Index nf = f_grid.nelem();
2766
2767 // Interplation grid positions and weights
2768 ArrayOfGridPos gp_lat(1), gp_lon(1);
2769 Matrix itw;
2771 gp_lat[0], gp_lon[0], atmosphere_dim, lat_grid, lon_grid, rtp_pos);
2772 interp_atmsurface_gp2itw(itw, atmosphere_dim, gp_lat, gp_lon);
2773
2774 // Skin temperature
2775 Vector surface_skin_t(1);
2776 surface_props_interp(surface_skin_t,
2777 "Skin temperature",
2778 atmosphere_dim,
2779 gp_lat,
2780 gp_lon,
2781 itw,
2782 surface_props_data,
2783 surface_props_names);
2784
2785 // Reflectivities
2786 Vector reflectivities(nr);
2787 //
2788 for (Index i=0; i<nr; i++) {
2789 Vector rv(1);
2790 ostringstream sstr;
2791 sstr << "Scalar reflectivity " << i;
2793 sstr.str(),
2794 atmosphere_dim,
2795 gp_lat,
2796 gp_lon,
2797 itw,
2798 surface_props_data,
2799 surface_props_names);
2800 reflectivities[i] = rv[0];
2801 }
2802
2803 // Create surface_scalar_reflectivity
2804 Vector surface_scalar_reflectivity(nf);
2805 ArrayOfGridPos gp(nf);
2806 if (nr==1) {
2807 gp4length1grid(gp);
2808 } else {
2809 gridpos(gp, f_reflectivities, f_grid, 1e99);
2811 }
2812 Matrix itw2(nf, 2);
2813 interpweights(itw2, gp);
2814 interp( surface_scalar_reflectivity, itw2, reflectivities, gp);
2815
2816 // Call non-Jacobian method
2818 surface_rmatrix,
2819 surface_emission,
2820 f_grid,
2821 stokes_dim,
2822 atmosphere_dim,
2823 rtp_pos,
2824 rtp_los,
2825 specular_los,
2826 surface_skin_t[0],
2827 surface_scalar_reflectivity,
2828 verbosity);
2829
2830 // Jacobian part
2831 if (jacobian_do) {
2832 dsurface_check(surface_props_names,
2833 dsurface_names,
2834 dsurface_rmatrix_dx,
2835 dsurface_emission_dx);
2836
2837 Index irq;
2838
2839 // Skin temperature
2840 irq = find_first(dsurface_names, String("Skin temperature"));
2841 if (irq >= 0) {
2842 const Numeric dd = 0.1;
2843 Matrix surface_los2;
2844 surfaceFlatScalarReflectivity(surface_los2,
2845 dsurface_rmatrix_dx[irq],
2846 dsurface_emission_dx[irq],
2847 f_grid,
2848 stokes_dim,
2849 atmosphere_dim,
2850 rtp_pos,
2851 rtp_los,
2852 specular_los,
2853 surface_skin_t[0]+dd,
2854 surface_scalar_reflectivity,
2855 verbosity);
2856 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
2857 dsurface_rmatrix_dx[irq] /= dd;
2858 dsurface_emission_dx[irq] -= surface_emission;
2859 dsurface_emission_dx[irq] /= dd;
2860 }
2861
2862 // Reflectivities
2863 for (Index i=0; i<nr; i++) {
2864 ostringstream sstr;
2865 sstr << "Scalar reflectivity " << i;
2866 irq = find_first(dsurface_names, String(sstr.str()));
2867 if (irq >= 0) {
2868 const Numeric dd = 1e-4;
2869 Vector rpert = reflectivities;
2870 rpert[i] += dd;
2871 interp( surface_scalar_reflectivity, itw2, rpert, gp);
2872 //
2873 Matrix surface_los2;
2874 surfaceFlatScalarReflectivity(surface_los2,
2875 dsurface_rmatrix_dx[irq],
2876 dsurface_emission_dx[irq],
2877 f_grid,
2878 stokes_dim,
2879 atmosphere_dim,
2880 rtp_pos,
2881 rtp_los,
2882 specular_los,
2883 surface_skin_t[0],
2884 surface_scalar_reflectivity,
2885 verbosity);
2886 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
2887 dsurface_rmatrix_dx[irq] /= dd;
2888 dsurface_emission_dx[irq] -= surface_emission;
2889 dsurface_emission_dx[irq] /= dd;
2890 }
2891 }
2892 }
2893}
2894
2895/* Workspace method: Doxygen documentation will be auto-generated */
2896void SurfaceTessem(Matrix& surface_los,
2897 Tensor4& surface_rmatrix,
2898 ArrayOfTensor4& dsurface_rmatrix_dx,
2899 Matrix& surface_emission,
2900 ArrayOfMatrix& dsurface_emission_dx,
2901 const Index& stokes_dim,
2902 const Index& atmosphere_dim,
2903 const Vector& lat_grid,
2904 const Vector& lon_grid,
2905 const Vector& f_grid,
2906 const Vector& rtp_pos,
2907 const Vector& rtp_los,
2908 const TessemNN& net_h,
2909 const TessemNN& net_v,
2910 const Tensor3& surface_props_data,
2911 const ArrayOfString& surface_props_names,
2912 const ArrayOfString& dsurface_names,
2913 const Index& jacobian_do,
2914 const Verbosity& verbosity) {
2915 // Check surface_data
2916 surface_props_check(atmosphere_dim,
2917 lat_grid,
2918 lon_grid,
2919 surface_props_data,
2920 surface_props_names);
2921
2922 // Interplation grid positions and weights
2923 ArrayOfGridPos gp_lat(1), gp_lon(1);
2924 Matrix itw;
2926 gp_lat[0], gp_lon[0], atmosphere_dim, lat_grid, lon_grid, rtp_pos);
2927 interp_atmsurface_gp2itw(itw, atmosphere_dim, gp_lat, gp_lon);
2928
2929 // Skin temperature
2930 Vector skin_t(1);
2931 surface_props_interp(skin_t,
2932 "Water skin temperature",
2933 atmosphere_dim,
2934 gp_lat,
2935 gp_lon,
2936 itw,
2937 surface_props_data,
2938 surface_props_names);
2939
2940 // Wind speed
2941 Vector wind_speed(1);
2942 surface_props_interp(wind_speed,
2943 "Wind speed",
2944 atmosphere_dim,
2945 gp_lat,
2946 gp_lon,
2947 itw,
2948 surface_props_data,
2949 surface_props_names);
2950
2951 // Salinity
2952 Vector salinity(1);
2953 surface_props_interp(salinity,
2954 "Salinity",
2955 atmosphere_dim,
2956 gp_lat,
2957 gp_lon,
2958 itw,
2959 surface_props_data,
2960 surface_props_names);
2961
2962 // Call TESSEM
2963 surfaceTessem(surface_los,
2964 surface_rmatrix,
2965 surface_emission,
2966 atmosphere_dim,
2967 stokes_dim,
2968 f_grid,
2969 rtp_pos,
2970 rtp_los,
2971 skin_t[0],
2972 net_h,
2973 net_v,
2974 salinity[0],
2975 wind_speed[0],
2976 verbosity);
2977
2978 // Jacobian part
2979 if (jacobian_do) {
2980 dsurface_check(surface_props_names,
2981 dsurface_names,
2982 dsurface_rmatrix_dx,
2983 dsurface_emission_dx);
2984
2985 Index irq;
2986
2987 // Skin temperature
2988 irq = find_first(dsurface_names, String("Water skin temperature"));
2989 if (irq >= 0) {
2990 const Numeric dd = 0.1;
2991 Matrix surface_los2;
2992 surfaceTessem(surface_los2,
2993 dsurface_rmatrix_dx[irq],
2994 dsurface_emission_dx[irq],
2995 atmosphere_dim,
2996 stokes_dim,
2997 f_grid,
2998 rtp_pos,
2999 rtp_los,
3000 skin_t[0] + dd,
3001 net_h,
3002 net_v,
3003 salinity[0],
3004 wind_speed[0],
3005 verbosity);
3006 //
3007 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
3008 dsurface_rmatrix_dx[irq] /= dd;
3009 dsurface_emission_dx[irq] -= surface_emission;
3010 dsurface_emission_dx[irq] /= dd;
3011 }
3012
3013 // Wind speed
3014 irq = find_first(dsurface_names, String("Wind speed"));
3015 if (irq >= 0) {
3016 const Numeric dd = 0.1;
3017 Matrix surface_los2;
3018 surfaceTessem(surface_los2,
3019 dsurface_rmatrix_dx[irq],
3020 dsurface_emission_dx[irq],
3021 atmosphere_dim,
3022 stokes_dim,
3023 f_grid,
3024 rtp_pos,
3025 rtp_los,
3026 skin_t[0],
3027 net_h,
3028 net_v,
3029 salinity[0],
3030 wind_speed[0] + dd,
3031 verbosity);
3032 //
3033 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
3034 dsurface_rmatrix_dx[irq] /= dd;
3035 dsurface_emission_dx[irq] -= surface_emission;
3036 dsurface_emission_dx[irq] /= dd;
3037 }
3038
3039 // Salinity
3040 irq = find_first(dsurface_names, String("Salinity"));
3041 if (irq >= 0) {
3042 const Numeric dd = 0.0005;
3043 Matrix surface_los2;
3044 surfaceTessem(surface_los2,
3045 dsurface_rmatrix_dx[irq],
3046 dsurface_emission_dx[irq],
3047 atmosphere_dim,
3048 stokes_dim,
3049 f_grid,
3050 rtp_pos,
3051 rtp_los,
3052 skin_t[0],
3053 net_h,
3054 net_v,
3055 salinity[0] + dd,
3056 wind_speed[0],
3057 verbosity);
3058 //
3059 dsurface_rmatrix_dx[irq] -= surface_rmatrix;
3060 dsurface_rmatrix_dx[irq] /= dd;
3061 dsurface_emission_dx[irq] -= surface_emission;
3062 dsurface_emission_dx[irq] /= dd;
3063 }
3064 }
3065}
3066
3067/* Workspace method: Doxygen documentation will be auto-generated */
3069 const ArrayOfString& iy_aux_vars,
3070 const ArrayOfMatrix& iy_aux,
3071 const Verbosity&) {
3072 Index ihit = -1;
3073
3074 for (Index i = 0; i < iy_aux_vars.nelem(); i++) {
3075 if (iy_aux_vars[i] == "Optical depth") {
3076 ihit = i;
3077 break;
3078 }
3079 }
3080
3081 if (ihit < 0)
3082 throw runtime_error("No element in *iy_aux* holds optical depths.");
3083
3084 const Index n = iy_aux[ihit].nrows();
3085
3086 transmittance.resize(n);
3087
3088 for (Index i = 0; i < n; i++) {
3089 transmittance[i] = exp(-iy_aux[ihit](i, 0));
3090 }
3091}
This file contains the definition of Array.
Index find_first(const Array< base > &x, const base &w)
Find first occurance.
Definition: array.h:292
The global header file for ARTS.
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 Numeric surface_type_aux, const ArrayOfAgenda &input_agenda_array)
Definition: auto_md.cc:24443
void surface_rtprop_sub_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:24510
void iy_main_agendaExecute(Workspace &ws, Matrix &iy, ArrayOfMatrix &iy_aux, Ppath &ppath, ArrayOfTensor3 &diy_dx, 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:23520
void iy_surface_agenda_arrayExecute(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, const Index agenda_array_index, const String &iy_unit, const Tensor3 &iy_transmittance, const Index iy_id, const Index cloudbox_on, const Index jacobian_do, const Agenda &iy_main_agenda, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Vector &rte_pos2, const Numeric surface_type_aux, const ArrayOfAgenda &input_agenda_array)
Definition: auto_md.cc:23781
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:24392
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_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:44
String name() const
Agenda name.
This can be used to make arrays out of anything.
Definition: array.h:107
Index nelem() const ARTS_NOEXCEPT
Number of elements.
Definition: array.h:195
Index nrows() const ARTS_NOEXCEPT
Returns the number of rows.
Definition: matpackI.cc:433
Index ncols() const ARTS_NOEXCEPT
Returns the number of columns.
Definition: matpackI.cc:436
Index npages() const
Returns the number of pages.
Definition: matpackIII.h:144
Index nrows() const
Returns the number of rows.
Definition: matpackIII.h:147
Index ncols() const
Returns the number of columns.
Definition: matpackIII.h:150
Index npages() const
Returns the number of pages.
Definition: matpackIV.cc:58
Index nbooks() const
Returns the number of books.
Definition: matpackIV.cc:55
Index ncols() const
Returns the number of columns.
Definition: matpackIV.cc:64
Index nrows() const
Returns the number of rows.
Definition: matpackIV.cc:61
Index ncols() const
Returns the number of columns.
Definition: matpackV.cc:56
Index nshelves() const
Returns the number of shelves.
Definition: matpackV.cc:44
Index npages() const
Returns the number of pages.
Definition: matpackV.cc:50
Index nrows() const
Returns the number of rows.
Definition: matpackV.cc:53
Index nbooks() const
Returns the number of books.
Definition: matpackV.cc:47
Index ncols() const
Returns the number of columns.
Definition: matpackVI.cc:57
Index npages() const
Returns the number of pages.
Definition: matpackVI.cc:51
Index nrows() const
Returns the number of rows.
Definition: matpackVI.cc:54
Index nshelves() const
Returns the number of shelves.
Definition: matpackVI.cc:45
Index nbooks() const
Returns the number of books.
Definition: matpackVI.cc:48
Index nvitrines() const
Returns the number of vitrines.
Definition: matpackVI.cc:42
Index nelem() const ARTS_NOEXCEPT
Returns the number of elements.
Definition: matpackI.cc:48
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:1225
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1056
The range class.
Definition: matpackI.h:165
A telsem atlas.
Definition: telsem.h:59
Index calc_cellnum(Numeric lat, Numeric lon) const
Definition: telsem.cc:142
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:170
Vector get_emis_h(Index cellnum) const
Definition: telsem.h:159
std::pair< Numeric, Numeric > get_coordinates(Index cellnum) const
Definition: telsem.cc:225
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:287
Vector get_emis_v(Index i) const
Definition: telsem.h:137
The Tensor3 class.
Definition: matpackIII.h:339
void resize(Index p, Index r, Index c)
Resize function.
Definition: matpackIII.cc:658
The Tensor4 class.
Definition: matpackIV.h:421
void resize(Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackIV.cc:1058
The VectorView class.
Definition: matpackI.h:626
The Vector class.
Definition: matpackI.h:876
void resize(Index n)
Resize function.
Definition: matpackI.cc:408
Workspace class.
Definition: workspace_ng.h:40
void mult(MatrixView C, ConstMatrixView A, const Block &B)
#define ARTS_ASSERT(condition,...)
Definition: debug.h:83
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:107
This file contains functions that are adapted from FASTEM code which is used to calculate surface emi...
Numeric sphdist(const Numeric &lat1, const Numeric &lon1, const Numeric &lat2, const Numeric &lon2)
sphdist
Definition: geodetic.cc:1335
void lon_shiftgrid(Vector &longrid_out, ConstVectorView longrid_in, const Numeric lon)
lon_shiftgrid
Definition: geodetic.cc:1408
Numeric refell2d(ConstVectorView refellipsoid, ConstVectorView lat_grid, const GridPos gp)
refell2d
Definition: geodetic.cc:1304
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:474
void jacobian_type_extrapol(ArrayOfGridPos &gp)
Adopts grid positions to extrapolation used for jacobians.
Definition: jacobian.cc:835
#define abs(x)
#define q
#define min(a, b)
#define max(a, b)
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:1024
void surfaceSplitSpecularTo3beams(Matrix &surface_los, Tensor4 &surface_rmatrix, const Index &atmosphere_dim, const Vector &rtp_pos, const Vector &rtp_los, const Numeric &specular_factor, const Numeric &dza, const Verbosity &)
WORKSPACE METHOD: surfaceSplitSpecularTo3beams.
Definition: m_surface.cc:1759
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:2098
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:1868
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:200
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:1131
void surface_rtpropCallAgendaX(Workspace &ws, 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 ArrayOfIndex &surface_types, const Vector &surface_types_aux, const Vector &surface_types_weights, const Verbosity &)
WORKSPACE METHOD: surface_rtpropCallAgendaX.
Definition: m_surface.cc:2363
void surfaceSemiSpecularBy3beams(Workspace &ws, Numeric &surface_skin_t, Matrix &surface_los, Tensor4 &surface_rmatrix, Matrix &surface_emission, const Index &atmosphere_dim, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Agenda &surface_rtprop_sub_agenda, const Numeric &specular_factor, const Numeric &dza, const Verbosity &)
WORKSPACE METHOD: surfaceSemiSpecularBy3beams.
Definition: m_surface.cc:1612
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:767
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:915
Numeric EARTH_RADIUS
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:1210
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:2409
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:730
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:876
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 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:571
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:2519
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:1444
void surface_typeInterpTypeMask(ArrayOfIndex &surface_types, Vector &surface_types_aux, Vector &surface_types_weights, 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 String &method, const Verbosity &)
WORKSPACE METHOD: surface_typeInterpTypeMask.
Definition: m_surface.cc:2223
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:1964
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:1280
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:1364
void iySurfaceCallAgendaX(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, const String &iy_unit, const Tensor3 &iy_transmittance, const Index &iy_id, const Index &cloudbox_on, const Index &jacobian_do, const Vector &f_grid, const Agenda &iy_main_agenda, const Vector &rtp_pos, const Vector &rtp_los, const Vector &rte_pos2, const ArrayOfAgenda &iy_surface_agenda_array, const ArrayOfIndex &surface_types, const Vector &surface_types_aux, const Vector &surface_types_weights, const Verbosity &)
WORKSPACE METHOD: iySurfaceCallAgendaX.
Definition: m_surface.cc:247
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 &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:443
void transmittanceFromIy_aux(Vector &transmittance, const ArrayOfString &iy_aux_vars, const ArrayOfMatrix &iy_aux, const Verbosity &)
WORKSPACE METHOD: transmittanceFromIy_aux.
Definition: m_surface.cc:3068
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:1510
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:326
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:2205
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:2896
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:2490
const Numeric DEG2RAD
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:141
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:66
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:2734
Numeric fac(const Index n)
fac
Definition: math_funcs.cc:63
Numeric sign(const Numeric &x)
sign
Definition: math_funcs.cc:417
void cross3(VectorView c, const ConstVectorView &a, const ConstVectorView &b) ARTS_NOEXCEPT
cross3
Definition: matpackI.cc:1393
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
const Joker joker
std::complex< Numeric > Complex
Declarations having to do with the four output streams.
#define CREATE_OUT3
Definition: messages.h:207
#define CREATE_OUT2
Definition: messages.h:206
Array< String > ArrayOfString
An array of Strings.
Definition: mystring.h:290
my_basic_string< char > String
The String type for ARTS.
Definition: mystring.h:287
constexpr Numeric l(const Index p0, const Index n, const Numeric x, const SortedVectorType &xi, const Index j, const std::pair< Numeric, Numeric > cycle={ -180, 180}) noexcept
invlib::Vector< ArtsVector > Vector
invlib wrapper type for ARTS vectors.
Definition: oem.h:31
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:312
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:299
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:1067
void resolve_lon(Numeric &lon, const Numeric &lon5, const Numeric &lon6)
Resolves which longitude angle that shall be used.
Definition: ppath.cc:500
Numeric plevel_angletilt(const Numeric &r, const Numeric &c1)
Calculates the angular tilt of the surface or a pressure level.
Definition: ppath.cc:616
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:579
Propagation path structure and functions.
#define d
#define w
#define c
#define b
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:2034
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:1973
void adjust_los(VectorView los, const Index &atmosphere_dim)
Ensures that the zenith and azimuth angles of a line-of-sight vector are inside defined ranges.
Definition: rte.cc:83
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
Numeric fd[2]
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.h:48
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...