ARTS 2.5.0 (git: 9ee3ac6c)
m_cloudbox.cc
Go to the documentation of this file.
1/* Copyright (C) 2002-2017
2 Patrick Eriksson <patrick.eriksson@chalmers.se>
3 Stefan Buehler <sbuehler@ltu.se>
4 Claudia Emde <claudia.emde@dlr.de>
5 Cory Davis <cory.davis@metservice.com>
6 Jana Mendrok <jana.mendrok@gmail.com>
7 Daniel Kreyling <daniel.kreyling@nict.go.jp>
8 Manfred Brath <manfred.brath@uni-hamburg.de>
9
10 This program is free software; you can redistribute it and/or modify it
11 under the terms of the GNU General Public License as published by the
12 Free Software Foundation; either version 2, or (at your option) any
13 later version.
14
15 This program is distributed in the hope that it will be useful,
16 but WITHOUT ANY WARRANTY; without even the implied warranty of
17 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 GNU General Public License for more details.
19
20 You should have received a copy of the GNU General Public License
21 along with this program; if not, write to the Free Software
22 Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
23 USA. */
24
25/*===========================================================================
26 === File description
27 ===========================================================================*/
28
39/*===========================================================================
40 === External declarations
41 ===========================================================================*/
42#include <cmath>
43#include <cstdlib>
44#include <stdexcept>
45
46#include "array.h"
47#include "arts.h"
48#include "auto_md.h"
49#include "check_input.h"
50#include "cloudbox.h"
51#include "file.h"
52#include "gridded_fields.h"
53#include "interpolation.h"
55#include "lin_alg.h"
56#include "logic.h"
57#include "math_funcs.h"
58#include "messages.h"
59#include "microphysics.h"
60#include "optproperties.h"
61#include "parameters.h"
62#include "physics_funcs.h"
63#include "rte.h"
64#include "sorting.h"
65#include "special_interp.h"
66#include "xml_io.h"
67
68extern const Index GFIELD3_P_GRID;
69extern const Index GFIELD3_LAT_GRID;
70extern const Index GFIELD3_LON_GRID;
71extern const Numeric PI;
72extern const Numeric DEG2RAD;
73extern const Numeric RAD2DEG;
74extern const Numeric LAT_LON_MIN;
75extern const Numeric DENSITY_OF_ICE;
76
77/*===========================================================================
78 === The functions (in alphabetical order)
79 ===========================================================================*/
80
81/* Workspace method: Doxygen documentation will be auto-generated */
82void cloudboxOff(Index& cloudbox_on,
83 Index& ppath_inside_cloudbox_do,
84 ArrayOfIndex& cloudbox_limits,
85 Agenda& iy_cloudbox_agenda,
86 Tensor4& pnd_field,
87 ArrayOfTensor4& dpnd_field_dx,
88 ArrayOfString& scat_species,
91 Index& scat_data_checked,
92 Matrix& particle_masses,
93 const ArrayOfRetrievalQuantity& jacobian_quantities,
94 const Verbosity&) {
95 cloudbox_on = 0;
96 ppath_inside_cloudbox_do = 0;
97 cloudbox_limits.resize(0);
98 iy_cloudbox_agenda = Agenda();
99 iy_cloudbox_agenda.set_name("iy_cloudbox_agenda");
100 pnd_field.resize(0, 0, 0, 0);
101 // we need to size dpnd_field to be consistent with jacobian_quantities.
102 dpnd_field_dx.resize(jacobian_quantities.nelem());
103 scat_data.resize(0);
104 scat_species.resize(0);
105 // remove scat_data_raw resizing once all scat solvers have been convert to
106 // use of (new-type) scat_data
107 scat_data_raw.resize(0);
108 scat_data_checked = 0;
109 particle_masses.resize(0, 0);
110}
111
112/* Workspace method: Doxygen documentation will be auto-generated */
113void cloudboxSetAutomatically( // WS Output:
114 //Workspace& /* ws */,
115 Index& cloudbox_on,
116 ArrayOfIndex& cloudbox_limits,
117 //Agenda& iy_cloudbox_agenda,
118 // WS Input:
119 const Index& atmosphere_dim,
120 const Vector& p_grid,
121 const Vector& lat_grid,
122 const Vector& lon_grid,
123 const Tensor4& particle_field,
124 // Control Parameters
125 const Numeric& cloudbox_margin,
126 const Verbosity& verbosity) {
127 // Check existing WSV
128 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
129 // includes p_grid chk_if_decreasing
130 chk_atm_grids(atmosphere_dim, p_grid, lat_grid, lon_grid);
131 // Set cloudbox_on
132 cloudbox_on = 1;
133
134 ARTS_USER_ERROR_IF (atmosphere_dim > 1,
135 "cloudboxSetAutomatically not yet available for 2D and 3D cases.")
136
137 Index np = p_grid.nelem();
138
139 // Allocate cloudbox_limits
140 cloudbox_limits.resize(atmosphere_dim * 2);
141
142 // Initialize boundary counters
143 Index p1 = np - 1, p2 = 0;
144
145 // OLE: Commented out until code that uses it at the end of this function is commented back in
146 /*
147 if ( atmosphere_dim > 1 )
148 {
149 Index lat1 = particle_field.nrows()-1;
150 Index lat2 = 0;
151 }
152 if ( atmosphere_dim > 2 )
153 {
154 Index lon1 = particle_field.ncols()-1;
155 Index lon2 = 0;
156 }
157*/
158
159 bool any_not_empty = false;
160
161 if (!particle_field.empty()) {
162 bool one_not_empty = false;
163 Index nss = particle_field.nbooks();
164
165 //--------- Start loop over fields in particle_field------------------------
166 for (Index l = 0; l < nss; l++) {
167 //not_empty is set to true, if a single value of particle_field
168 //is unequal zero (and not NaN), i.e. if we actually have some amount of
169 //these scattering species in the atmosphere.
170 chk_scat_species_field(one_not_empty,
171 particle_field(l, joker, joker, joker),
172 "particle_field",
173 atmosphere_dim,
174 p_grid,
175 lat_grid,
176 lon_grid);
177
178 //if particles found, enter detailed search
179 if (one_not_empty) {
180 any_not_empty = true;
182 p2,
183 particle_field(l, joker, joker, joker),
184 atmosphere_dim,
185 cloudbox_margin);
186 }
187 }
188 }
189
190 if (any_not_empty) {
191 // decrease lower cb limit by one to ensure that linear interpolation of
192 // particle number densities is possible.
193 p1 = max(p1 - 1, Index(0));
194
195 Numeric p_margin1;
196
197 // alter lower cloudbox_limit by cloudbox_margin, using barometric
198 // height formula
199 p_margin1 = barometric_heightformula(p_grid[p1], cloudbox_margin);
200 while ((p_grid[p1] < p_margin1) && (p1 > 0)) p1--;
201 cloudbox_limits[0] = p1;
202
203 // increase upper cb limit by one to ensure that linear interpolation of
204 // particle number densities is possible.
205 p2 = min(p2 + 1, np - 1);
206 // set upper cloudbox_limit
207 // if cloudbox reaches to the upper most pressure level
208 if (p2 >= np - 1) {
210 out2 << "The cloud reaches to TOA!\n"
211 << "Check your *particle_field* data, if realistic!\n";
212 }
213 cloudbox_limits[1] = p2;
214
215 // ARTS_ASSERT keyword arguments
216
217 // The pressure in *p1* must be bigger than the pressure in *p2*.
218 ARTS_ASSERT(p_grid[p1] > p_grid[p2]);
219 // The pressure in *p1* must be larger than the last value in *p_grid*.
220 ARTS_ASSERT(p_grid[p1] > p_grid[p_grid.nelem() - 1]);
221 // The pressure in *p2* must be smaller than the first value in *p_grid*."
222 ARTS_ASSERT(p_grid[p2] < p_grid[0]);
223
224 /*
225 if ( atmosphere_dim >= 2 )
226 {
227 // The latitude in *lat2* must be bigger than the latitude in *lat1*.
228 ARTS_ASSERT ( lat_grid[lat2] > lat_grid[lat1] );
229 // The latitude in *lat1* must be >= the second value in *lat_grid*.
230 ARTS_ASSERT ( lat_grid[lat1] >= lat_grid[1] );
231 // The latitude in *lat2* must be <= the next to last value in *lat_grid*.
232 ARTS_ASSERT ( lat_grid[lat2] <= lat_grid[lat_grid.nelem()-2] );
233 }
234 if ( atmosphere_dim == 3 )
235 {
236 // The longitude in *lon2* must be bigger than the longitude in *lon1*.
237 ARTS_ASSERT ( lon_grid[lon2] > lon_grid[lon1] );
238 // The longitude in *lon1* must be >= the second value in *lon_grid*.
239 ARTS_ASSERT ( lon_grid[lon1] >= lon_grid[1] );
240 // The longitude in *lon2* must be <= the next to last value in *lon_grid*.
241 ARTS_ASSERT ( lon_grid[lon2] <= lon_grid[lon_grid.nelem()-2] );
242 }
243 */
244 }
245
246 else
247 // if all particle fields are zero at each level and cloudbox was not preset,
248 // switch cloudbox off.
249 {
251 cloudbox_on = 0;
252 cloudbox_limits[1] =
253 -1; // just for consistency with cloudboxSetAutomatically
254 out0 << "Cloudbox is switched off!\n";
255 }
256}
257
258/* Workspace method: Doxygen documentation will be auto-generated */
259void cloudboxSetFullAtm( //WS Output
260 Index& cloudbox_on,
261 ArrayOfIndex& cloudbox_limits,
262 // WS Input
263 const Index& atmosphere_dim,
264 const Vector& p_grid,
265 const Vector& lat_grid,
266 const Vector& lon_grid,
267 const Verbosity&) {
268 cloudbox_on = 1;
269 cloudbox_limits.resize(2 * atmosphere_dim);
270
271 cloudbox_limits[0] = 0;
272 cloudbox_limits[1] = p_grid.nelem() - 1;
273
274 if (atmosphere_dim > 1) {
275 Index last_lat = lat_grid.nelem() - 1;
276
277 // find minimum lat_grid point i with lat_grid[i]-lat_grid[0]>=LAT_LON_MIN
278 Index i = 1;
279 while ((i < last_lat - 1) && (lat_grid[i] - lat_grid[0] < LAT_LON_MIN)) i++;
280 ARTS_USER_ERROR_IF (i == last_lat - 1,
281 "Can not define lower latitude edge of cloudbox:\n"
282 "Extend of atmosphere too small. Distance to minimum latitude\n"
283 "has to be at least ", LAT_LON_MIN, "deg, but only ",
284 lat_grid[i - 1] - lat_grid[0], " available here.")
285 cloudbox_limits[2] = i;
286
287 // find maximum lat_grid point j with lat_grid[-1]-lat_grid[j]>=LAT_LON_MIN
288 // and j>i
289 Index j = last_lat - 1;
290 while ((j > i) && (lat_grid[last_lat] - lat_grid[j] < LAT_LON_MIN)) j--;
291 ARTS_USER_ERROR_IF (j == i,
292 "Can not define upper latitude edge of cloudbox:\n"
293 "Extend of atmosphere too small. Distance to maximum latitude\n"
294 "has to be at least ", LAT_LON_MIN, "deg, but only ",
295 lat_grid[last_lat] - lat_grid[j + 1], " available here.")
296 cloudbox_limits[3] = j;
297 }
298
299 if (atmosphere_dim > 2) {
300 const Numeric latmax = max(abs(lat_grid[cloudbox_limits[2]]),
301 abs(lat_grid[cloudbox_limits[3]]));
302 const Numeric lfac = 1 / cos(DEG2RAD * latmax);
303 Index last_lon = lon_grid.nelem() - 1;
304
305 // find minimum lon_grid point i with lon_grid[i]-lon_grid[0]>=LAT_LON_MIN/lfac
306 Index i = 1;
307 while ((i < last_lon - 1) &&
308 (lon_grid[i] - lon_grid[0] < LAT_LON_MIN / lfac))
309 i++;
310 ARTS_USER_ERROR_IF (i == last_lon - 1,
311 "Can not define lower longitude edge of cloudbox:\n"
312 "Extend of atmosphere too small. Distance to minimum longitude\n"
313 "has to be at least ", LAT_LON_MIN / lfac, "deg, but only ",
314 lon_grid[i - 1] - lon_grid[0], " available here.")
315 cloudbox_limits[4] = i;
316
317 // find maximum lon_grid point j with lon_grid[-1]-lon_grid[j]>=LAT_LON_MIN/lfac
318 // and j>i
319 Index j = last_lon - 1;
320 while ((j > i) && (lon_grid[last_lon] - lon_grid[j] < LAT_LON_MIN / lfac))
321 j--;
322 ARTS_USER_ERROR_IF (j == i,
323 "Can not define upper longitude edge of cloudbox:\n"
324 "Extend of atmosphere too small. Distance to maximum longitude\n"
325 "has to be at least ", LAT_LON_MIN / lfac, "deg, but only ",
326 lon_grid[last_lon] - lon_grid[j + 1], " available here.")
327 cloudbox_limits[5] = j;
328 }
329}
330
331/* Workspace method: Doxygen documentation will be auto-generated */
332void cloudboxSetManually( // WS Output:
333 Index& cloudbox_on,
334 ArrayOfIndex& cloudbox_limits,
335 // WS Input:
336 const Index& atmosphere_dim,
337 const Vector& p_grid,
338 const Vector& lat_grid,
339 const Vector& lon_grid,
340 // Control Parameters:
341 const Numeric& p1,
342 const Numeric& p2,
343 const Numeric& lat1,
344 const Numeric& lat2,
345 const Numeric& lon1,
346 const Numeric& lon2,
347 const Verbosity&) {
348 // Check existing WSV
349 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
350 chk_atm_grids(atmosphere_dim, p_grid, lat_grid, lon_grid);
351
352 // Check keyword arguments
353 ARTS_USER_ERROR_IF (p1 <= p2,
354 "The pressure in *p1* must be bigger than the "
355 "pressure in *p2*.");
356 ARTS_USER_ERROR_IF (p1 <= p_grid[p_grid.nelem() - 1],
357 "The pressure in *p1* must be larger than the "
358 "last value in *p_grid*.");
359 ARTS_USER_ERROR_IF (p2 >= p_grid[0],
360 "The pressure in *p2* must be smaller than the "
361 "first value in *p_grid*.");
362 if (atmosphere_dim >= 2) {
363 ARTS_USER_ERROR_IF (lat2 <= lat1,
364 "The latitude in *lat2* must be bigger than the "
365 "latitude in *lat1*.");
366 ARTS_USER_ERROR_IF (lat1 < lat_grid[1],
367 "The latitude in *lat1* must be >= the "
368 "second value in *lat_grid*.");
369 ARTS_USER_ERROR_IF (lat2 > lat_grid[lat_grid.nelem() - 2],
370 "The latitude in *lat2* must be <= the "
371 "next to last value in *lat_grid*.");
372 }
373 if (atmosphere_dim == 3) {
374 ARTS_USER_ERROR_IF (lon2 <= lon1,
375 "The longitude in *lon2* must be bigger than the "
376 "longitude in *lon1*.");
377 ARTS_USER_ERROR_IF (lon1 < lon_grid[1],
378 "The longitude in *lon1* must be >= the "
379 "second value in *lon_grid*.");
380 ARTS_USER_ERROR_IF (lon2 > lon_grid[lon_grid.nelem() - 2],
381 "The longitude in *lon2* must be <= the "
382 "next to last value in *lon_grid*.");
383 }
384
385 // Set cloudbox_on
386 cloudbox_on = 1;
387
388 // Allocate cloudbox_limits
389 cloudbox_limits.resize(atmosphere_dim * 2);
390
391 // Pressure limits
392 if (p1 > p_grid[1]) {
393 cloudbox_limits[0] = 0;
394 } else {
395 for (cloudbox_limits[0] = 1; p_grid[cloudbox_limits[0] + 1] >= p1;
396 cloudbox_limits[0]++) {
397 }
398 }
399 if (p2 < p_grid[p_grid.nelem() - 2]) {
400 cloudbox_limits[1] = p_grid.nelem() - 1;
401 } else {
402 for (cloudbox_limits[1] = p_grid.nelem() - 2;
403 p_grid[cloudbox_limits[1] - 1] <= p2;
404 cloudbox_limits[1]--) {
405 }
406 }
407
408 // Latitude limits
409 if (atmosphere_dim >= 2) {
410 for (cloudbox_limits[2] = 1; lat_grid[cloudbox_limits[2] + 1] <= lat1;
411 cloudbox_limits[2]++) {
412 }
413 for (cloudbox_limits[3] = lat_grid.nelem() - 2;
414 lat_grid[cloudbox_limits[3] - 1] >= lat2;
415 cloudbox_limits[3]--) {
416 }
417 }
418
419 // Longitude limits
420 if (atmosphere_dim == 3) {
421 for (cloudbox_limits[4] = 1; lon_grid[cloudbox_limits[4] + 1] <= lon1;
422 cloudbox_limits[4]++) {
423 }
424 for (cloudbox_limits[5] = lon_grid.nelem() - 2;
425 lon_grid[cloudbox_limits[5] - 1] >= lon2;
426 cloudbox_limits[5]--) {
427 }
428 }
429}
430
431/* Workspace method: Doxygen documentation will be auto-generated */
433 Index& cloudbox_on,
434 ArrayOfIndex& cloudbox_limits,
435 // WS Input:
436 const Index& atmosphere_dim,
437 const Tensor3& z_field,
438 const Vector& lat_grid,
439 const Vector& lon_grid,
440 // Control Parameters:
441 const Numeric& z1,
442 const Numeric& z2,
443 const Numeric& lat1,
444 const Numeric& lat2,
445 const Numeric& lon1,
446 const Numeric& lon2,
447 const Verbosity&) {
448 // Check existing WSV
449 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
450
451 // Check keyword arguments
452 ARTS_USER_ERROR_IF (z1 >= z2,
453 "The altitude in *z1* must be smaller than the "
454 "altitude in *z2*.");
455 /* These cases are in fact handled
456 if( z1 <= z_field(0, 0, 0) )
457 throw runtime_error( "The altitude in *z1* must be larger than the "
458 "first value in *z_field*." );
459 if( z2 >= z_field(z_field.npages()-1, 0, 0) )
460 throw runtime_error( "The altitude in *z2* must be smaller than the "
461 "last value in *z_field*." );
462 */
463 if (atmosphere_dim == 3) {
464 ARTS_USER_ERROR_IF (lat2 <= lat1,
465 "The latitude in *lat2* must be bigger than the "
466 " latitude in *lat1*.");
467 ARTS_USER_ERROR_IF (lat1 < lat_grid[1],
468 "The latitude in *lat1* must be >= the "
469 "second value in *lat_grid*.");
470 ARTS_USER_ERROR_IF (lat2 > lat_grid[lat_grid.nelem() - 2],
471 "The latitude in *lat2* must be <= the "
472 "next to last value in *lat_grid*.");
473 ARTS_USER_ERROR_IF (lon2 <= lon1,
474 "The longitude in *lon2* must be bigger than the "
475 "longitude in *lon1*.");
476 ARTS_USER_ERROR_IF (lon1 < lon_grid[1],
477 "The longitude in *lon1* must be >= the "
478 "second value in *lon_grid*.");
479 ARTS_USER_ERROR_IF (lon2 > lon_grid[lon_grid.nelem() - 2],
480 "The longitude in *lon2* must be <= the "
481 "next to last value in *lon_grid*.");
482 }
483
484 // Set cloudbox_on
485 cloudbox_on = 1;
486
487 // Allocate cloudbox_limits
488 cloudbox_limits.resize(atmosphere_dim * 2);
489
490 // Pressure/altitude limits
491 if (z1 < z_field(1, 0, 0)) {
492 cloudbox_limits[0] = 0;
493 } else {
494 for (cloudbox_limits[0] = 1; z_field(cloudbox_limits[0] + 1, 0, 0) <= z1;
495 cloudbox_limits[0]++) {
496 }
497 }
498 if (z2 > z_field(z_field.npages() - 2, 0, 0)) {
499 cloudbox_limits[1] = z_field.npages() - 1;
500 } else {
501 for (cloudbox_limits[1] = z_field.npages() - 2;
502 z_field(cloudbox_limits[1] - 1, 0, 0) >= z2;
503 cloudbox_limits[1]--) {
504 }
505 }
506
507 // Latitude limits
508 if (atmosphere_dim >= 2) {
509 for (cloudbox_limits[2] = 1; lat_grid[cloudbox_limits[2] + 1] <= lat1;
510 cloudbox_limits[2]++) {
511 }
512 for (cloudbox_limits[3] = lat_grid.nelem() - 2;
513 lat_grid[cloudbox_limits[3] - 1] >= lat2;
514 cloudbox_limits[3]--) {
515 }
516 }
517
518 // Longitude limits
519 if (atmosphere_dim == 3) {
520 for (cloudbox_limits[4] = 1; lon_grid[cloudbox_limits[4] + 1] <= lon1;
521 cloudbox_limits[4]++) {
522 }
523 for (cloudbox_limits[5] = lon_grid.nelem() - 2;
524 lon_grid[cloudbox_limits[5] - 1] >= lon2;
525 cloudbox_limits[5]--) {
526 }
527 }
528}
529
530/* Workspace method: Doxygen documentation will be auto-generated */
532 const Tensor7& cloudbox_field,
533 const Vector& rte_pos,
534 const Vector& rte_los,
535 const Index& jacobian_do,
536 const Index& cloudbox_on,
537 const ArrayOfIndex& cloudbox_limits,
538 const Index& atmosphere_dim,
539 const Vector& p_grid,
540 const Vector& lat_grid,
541 const Vector& lon_grid,
542 const Tensor3& z_field,
543 const Matrix& z_surface,
544 const Index& stokes_dim,
545 const Vector& za_grid,
546 const Vector& aa_grid,
547 const Vector& f_grid,
548 //const Index& p_interp_order,
549 const Index& za_interp_order,
550 const Index& za_restrict,
551 const Index& cos_za_interp,
552 const Numeric& za_extpolfac,
553 const Index& aa_interp_order,
554 const Verbosity&) {
555 //--- Check input -----------------------------------------------------------
556 ARTS_USER_ERROR_IF (!(atmosphere_dim == 1 || atmosphere_dim == 3),
557 "The atmospheric dimensionality must be 1 or 3.");
558 ARTS_USER_ERROR_IF (jacobian_do,
559 "This method does not support jacobians (jacobian_do must be 0)");
560 ARTS_USER_ERROR_IF (!cloudbox_on,
561 "The cloud box is not activated and no outgoing "
562 "field can be returned.");
563 ARTS_USER_ERROR_IF (cloudbox_limits.nelem() != 2 * atmosphere_dim,
564 "*cloudbox_limits* is a vector which contains the upper and lower\n"
565 "limit of the cloud for all atmospheric dimensions.\n"
566 "So its length must be 2 x *atmosphere_dim*");
567 ARTS_USER_ERROR_IF (za_grid.nelem() == 0,
568 "The variable *za_grid* is empty. Are dummy "
569 "values from *cloudboxOff used?");
570 ARTS_USER_ERROR_IF (!(za_interp_order < za_grid.nelem()),
571 "Zenith angle interpolation order *za_interp_order*"
572 " must be smaller\n"
573 "than number of angles in *za_grid*.");
574 ARTS_USER_ERROR_IF (atmosphere_dim > 1 && !(aa_interp_order < aa_grid.nelem()),
575 "Azimuth angle interpolation order *aa_interp_order*"
576 " must be smaller\n"
577 "than number of angles in *aa_grid*.");
578 ARTS_USER_ERROR_IF (cloudbox_field.nlibraries() != f_grid.nelem(),
579 "Inconsistency in size between f_grid and cloudbox_field! "
580 "(This method does not yet handle dispersion type calculations.)");
581 //---------------------------------------------------------------------------
582
583 // At each lat/lon, one value of the intensity field is defined at the
584 // surface. Simplest way to handle this is to insert z_surface into z_field
585 Tensor3 z_with_surface = z_field;
586 for (Index ilat = 0; ilat < z_surface.nrows(); ilat++) {
587 for (Index ilon = 0; ilon < z_surface.ncols(); ilon++) {
588 Index ip = 0;
589 while (z_surface(ilat, ilon) >= z_field(ip + 1, ilat, ilon)) {
590 ip++;
591 }
592 z_with_surface(ip, ilat, ilon) = z_surface(ilat, ilon);
593 }
594 }
595
596 // Convert rte_pos to grid positions
597 GridPos gp_p, gp_lat, gp_lon;
598 rte_pos2gridpos(gp_p,
599 gp_lat,
600 gp_lon,
601 atmosphere_dim,
602 p_grid,
603 lat_grid,
604 lon_grid,
605 z_with_surface,
606 rte_pos);
607
608 //--- Determine if at border or inside of cloudbox (or outside!)
609 //
610 // Let us introduce a number coding for cloud box borders.
611 // Borders have the same number as position in *cloudbox_limits*.
612 // Inside cloud box is coded as 99, and outside as > 100.
613 Index border = 999;
614 //
615 //- Check if at any border
616 if (is_gridpos_at_index_i(gp_p, cloudbox_limits[0], false)) {
617 border = 0;
618 } else if (is_gridpos_at_index_i(gp_p, cloudbox_limits[1], false)) {
619 border = 1;
620 } else if (atmosphere_dim > 1) {
621 if (is_gridpos_at_index_i(gp_lat, cloudbox_limits[2], false)) {
622 border = 2;
623 } else if (is_gridpos_at_index_i(gp_lat, cloudbox_limits[3], false)) {
624 border = 3;
625 } else if (atmosphere_dim > 2) {
626 if (is_gridpos_at_index_i(gp_lon, cloudbox_limits[4], false)) {
627 border = 4;
628 } else if (is_gridpos_at_index_i(gp_lon, cloudbox_limits[5], false)) {
629 border = 5;
630 }
631 }
632 }
633
634 //
635 //- Check if inside (when border<100 here, it means we are on a cloudbox border)
636 if (border > 100) {
637 // Assume inside as it is easiest to detect if outside (can be detected
638 // check in one dimension at the time)
639 bool inside = true;
640 Numeric fgp;
641
642 // Check in pressure dimension
643 fgp = fractional_gp(gp_p);
644 if (fgp < Numeric(cloudbox_limits[0]) ||
645 fgp > Numeric(cloudbox_limits[1])) {
646 inside = false;
647 }
648
649 // Check in lat and lon dimensions
650 if (atmosphere_dim == 3 && inside) {
651 fgp = fractional_gp(gp_lat);
652 if (fgp < Numeric(cloudbox_limits[2]) ||
653 fgp > Numeric(cloudbox_limits[3])) {
654 inside = false;
655 }
656 fgp = fractional_gp(gp_lon);
657 if (fgp < Numeric(cloudbox_limits[4]) ||
658 fgp > Numeric(cloudbox_limits[5])) {
659 inside = false;
660 }
661 }
662
663 if (inside) {
664 border = 99;
665 }
666 }
667
668 // If outside, something is wrong
669 ARTS_USER_ERROR_IF (border > 100,
670 "Given position has been found to be outside the cloudbox.");
671
672 //- Sizes
673 const Index nf = f_grid.nelem();
674 DEBUG_ONLY(const Index np = cloudbox_limits[1] - cloudbox_limits[0] + 1);
675 const Index nza = za_grid.nelem();
676 const Index naa = cloudbox_field.nrows();
677
678 //- Resize *iy*
679 iy.resize(nf, stokes_dim);
680 //- Make space for spatially interpolated field (we perfrom angle
681 //interpolation separately on this then)
682 Tensor4 i_field_local(nf, nza, naa, stokes_dim);
683
684 // Index of the p/lat/lon slice (first or last element slice in this
685 // dimension), where rte_pos is located. if located on a border.
686 Index border_index;
687 if (border < 99) {
688 if (border % 2) // odd border id, ie top or north or east
689 border_index = cloudbox_limits[border] - cloudbox_limits[border - 1];
690 else // even border id, ie bottom or south or west
691 border_index = 0;
692 }
693
694 // Sensor points inside the cloudbox
695 if (border == 99) {
696 ARTS_USER_ERROR_IF (atmosphere_dim == 3,
697 "Radiation extraction for a position inside cloudbox\n"
698 "is not yet implemented for 3D cases.\n");
699 ARTS_ASSERT(atmosphere_dim == 1);
700
701 ARTS_ASSERT(is_size(cloudbox_field, nf, np, 1, 1, nza, 1, stokes_dim));
702
703 // Grid position in *p_grid*
704 gp_p.idx = gp_p.idx - cloudbox_limits[0];
705 gridpos_upperend_check(gp_p, cloudbox_limits[1] - cloudbox_limits[0]);
706 Vector itw_p(2);
707 interpweights(itw_p, gp_p);
708
709 for (Index is = 0; is < stokes_dim; is++)
710 for (Index iv = 0; iv < nf; iv++)
711 for (Index i_za = 0; i_za < nza; i_za++)
712 i_field_local(iv, i_za, 0, is) = interp(
713 itw_p, cloudbox_field(iv, joker, 0, 0, i_za, 0, is), gp_p);
714 }
715
716 // Sensor outside the cloudbox
717
718 // --- 1D ------------------------------------------------------------------
719 else if (atmosphere_dim == 1) {
720 i_field_local =
721 cloudbox_field(joker, border_index, 0, 0, joker, joker, joker);
722 }
723
724 // --- 3D ------------------------------------------------------------------
725 else {
726 ARTS_ASSERT(is_size(cloudbox_field,
727 nf,
728 cloudbox_field.nvitrines(),
729 cloudbox_field.nshelves(),
730 cloudbox_field.nbooks(),
731 za_grid.nelem(),
732 aa_grid.nelem(),
733 stokes_dim));
734
735 // Interpolation weights (for 2D "red" interpolation)
736 Vector itw(4);
737
738 // Outgoing from cloudbox top or bottom, i.e. from a pressure level
739 if (border <= 1) {
740 // Lat and lon grid positions with respect to cloud box
741 GridPos cb_gp_lat, cb_gp_lon;
742 cb_gp_lat = gp_lat;
743 cb_gp_lon = gp_lon;
744 cb_gp_lat.idx -= cloudbox_limits[2];
745 cb_gp_lon.idx -= cloudbox_limits[4];
746 //
747 gridpos_upperend_check(cb_gp_lat,
748 cloudbox_limits[3] - cloudbox_limits[2]);
749 gridpos_upperend_check(cb_gp_lon,
750 cloudbox_limits[5] - cloudbox_limits[4]);
751
752 interpweights(itw, cb_gp_lat, cb_gp_lon);
753
754 for (Index is = 0; is < stokes_dim; is++)
755 for (Index iv = 0; iv < nf; iv++)
756 for (Index i_za = 0; i_za < nza; i_za++)
757 for (Index i_aa = 0; i_aa < naa; i_aa++)
758 i_field_local(iv, i_za, i_aa, is) =
759 interp(itw,
760 cloudbox_field(
761 iv, border_index, joker, joker, i_za, i_aa, is),
762 cb_gp_lat,
763 cb_gp_lon);
764 }
765
766 // Outgoing from cloudbox north or south border, i.e. from a latitude level
767 else if (border <= 3) {
768 // Pressure and lon grid positions with respect to cloud box
769 GridPos cb_gp_p, cb_gp_lon;
770 cb_gp_p = gp_p;
771 cb_gp_lon = gp_lon;
772 cb_gp_p.idx -= cloudbox_limits[0];
773 cb_gp_lon.idx -= cloudbox_limits[4];
774 //
775 gridpos_upperend_check(cb_gp_p, cloudbox_limits[1] - cloudbox_limits[0]);
776 gridpos_upperend_check(cb_gp_lon,
777 cloudbox_limits[5] - cloudbox_limits[4]);
778
779 interpweights(itw, cb_gp_p, cb_gp_lon);
780
781 for (Index is = 0; is < stokes_dim; is++)
782 for (Index iv = 0; iv < nf; iv++)
783 for (Index i_za = 0; i_za < nza; i_za++)
784 for (Index i_aa = 0; i_aa < naa; i_aa++)
785 i_field_local(iv, i_za, i_aa, is) =
786 interp(itw,
787 cloudbox_field(
788 iv, joker, border_index, joker, i_za, i_aa, is),
789 cb_gp_p,
790 cb_gp_lon);
791 }
792
793 // Outgoing from cloudbox east or west border, i.e. from a longitude level
794 else {
795 // Pressure and lat grid positions with respect to cloud box
796 GridPos cb_gp_p, cb_gp_lat;
797 cb_gp_p = gp_p;
798 cb_gp_lat = gp_lat;
799 cb_gp_p.idx -= cloudbox_limits[0];
800 cb_gp_lat.idx -= cloudbox_limits[2];
801 //
802 gridpos_upperend_check(cb_gp_p, cloudbox_limits[1] - cloudbox_limits[0]);
803 gridpos_upperend_check(cb_gp_lat,
804 cloudbox_limits[3] - cloudbox_limits[2]);
805
806 interpweights(itw, cb_gp_p, cb_gp_lat);
807
808 for (Index is = 0; is < stokes_dim; is++)
809 for (Index iv = 0; iv < nf; iv++)
810 for (Index i_za = 0; i_za < nza; i_za++)
811 for (Index i_aa = 0; i_aa < naa; i_aa++)
812 i_field_local(iv, i_za, i_aa, is) =
813 interp(itw,
814 cloudbox_field(
815 iv, joker, joker, border_index, i_za, i_aa, is),
816 cb_gp_p,
817 cb_gp_lat);
818 }
819 }
820
821 // now, do Nth-oder polynomial interpoaltion of angle(s)
822 //
823 // a bunch of options for testing:
824 // a) interpolation over full zenith angle grid vs over hemispheres
825 // separately.
826 // b) interpolation in plain zenith angles vs. in cosines of zenith angle.
827
828 // find range of za_grid that we will do interpolation over.
829 Index za_start = 0;
830 Index za_extend = za_grid.nelem();
831 if (za_restrict) {
832 // which hemisphere do we need?
833 ARTS_USER_ERROR_IF (is_same_within_epsilon(rte_los[0], 90., 1e-6),
834 //FIXME: we should allow this though, if za_grid has a grid point
835 //at 90deg.
836 "Hemisphere-restricted zenith angle interpolation not allowed\n"
837 "for 90degree views.");
838 if (rte_los[0] > 90) {
839 // upwelling, i.e. second part of za_grid. that is, we need to find
840 // the first point in za_grid where za>90. and update za_start
841 // accordingly.
842 while (za_start < za_grid.nelem() && za_grid[za_start] < 90.) {
843 za_start++;
844 }
845 ARTS_USER_ERROR_IF (za_start == za_grid.nelem(),
846 "No za_grid grid point found in 90-180deg hemisphere.\n"
847 "No hemispheric interpolation possible.");
848 za_extend -= za_start;
849 } else {
850 // downwelling, i.e. first part of za_grid. that is, we need to
851 // find the last point in za_grid where za<90. and update za_extend
852 // accordingly.
853 while (za_extend > 0 && za_grid[za_extend - 1] > 90.) {
854 za_extend--;
855 }
856 ARTS_USER_ERROR_IF (za_extend == 0,
857 "No za_grid grid point found in 0-90deg hemisphere.\n"
858 "No hemispheric interpolation possible.");
859 }
860 ARTS_USER_ERROR_IF (!(za_interp_order < za_extend),
861 "Zenith angle interpolation order *za_interp_order* (",
862 za_interp_order, ") must be smaller\n"
863 "than number of angles in respective hemisphere (", za_extend,
864 ").")
865 }
866
867 const auto range = Range(za_start, za_extend);
868 const ConstVectorView za_g(za_grid[range]);
869
870 // Check the zenith-grid
871 Interpolation::check_lagrange_interpolation(za_g, za_interp_order, rte_los[0], za_extpolfac);
872
873 // Find position of zenith, either by cosine or linear weights
874 const auto lag_za = cos_za_interp ?
875 LagrangeInterpolation(0, rte_los[0], za_g, za_interp_order, false, Interpolation::GridType::CosDeg) :
876 LagrangeInterpolation(0, rte_los[0], za_g, za_interp_order);
877
878 // First position if 1D atmosphere, otherwise compute cyclic for a azimuth grid [-180, 180]
879 const auto lag_aa = (atmosphere_dim > 1) ?
880 LagrangeInterpolation(0, rte_los[1], aa_grid, aa_interp_order, false, Interpolation::GridType::Cyclic, {-180, 180}) :
882
883 // Corresponding interpolation weights
884 const auto itw_angs=interpweights(lag_za, lag_aa);
885
886 for (Index is = 0; is < stokes_dim; is++) {
887 for (Index iv = 0; iv < nf; iv++) {
888 iy(iv, is) =
889 interp(i_field_local(iv, range, joker, is),
890 itw_angs,
891 lag_za,
892 lag_aa);
893 }
894 }
895}
896
897/* Workspace method: Doxygen documentation will be auto-generated */
898void cloudbox_fieldCrop(Tensor7& cloudbox_field,
899 ArrayOfIndex& cloudbox_limits,
900 const Index& atmosphere_dim,
901 const Index& cloudbox_on,
902 const Index& new_limit0,
903 const Index& new_limit1,
904 const Index& new_limit2,
905 const Index& new_limit3,
906 const Index& new_limit4,
907 const Index& new_limit5,
908 const Verbosity&) {
909 ARTS_USER_ERROR_IF (!cloudbox_on,
910 "No need to use this method with cloudbox=0.");
911 ARTS_USER_ERROR_IF (new_limit0 < cloudbox_limits[0],
912 "new_limits0 < cloudbox_limits[0], not allowed!");
913 ARTS_USER_ERROR_IF (new_limit1 > cloudbox_limits[1],
914 "new_limits1 > cloudbox_limits[1], not allowed!");
915
916 Tensor7 fcopy = cloudbox_field;
917
918 if (atmosphere_dim == 1) {
919 cloudbox_field = fcopy(
920 joker,
921 Range(new_limit0 - cloudbox_limits[0], new_limit1 - new_limit0 + 1),
922 joker,
923 joker,
924 joker,
925 joker,
926 joker);
927 cloudbox_limits[0] = new_limit0;
928 cloudbox_limits[1] = new_limit1;
929 } else {
930 ARTS_USER_ERROR_IF (new_limit2 < cloudbox_limits[2],
931 "new_limits2 < cloudbox_limits[2], not allowed!");
932 ARTS_USER_ERROR_IF (new_limit3 > cloudbox_limits[3],
933 "new_limits3 > cloudbox_limits[3], not allowed!");
934
935 if (atmosphere_dim == 2) {
936 cloudbox_field = fcopy(
937 joker,
938 Range(new_limit0 - cloudbox_limits[0], new_limit1 - new_limit0 + 1),
939 Range(new_limit2 - cloudbox_limits[2], new_limit3 - new_limit2 - 1),
940 joker,
941 joker,
942 joker,
943 joker);
944 cloudbox_limits[0] = new_limit0;
945 cloudbox_limits[1] = new_limit1;
946 cloudbox_limits[2] = new_limit2;
947 cloudbox_limits[3] = new_limit3;
948 } else {
949 ARTS_USER_ERROR_IF (new_limit4 < cloudbox_limits[4],
950 "new_limits4 < cloudbox_limits[4], not allowed!");
951 ARTS_USER_ERROR_IF (new_limit5 > cloudbox_limits[5],
952 "new_limits5 > cloudbox_limits[5], not allowed!");
953 cloudbox_field = fcopy(
954 joker,
955 Range(new_limit0 - cloudbox_limits[0], new_limit1 - new_limit0 + 1),
956 Range(new_limit2 - cloudbox_limits[2], new_limit3 - new_limit2 + 1),
957 Range(new_limit4 - cloudbox_limits[4], new_limit5 - new_limit4 + 1),
958 joker,
959 joker,
960 joker);
961 cloudbox_limits[0] = new_limit0;
962 cloudbox_limits[1] = new_limit1;
963 cloudbox_limits[2] = new_limit2;
964 cloudbox_limits[3] = new_limit3;
965 cloudbox_limits[4] = new_limit4;
966 cloudbox_limits[5] = new_limit5;
967 }
968 }
969}
970
971/* Workspace method: Doxygen documentation will be auto-generated */
972void particle_fieldCleanup( //WS Output:
973 Tensor4& particle_field_out,
974 //WS Input:
975 const Tensor4& particle_field_in,
976 const Numeric& threshold,
977 const Verbosity&) {
978 if (&particle_field_out != &particle_field_in) {
979 particle_field_out = particle_field_in;
980 }
981
982 // Check that particle_field contains realistic values
983 //(values smaller than the threshold will be set to 0)
984 for (Index i = 0; i < particle_field_out.nbooks(); i++) {
985 for (Index j = 0; j < particle_field_out.npages(); j++) {
986 for (Index k = 0; k < particle_field_out.nrows(); k++) {
987 for (Index l = 0; l < particle_field_out.ncols(); l++) {
988 if (particle_field_out(i, j, k, l) < threshold) {
989 particle_field_out(i, j, k, l) = 0.0;
990 }
991 }
992 }
993 }
994 }
995}
996
997/* Workspace method: Doxygen documentation will be auto-generated */
998void ScatSpeciesInit( //WS Output:
999 ArrayOfString& scat_species,
1002 Index& scat_data_checked,
1003 ArrayOfGriddedField3& pnd_field_raw,
1004 const Verbosity&) {
1005 scat_species.resize(0);
1006 scat_data_raw.resize(0);
1007 scat_meta.resize(0);
1008 pnd_field_raw.resize(0);
1009 scat_data_checked = 0;
1010}
1011
1012/* Workspace method: Doxygen documentation will be auto-generated */
1015 ArrayOfGriddedField3& pnd_field_raw,
1016 // WS Input (needed for checking the datafiles):
1017 const Index& atmosphere_dim,
1018 // Keywords:
1019 const ArrayOfString& scat_data_files,
1020 const ArrayOfString& pnd_field_files,
1021 const Verbosity& verbosity) {
1023
1024 //--- Check input ---------------------------------------------------------
1025
1026 // Atmosphere
1027 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1028 //chk_atm_grids( atmosphere_dim, p_grid, lat_grid, lon_grid );
1029
1030 //--- Reading the data ---------------------------------------------------
1031
1032 ARTS_USER_ERROR_IF (scat_data_files.nelem() != pnd_field_files.nelem(),
1033 "Number of elements in scat_data and pnd_field filelists is"
1034 "inconsistent.")
1035
1036 Index last_species = scat_data_raw.nelem() - 1;
1037 if (last_species == -1) {
1038 scat_data_raw.resize(1);
1039 last_species = 0;
1040 }
1041
1042 // Create empty dummy elements to append to *scat_data_raw* and *pnd_field_raw*.
1043 SingleScatteringData scat_data_single;
1044 GriddedField3 pnd_field_data;
1045
1046 for (Index i = 0; i < scat_data_files.nelem(); i++) {
1047 // Append *scat_data_raw* and *pnd_field_raw* with empty Arrays of Tensors.
1048 scat_data_raw[last_species].push_back(scat_data_single);
1049 pnd_field_raw.push_back(pnd_field_data);
1050
1051 out2 << " Read single scattering data file " << scat_data_files[i] << "\n";
1053 scat_data_files[i],
1054 scat_data_raw[last_species][scat_data_raw[last_species].nelem() - 1],
1055 verbosity);
1056
1057 out2 << " Read particle number density field\n";
1058 if (pnd_field_files[i].nelem() < 1) {
1060 out1 << "Warning: No pnd_field_file specified. Ignored here,\n"
1061 << "but user HAS TO add that later on!. \n";
1062 } else {
1063 xml_read_from_file(pnd_field_files[i],
1064 pnd_field_raw[pnd_field_raw.nelem() - 1],
1065 verbosity);
1066
1067 chk_pnd_data(pnd_field_raw[pnd_field_raw.nelem() - 1],
1068 pnd_field_files[i],
1069 atmosphere_dim,
1070 verbosity);
1071 }
1072 }
1073}
1074
1075/* Workspace method: Doxygen documentation will be auto-generated */
1078 ArrayOfGriddedField3& pnd_field_raw,
1079 // WS Input(needed for checking the datafiles):
1080 const Index& atmosphere_dim,
1081 // Keywords:
1082 const ArrayOfString& scat_data_files,
1083 const String& pnd_fieldarray_file,
1084 const Verbosity& verbosity) {
1086
1087 //--- Check input ---------------------------------------------------------
1088
1089 // Atmosphere
1090 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1091 //chk_atm_grids ( atmosphere_dim, p_grid, lat_grid, lon_grid );
1092
1093 //--- Reading the data ---------------------------------------------------
1095 arr_ssd.resize(scat_data_files.nelem());
1096
1097 for (Index i = 0; i < scat_data_files.nelem(); i++) {
1098 out2 << " Read single scattering data file " << scat_data_files[i] << "\n";
1099 xml_read_from_file(scat_data_files[i], arr_ssd[i], verbosity);
1100 }
1101
1102 // append as new scattering species
1103 if (scat_data_raw.nelem() == 0) {
1104 scat_data_raw.resize(1);
1105 scat_data_raw[0] = arr_ssd;
1106 } else
1107 scat_data_raw.push_back(arr_ssd);
1108
1109 out2 << " Read particle number density data \n";
1110 ArrayOfGriddedField3 pnd_tmp;
1111 xml_read_from_file(pnd_fieldarray_file, pnd_tmp, verbosity);
1112
1113 chk_pnd_raw_data(pnd_tmp, pnd_fieldarray_file, atmosphere_dim, verbosity);
1114
1115 // append to pnd_field_raw
1116 for (Index i = 0; i < pnd_tmp.nelem(); ++i)
1117 pnd_field_raw.push_back(pnd_tmp[i]);
1118}
1119
1120/* Workspace method: Doxygen documentation will be auto-generated */
1123 ArrayOfGriddedField3& vmr_field_raw,
1124 ArrayOfArrayOfSpeciesTag& abs_species,
1125 Index& propmat_clearsky_agenda_checked,
1126 Index& abs_xsec_agenda_checked,
1127 // WS Input (needed for checking the datafiles):
1128 const Index& atmosphere_dim,
1129 const Vector& f_grid,
1130 // Keywords:
1131 const ArrayOfString& scat_data_files,
1132 const ArrayOfString& pnd_field_files,
1133 const Verbosity& verbosity) {
1135
1136 //--- Check input ---------------------------------------------------------
1137
1138 // Atmosphere
1139 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1140 //chk_atm_grids( atmosphere_dim, p_grid, lat_grid, lon_grid );
1141
1142 // Frequency grid
1143 ARTS_USER_ERROR_IF (f_grid.empty(), "The frequency grid is empty.");
1144 chk_if_increasing("f_grid", f_grid);
1145
1146 //--- Reading the data ---------------------------------------------------
1147
1148 ARTS_USER_ERROR_IF (scat_data_files.nelem() != pnd_field_files.nelem(),
1149 "Number of elements in scat_data and pnd_field filelists is"
1150 "inconsistent.")
1151
1152 Index last_species = scat_data_raw.nelem() - 1;
1153 if (last_species == -1) {
1154 scat_data_raw.resize(1);
1155 last_species = 0;
1156 }
1157
1158 // Create empty dummy elements to append to *scat_data_raw* and *pnd_field_raw*.
1159 SingleScatteringData scat_data_single;
1160 GriddedField3 pnd_field_data;
1161 ArrayOfString species(1);
1162 species[0] = "particles";
1163
1164 for (Index i = 0; i < scat_data_files.nelem(); i++) {
1165 // Append *scat_data_raw* and *pnd_field_raw* with empty Arrays of Tensors.
1166 scat_data_raw[last_species].push_back(scat_data_single);
1167 vmr_field_raw.push_back(pnd_field_data);
1168
1169 out2 << " Read single scattering data file " << scat_data_files[i] << "\n";
1171 scat_data_files[i],
1172 scat_data_raw[last_species][scat_data_raw[last_species].nelem() - 1],
1173 verbosity);
1174
1175 out2 << " Check single scattering properties\n";
1177 "scat_data_single.f_grid to f_grid",
1178 scat_data_raw[last_species][scat_data_raw[last_species].nelem() - 1]
1179 .f_grid,
1180 f_grid);
1181
1182 out2 << " Read particle number density field\n";
1183 if (pnd_field_files[i].nelem() < 1) {
1185 out1 << "Warning: No pnd_field_file specified. Ignored here,\n"
1186 << "but user HAS TO add that later on!\n";
1187 } else {
1188 try {
1189 xml_read_from_file(pnd_field_files[i],
1190 vmr_field_raw[vmr_field_raw.nelem() - 1],
1191 verbosity);
1192 } catch (...) {
1194 try {
1195 xml_read_from_file(pnd_field_files[i], tmp, verbosity);
1196 if (tmp.nelem() == 1) {
1197 vmr_field_raw[vmr_field_raw.nelem() - 1] = tmp[0];
1198 } else {
1200 "The file ", pnd_field_files[i], "\n"
1201 "is neither GriddedField3 nor a 1-long ArrayOfGriddedField3.\n")
1202 }
1203 } catch (...) {
1205 "The file ", pnd_field_files[i], " does not exist or\n"
1206 "its type is neither GriddedField3 nor a 1-long ArrayOfGriddedField3.\n")
1207 }
1208 }
1209
1210 chk_pnd_data(vmr_field_raw[vmr_field_raw.nelem() - 1],
1211 pnd_field_files[i],
1212 atmosphere_dim,
1213 verbosity);
1214 }
1215
1216 out2 << " Append 'particle' field to abs_species\n";
1217 abs_speciesAdd(abs_species,
1218 propmat_clearsky_agenda_checked,
1219 abs_xsec_agenda_checked,
1220 species,
1221 verbosity);
1222 }
1223 scat_dataCheck(scat_data_raw, "sane", 1e-2, verbosity);
1224}
1225
1226/* Workspace method: Doxygen documentation will be auto-generated */
1230 // Keywords:
1231 const ArrayOfString& scat_data_files,
1232 const Verbosity& verbosity) {
1235
1236 //--- Reading the data ---------------------------------------------------
1239
1240 arr_ssd.resize(scat_data_files.nelem());
1241 arr_smd.resize(scat_data_files.nelem());
1242
1243 Index meta_naming_conv = 0;
1244
1245 for (Index i = 0; i < 1 && i < scat_data_files.nelem(); i++) {
1246 out3 << " Read single scattering data file " << scat_data_files[i] << "\n";
1247 xml_read_from_file(scat_data_files[i], arr_ssd[i], verbosity);
1248
1249 // make meta data name from scat data name
1250 ArrayOfString strarr;
1251 String scat_meta_file;
1252
1253 if (i == 0) {
1254 scat_data_files[i].split(strarr, ".xml");
1255 scat_meta_file = strarr[0] + ".meta.xml";
1256
1257 try {
1258 find_xml_file(scat_meta_file, verbosity);
1259 } catch (const runtime_error&) {
1260 }
1261
1262 if (file_exists(scat_meta_file)) {
1263 out3 << " Read scattering meta data\n";
1264
1265 xml_read_from_file(scat_meta_file, arr_smd[i], verbosity);
1266
1267 meta_naming_conv = 1;
1268 } else {
1269 try {
1270 scat_data_files[i].split(strarr, "scat_data");
1271 ARTS_USER_ERROR_IF (strarr.nelem() < 2,
1272 "Splitting scattering data filename up at 'scat_data' also failed.");
1273 scat_meta_file = strarr[0] + "scat_meta" + strarr[1];
1274
1275 out3 << " Read scattering meta data\n";
1276 xml_read_from_file(scat_meta_file, arr_smd[i], verbosity);
1277
1278 meta_naming_conv = 2;
1279 } catch (const std::runtime_error& e) {
1281 "No meta data file following one of the allowed naming "
1282 "conventions was found.\n"
1283 "Allowed are "
1284 "*.meta.xml from *.xml and "
1285 "*scat_meta* from *scat_data*\n"
1286 "Scattering meta data file not found: ", scat_meta_file,
1287 "\n",
1288 e.what())
1289 }
1290 }
1291 }
1292 }
1293
1294 ArrayOfString fail_msg;
1295
1296#pragma omp parallel for if (!arts_omp_in_parallel() && \
1297 scat_data_files.nelem() > 1) \
1298 num_threads(arts_omp_get_max_threads() > 16 ? 16 \
1299 : arts_omp_get_max_threads()) \
1300 shared(out3, scat_data_files, arr_ssd, arr_smd)
1301 for (Index i = 1; i < scat_data_files.nelem(); i++) {
1302 // make meta data name from scat data name
1303 ArrayOfString strarr;
1304 String scat_meta_file;
1307
1308 try {
1309 out3 << " Read single scattering data file " << scat_data_files[i]
1310 << "\n";
1311 xml_read_from_file(scat_data_files[i], ssd, verbosity);
1312
1313 scat_data_files[i].split(strarr, ".xml");
1314 scat_meta_file = strarr[0] + ".meta.xml";
1315
1316 if (meta_naming_conv == 1) {
1317 scat_data_files[i].split(strarr, ".xml");
1318 scat_meta_file = strarr[0] + ".meta.xml";
1319
1320 out3 << " Read scattering meta data\n";
1321 xml_read_from_file(scat_meta_file, smd, verbosity);
1322 } else {
1323 scat_data_files[i].split(strarr, "scat_data");
1324 scat_meta_file = strarr[0] + "scat_meta" + strarr[1];
1325
1326 out3 << " Read scattering meta data\n";
1327 xml_read_from_file(scat_meta_file, smd, verbosity);
1328 }
1329 } catch (const std::exception& e) {
1330 ostringstream os;
1331 os << "Run-time error reading scattering data : \n" << e.what();
1332#pragma omp critical(ybatchCalc_push_fail_msg)
1333 fail_msg.push_back(os.str());
1334 }
1335
1336 //FIXME: currently nothing is done in chk_scattering_meta_data!
1337 chk_scattering_meta_data(smd, scat_meta_file, verbosity);
1338
1339#pragma omp critical(ScatSpeciesScatAndMetaRead_assign_ssd)
1340 arr_ssd[i] = std::move(ssd);
1341#pragma omp critical(ScatSpeciesScatAndMetaRead_assign_smd)
1342 arr_smd[i] = std::move(smd);
1343 }
1344
1345 if (fail_msg.nelem()) {
1346 std::ostringstream os;
1347 for (auto& msg : fail_msg) os << msg << '\n';
1348
1349 ARTS_USER_ERROR ( os.str());
1350 }
1351
1352 // check if arrays have same size
1353 chk_scattering_data(arr_ssd, arr_smd, verbosity);
1354
1355 // append as new scattering species
1356 scat_data_raw.push_back(std::move(arr_ssd));
1357 scat_meta.push_back(std::move(arr_smd));
1358}
1359
1360/* Workspace method: Doxygen documentation will be auto-generated */
1361void ScatElementsSelect( //WS Output:
1364 // WS Input:
1365 const ArrayOfString& scat_species,
1366 const String& species,
1367 const String& sizeparam,
1368 const Numeric& sizemin,
1369 const Numeric& sizemax,
1370 const Numeric& tolerance,
1371 const String& delim,
1372 const Verbosity&) {
1373 // first check that sizes of scat_species and scat_data_raw/scat_meta agree
1374 Index nspecies = scat_species.nelem();
1375 ARTS_USER_ERROR_IF (nspecies != scat_data_raw.nelem() || nspecies != scat_meta.nelem(),
1376 "Number of scattering species specified by scat_species does\n"
1377 "not agree with number of scattering species in\n"
1378 "scat_data_raw or scat_meta:\n"
1379 "scat_species has ", nspecies,
1380 " entries, while scat_data_raw has ", scat_data_raw.nelem(),
1381 " and scat_meta has ", scat_meta.nelem(), ".")
1382
1383 // create temporary containers for selected elements
1384 ArrayOfSingleScatteringData scat_data_raw_tmp;
1385 ArrayOfScatteringMetaData scat_meta_tmp;
1386
1387 String partfield_name;
1388 //find the species to handle: compare 'species' to 'partfield' part of
1389 //scat_species tags
1390 Index i_ss = -1;
1391 for (Index i = 0; i < scat_species.nelem(); i++) {
1392 parse_partfield_name(partfield_name, scat_species[i], delim);
1393 if (partfield_name == species) i_ss = i;
1394 }
1395 ARTS_USER_ERROR_IF (i_ss < 0,
1396 "Scattering species ", species, " not found among scat_species.")
1397
1398 // choosing the specified SingleScatteringData and ScatteringMetaData
1399 if (sizeparam == "diameter_max")
1400 for (Index i_se = 0; i_se < scat_meta[i_ss].nelem(); i_se++) {
1401 // scattering element diameter is extracted from the
1402 // scattering element's meta data and checked whether it's within size
1403 // selected range (sizemax < 0 check follows from wildcard usage and
1404 // means consider all sizes on the upper end)
1405 if (scat_meta[i_ss][i_se].diameter_max > sizemin - sizemin * tolerance &&
1406 (sizemax + sizemax * tolerance > scat_meta[i_ss][i_se].diameter_max ||
1407 sizemax < 0.)) {
1408 // copy selected scattering element to temp arrays
1409 scat_data_raw_tmp.push_back(scat_data_raw[i_ss][i_se]);
1410 scat_meta_tmp.push_back(scat_meta[i_ss][i_se]);
1411 }
1412 }
1413 else if (sizeparam == "diameter_volume_equ")
1414 for (Index i_se = 0; i_se < scat_meta[i_ss].nelem(); i_se++) {
1415 if (scat_meta[i_ss][i_se].diameter_volume_equ >
1416 sizemin - sizemin * tolerance &&
1417 (sizemax + sizemax * tolerance >
1418 scat_meta[i_ss][i_se].diameter_volume_equ ||
1419 sizemax < 0.)) {
1420 // copy selected scattering element to temp arrays
1421 scat_data_raw_tmp.push_back(scat_data_raw[i_ss][i_se]);
1422 scat_meta_tmp.push_back(scat_meta[i_ss][i_se]);
1423 }
1424 }
1425 else if (sizeparam == "diameter_area_equ_aerodynamical")
1426 for (Index i_se = 0; i_se < scat_meta[i_ss].nelem(); i_se++) {
1427 if (scat_meta[i_ss][i_se].diameter_area_equ_aerodynamical >
1428 sizemin - sizemin * tolerance &&
1429 (sizemax + sizemax * tolerance >
1430 scat_meta[i_ss][i_se].diameter_area_equ_aerodynamical ||
1431 sizemax < 0.)) {
1432 // copy selected scattering element to temp arrays
1433 scat_data_raw_tmp.push_back(scat_data_raw[i_ss][i_se]);
1434 scat_meta_tmp.push_back(scat_meta[i_ss][i_se]);
1435 }
1436 }
1437 else {
1439 "Size parameter ", sizeparam, "is unknown.")
1440 }
1441
1442 // To use a particle species field without associated scattering element
1443 // data poses a high risk of accidentially neglecting these species. That's
1444 // unlikely what the user intends. Hence throw error.
1445 ARTS_USER_ERROR_IF (scat_meta_tmp.nelem() < 1,
1446 "For scattering species ", species, " no scattering "
1447 "element matching the requested size range found.\n"
1448 "Check *scat_data_raw* and *scat_meta* input as well as your size limit "
1449 "selection!")
1450
1451 scat_meta[i_ss] = std::move(scat_meta_tmp);
1452 scat_data_raw[i_ss] = std::move(scat_data_raw_tmp);
1453
1454 // check if array is empty. should never apply (since we checked the re-worked
1455 // data before and that error should also catch cases that are empty from the
1456 // beginning).
1458}
1459
1460/* Workspace method: Doxygen documentation will be auto-generated */
1463 // Keywords:
1464 const ArrayOfString& scat_species,
1465 const String& species,
1466 const String& scat_species_delim,
1467 const Numeric& T_low,
1468 const Numeric& T_high,
1469 const Verbosity&) {
1470 const bool do_tl = (T_low >= 0.);
1471 const bool do_th = (T_high >= 0.);
1472
1473 if (do_tl || do_th) {
1474 Index i_ss = -1;
1475 if (species == "") {
1476 i_ss = scat_data_raw.nelem() - 1;
1477 ARTS_USER_ERROR_IF (i_ss == -1,
1478 "No *scat_data* available. Can not extend temperature range on "
1479 "inexistent data.")
1480 } else {
1481 // first check that sizes of scat_species and scat_data_raw agree
1482 Index nspecies = scat_species.nelem();
1483 ARTS_USER_ERROR_IF (nspecies != scat_data_raw.nelem(),
1484 "Number of scattering species specified by scat_species does\n"
1485 "not agree with number of scattering species in *scat_data*:\n"
1486 "scat_species has ", nspecies,
1487 " entries, while *scat_data* has ", scat_data_raw.nelem(),
1488 ".")
1489 String partfield_name;
1490 //find the species to handle: compare 'species' to 'partfield' part of
1491 //scat_species tags
1492 for (Index i = 0; i < scat_species.nelem(); i++) {
1494 partfield_name, scat_species[i], scat_species_delim);
1495 if (partfield_name == species) i_ss = i;
1496 }
1497 ARTS_USER_ERROR_IF (i_ss < 0,
1498 "Scattering species ", species, " not found among scat_species.")
1499 }
1500
1501 for (Index i_se = 0; i_se < scat_data_raw[i_ss].nelem(); i_se++) {
1502 const SingleScatteringData& ssdo = scat_data_raw[i_ss][i_se];
1503 const Index nTo = ssdo.T_grid.nelem();
1504 Index nTn = nTo;
1505 bool do_htl, do_hth;
1506 if (nTo > 1) {
1507 do_htl = (do_tl && (T_low < ssdo.T_grid[0]));
1508 do_hth = (do_th && (T_high > last(ssdo.T_grid)));
1509 } else {
1510 do_htl = false;
1511 do_hth = false;
1512 }
1513
1514 if (do_htl || do_hth) {
1515 // create new instance of SingleScatteringData
1517 Index iToff;
1518
1519 // determine new temperature grid
1520 if (do_htl) nTn += 1;
1521 if (do_hth) nTn += 1;
1522 Vector T_grid_new(nTn);
1523 if (do_htl) {
1524 T_grid_new[0] = T_low;
1525 iToff = 1;
1526 } else {
1527 iToff = 0;
1528 }
1529 for (Index iT = 0; iT < nTo; iT++)
1530 T_grid_new[iT + iToff] = scat_data_raw[i_ss][i_se].T_grid[iT];
1531 if (do_hth) T_grid_new[nTo + iToff] = T_high;
1532 ssdn.T_grid = std::move(T_grid_new);
1533
1534 // copy grids and other descriptive data that is to remain identical
1535 ssdn.ptype = ssdo.ptype;
1536 ostringstream description;
1537 description << ssdo.description; // here just copy. we append further
1538 // info below if applicable.
1539 ssdn.f_grid = ssdo.f_grid;
1540 ssdn.za_grid = ssdo.za_grid;
1541 ssdn.aa_grid = ssdo.aa_grid;
1542
1543 // determine size of current optical property data
1544 const Index nf = ssdo.f_grid.nelem();
1545 const Index nzas = ssdo.pha_mat_data.nshelves();
1546 const Index naas = ssdo.pha_mat_data.nbooks();
1547 const Index nzai = ssdo.pha_mat_data.npages();
1548 const Index naai = ssdo.pha_mat_data.nrows();
1549 const Index nmep = ssdo.pha_mat_data.ncols();
1550 const Index nmee = ssdo.ext_mat_data.ncols();
1551 const Index nvea = ssdo.abs_vec_data.ncols();
1552
1553 // create containers for extended optical property data
1554 ssdn.pha_mat_data.resize(nf, nTn, nzas, naas, nzai, naai, nmep);
1555 ssdn.ext_mat_data.resize(nf, nTn, nzai, naai, nmee);
1556 ssdn.abs_vec_data.resize(nf, nTn, nzai, naai, nvea);
1557
1558 // copy optical property data
1559 for (Index iT = 0; iT < nTo; iT++) {
1560 ssdn.pha_mat_data(
1561 joker, iT + iToff, joker, joker, joker, joker, joker) =
1562 ssdo.pha_mat_data(joker, iT, joker, joker, joker, joker, joker);
1563 ssdn.ext_mat_data(joker, iT + iToff, joker, joker, joker) =
1564 ssdo.ext_mat_data(joker, iT, joker, joker, joker);
1565 ssdn.abs_vec_data(joker, iT + iToff, joker, joker, joker) =
1566 ssdo.abs_vec_data(joker, iT, joker, joker, joker);
1567 }
1568
1569 // duplicate optical property data on T-edges if applicable
1570 if (do_htl) {
1571 ssdn.pha_mat_data(joker, 0, joker, joker, joker, joker, joker) =
1573 ssdn.ext_mat_data(joker, 0, joker, joker, joker) =
1574 ssdn.ext_mat_data(joker, 1, joker, joker, joker);
1575 ssdn.abs_vec_data(joker, 0, joker, joker, joker) =
1576 ssdn.abs_vec_data(joker, 1, joker, joker, joker);
1577 description << "\n"
1578 << "Low temperature limit extended by"
1579 << " duplicating previous low temperature limit"
1580 << " single scattering properties.";
1581 }
1582 if (do_hth) {
1583 ssdn.pha_mat_data(joker, nTn - 1, joker, joker, joker, joker, joker) =
1584 ssdn.pha_mat_data(
1585 joker, nTn - 2, joker, joker, joker, joker, joker);
1586 ssdn.ext_mat_data(joker, nTn - 1, joker, joker, joker) =
1587 ssdn.ext_mat_data(joker, nTn - 2, joker, joker, joker);
1588 ssdn.abs_vec_data(joker, nTn - 1, joker, joker, joker) =
1589 ssdn.abs_vec_data(joker, nTn - 2, joker, joker, joker);
1590 description << "\n"
1591 << "High temperature limit extended by"
1592 << " duplicating previous high temperature limit"
1593 << " single scattering properties.";
1594 }
1595 ssdn.description = description.str();
1596 scat_data_raw[i_ss][i_se] = std::move(ssdn);
1597 }
1598 }
1599 }
1600}
1601
1602/* Workspace method: Doxygen documentation will be auto-generated */
1604 Tensor4& pnd_field,
1605 ArrayOfTensor4& dpnd_field_dx,
1606 //WS Input
1607 const Vector& p_grid,
1608 const Vector& lat_grid,
1609 const Vector& lon_grid,
1610 const ArrayOfGriddedField3& pnd_field_raw,
1611 const Index& atmosphere_dim,
1612 const ArrayOfIndex& cloudbox_limits,
1613 const ArrayOfRetrievalQuantity& jacobian_quantities,
1614 const Index& zeropadding,
1615 const Verbosity& verbosity) {
1616 // Basic checks of input variables
1617 //
1618 // Particle number density data
1619 //
1620 ARTS_USER_ERROR_IF (pnd_field_raw.empty(),
1621 "No particle number density data given. Please use WSMs\n"
1622 "*ParticleTypeInit* and *ParticleTypeAdd(All)* for reading\n"
1623 "scattering element data.\n")
1624
1625 chk_atm_grids(atmosphere_dim, p_grid, lat_grid, lon_grid);
1626 ArrayOfIndex cloudbox_limits_tmp;
1627 /*if ( cloudbox_limits.empty() )
1628 {
1629 //If no limits set, the cloud(box) is supposed to cover the
1630 //complete atmosphere. This particularly to facilitate use of
1631 //scat_data_single&pnd_fields for non-scatt, greybody cloud calculations.
1632 //Hence, set the limits to first & last elements of the different grids.
1633 //Note: no margin left at lat/lon_grid edges!.
1634 cloudbox_limits_tmp.resize(2*atmosphere_dim);
1635
1636 // Pressure limits
1637 cloudbox_limits_tmp[0] = 0;
1638 cloudbox_limits_tmp[1] = p_grid.nelem() - 1;
1639 // Latitude limits
1640 if( atmosphere_dim >= 2 )
1641 {
1642 cloudbox_limits_tmp[2] = 0;
1643 cloudbox_limits_tmp[3] = lat_grid.nelem() - 1;
1644 }
1645 // Longitude limits
1646 if( atmosphere_dim == 3 )
1647 {
1648 cloudbox_limits_tmp[4] = 0;
1649 cloudbox_limits_tmp[5] = lon_grid.nelem() - 1;
1650 }
1651 }
1652 else */
1653 ARTS_USER_ERROR_IF (cloudbox_limits.nelem() != 2 * atmosphere_dim,
1654 "*cloudbox_limits* is a vector which contains the"
1655 "upper and lower limit of the cloud for all "
1656 "atmospheric dimensions. So its dimension must"
1657 "be 2 x *atmosphere_dim*");
1658 cloudbox_limits_tmp = cloudbox_limits;
1659
1660 // Check that pnd_field_raw has at least 2 grid-points in each dimension.
1661 // Otherwise, interpolation further down will fail with assertion.
1662
1663 for (Index d = 0; d < atmosphere_dim; d++) {
1664 for (Index i = 0; i < pnd_field_raw.nelem(); i++) {
1665 ARTS_USER_ERROR_IF (pnd_field_raw[i].get_grid_size(d) < 2,
1666 "Error in pnd_field_raw data. "
1667 "Dimension ", d, " (name: \"",
1668 pnd_field_raw[i].get_grid_name(d),
1669 "\") has only ",
1670 pnd_field_raw[i].get_grid_size(d),
1671 " element",
1672 ((pnd_field_raw[i].get_grid_size(d) == 1) ? "" : "s"),
1673 ". Must be at least 2.")
1674 }
1675 }
1676 const Index Np_cloud = cloudbox_limits_tmp[1] - cloudbox_limits_tmp[0] + 1;
1677
1678 ConstVectorView p_grid_cloud =
1679 p_grid[Range(cloudbox_limits_tmp[0], Np_cloud)];
1680
1681 // Check that no scatterers exist outside the cloudbox
1683 pnd_field_raw,
1684 p_grid,
1685 lat_grid,
1686 lon_grid,
1687 cloudbox_limits_tmp);
1688
1689 //==========================================================================
1690 if (atmosphere_dim == 1) {
1691 ArrayOfGriddedField3 pnd_field_tmp;
1692
1694 pnd_field_tmp, p_grid_cloud, pnd_field_raw, 1, zeropadding, verbosity);
1695
1696 FieldFromGriddedField(pnd_field,
1697 p_grid_cloud,
1698 pnd_field_tmp[0].get_numeric_grid(1),
1699 pnd_field_tmp[0].get_numeric_grid(2),
1700 pnd_field_tmp,
1701 verbosity);
1702 } else if (atmosphere_dim == 2) {
1703 const Index Nlat_cloud =
1704 cloudbox_limits_tmp[3] - cloudbox_limits_tmp[2] + 1;
1705
1706 ConstVectorView lat_grid_cloud =
1707 lat_grid[Range(cloudbox_limits_tmp[2], Nlat_cloud)];
1708
1709 if (zeropadding) {
1710 // FIXME: OLE: Implement this
1712 out0 << "WARNING: zeropadding currently only supported for 1D.";
1713 }
1714
1715 //Resize variables
1716 pnd_field.resize(pnd_field_raw.nelem(), Np_cloud, Nlat_cloud, 1);
1717
1718 // Gridpositions:
1719 ArrayOfGridPos gp_p(Np_cloud);
1720 ArrayOfGridPos gp_lat(Nlat_cloud);
1721
1722 // Interpolate pnd_field.
1723 // Loop over the scattering element:
1724 for (Index i = 0; i < pnd_field_raw.nelem(); ++i) {
1725 // Calculate grid positions:
1726 p2gridpos(gp_p,
1727 pnd_field_raw[i].get_numeric_grid(GFIELD3_P_GRID),
1728 p_grid_cloud);
1729 gridpos(gp_lat,
1730 pnd_field_raw[i].get_numeric_grid(GFIELD3_LAT_GRID),
1731 lat_grid_cloud);
1732
1733 // Interpolation weights:
1734 Tensor3 itw(Np_cloud, Nlat_cloud, 4);
1735 interpweights(itw, gp_p, gp_lat);
1736
1737 // Interpolate:
1738 interp(pnd_field(i, joker, joker, 0),
1739 itw,
1740 pnd_field_raw[i].data(joker, joker, 0),
1741 gp_p,
1742 gp_lat);
1743 }
1744 } else {
1745 const Index Nlat_cloud =
1746 cloudbox_limits_tmp[3] - cloudbox_limits_tmp[2] + 1;
1747 const Index Nlon_cloud =
1748 cloudbox_limits_tmp[5] - cloudbox_limits_tmp[4] + 1;
1749
1750 if (zeropadding) {
1751 // FIXME: OLE: Implement this
1753 out0 << "WARNING: zeropadding currently only supported for 1D.";
1754 }
1755
1756 ConstVectorView lat_grid_cloud =
1757 lat_grid[Range(cloudbox_limits_tmp[2], Nlat_cloud)];
1758 ConstVectorView lon_grid_cloud =
1759 lon_grid[Range(cloudbox_limits_tmp[4], Nlon_cloud)];
1760
1761 //Resize variables
1762 pnd_field.resize(pnd_field_raw.nelem(), Np_cloud, Nlat_cloud, Nlon_cloud);
1763
1764 // Gridpositions:
1765 ArrayOfGridPos gp_p(Np_cloud);
1766 ArrayOfGridPos gp_lat(Nlat_cloud);
1767 ArrayOfGridPos gp_lon(Nlon_cloud);
1768
1769 // Interpolate pnd_field.
1770 // Loop over the scattering element types:
1771 for (Index i = 0; i < pnd_field_raw.nelem(); ++i) {
1772 // Calculate grid positions:
1773 p2gridpos(gp_p,
1774 pnd_field_raw[i].get_numeric_grid(GFIELD3_P_GRID),
1775 p_grid_cloud);
1776 gridpos(gp_lat,
1777 pnd_field_raw[i].get_numeric_grid(GFIELD3_LAT_GRID),
1778 lat_grid_cloud);
1779 gridpos(gp_lon,
1780 pnd_field_raw[i].get_numeric_grid(GFIELD3_LON_GRID),
1781 lon_grid_cloud);
1782
1783 // Interpolation weights:
1784 Tensor4 itw(Np_cloud, Nlat_cloud, Nlon_cloud, 8);
1785 interpweights(itw, gp_p, gp_lat, gp_lon);
1786
1787 // Interpolate:
1788 interp(pnd_field(i, joker, joker, joker),
1789 itw,
1790 pnd_field_raw[i].data,
1791 gp_p,
1792 gp_lat,
1793 gp_lon);
1794 }
1795 }
1796
1797 // no (cloudy) Jacobians with this WSM, hence no calc.
1798 // but we need to size dpnd_field to be consistent with jacobian_quantities.
1799 dpnd_field_dx.resize(jacobian_quantities.nelem());
1800}
1801
1802/* Workspace method: Doxygen documentation will be auto-generated */
1804 const Index& atmosphere_dim,
1805 const Index& cloudbox_on,
1806 const ArrayOfIndex& cloudbox_limits,
1807 const Index& nzero,
1808 const Verbosity&) {
1809 /* if( !cloudbox_checked )
1810 throw runtime_error( "The cloudbox must be flagged to have passed a "
1811 "consistency check (cloudbox_checked=1)." );
1812*/
1813 ARTS_USER_ERROR_IF (atmosphere_dim == 1,
1814 "No use in calling this method for 1D.");
1815 ARTS_USER_ERROR_IF (!cloudbox_on,
1816 "No use in calling this method with cloudbox off.");
1817 ARTS_USER_ERROR_IF (nzero < 1,
1818 "The argument *nzero* must be > 0.");
1819
1820 // Sizes
1821 const Index npart = pnd_field.nbooks();
1822 const Index np = cloudbox_limits[1] - cloudbox_limits[0] + 1;
1823 const Index nlat = cloudbox_limits[3] - cloudbox_limits[2] + 1;
1824 Index nlon = 1;
1825 if (atmosphere_dim == 3) {
1826 nlon = cloudbox_limits[5] - cloudbox_limits[4] + 1;
1827 }
1828
1829 ARTS_USER_ERROR_IF (pnd_field.npages() != np ||
1830 pnd_field.nrows() != 1 ||
1831 pnd_field.ncols() != 1,
1832 "The input *pnd_field* is either not 1D or does not "
1833 "match pressure size of cloudbox.");
1834
1835 // Temporary container
1836 Tensor4 pnd_temp = pnd_field;
1837
1838 // Resize and fill
1839 pnd_field.resize(npart, np, nlat, nlon);
1840 pnd_field = 0;
1841 //
1842 for (Index ilon = nzero; ilon < nlon - nzero; ilon++) {
1843 for (Index ilat = nzero; ilat < nlat - nzero; ilat++) {
1844 for (Index ip = 0; ip < np; ip++) {
1845 for (Index is = 0; is < npart; is++) {
1846 pnd_field(is, ip, ilat, ilon) = pnd_temp(is, ip, 0, 0);
1847 }
1848 }
1849 }
1850 }
1851}
1852
1853/* Workspace method: Doxygen documentation will be auto-generated */
1854void pnd_fieldZero( //WS Output:
1855 Tensor4& pnd_field,
1856 ArrayOfTensor4& dpnd_field_dx,
1858 //WS Input:
1859 const Index& atmosphere_dim,
1860 const Vector& f_grid,
1861 const ArrayOfIndex& cloudbox_limits,
1862 const ArrayOfRetrievalQuantity& jacobian_quantities,
1863 const Verbosity&) {
1864 chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1865
1866 ARTS_USER_ERROR_IF (cloudbox_limits.nelem() != 2 * atmosphere_dim,
1867 "*cloudbox_limits* is a vector which contains the"
1868 "upper and lower limit of the cloud for all "
1869 "atmospheric dimensions. So its dimension must"
1870 "be 2 x *atmosphere_dim*");
1871
1872 // Resize pnd_field and set it to 0:
1873 Index np = cloudbox_limits[1] - cloudbox_limits[0] + 1;
1874 Index nlat = 1, nlon = 1;
1875 if (atmosphere_dim > 1) {
1876 nlat = cloudbox_limits[3] - cloudbox_limits[2] + 1;
1877 if (atmosphere_dim > 2) {
1878 nlon = cloudbox_limits[5] - cloudbox_limits[4] + 1;
1879 }
1880 }
1881
1882 // no (cloudy) Jacobians with this WSM, hence no setting.
1883 // but we need to size dpnd_field to be consistent with jacobian_quantities.
1884 dpnd_field_dx.resize(jacobian_quantities.nelem());
1885
1886 // Do only reset scat_data if it has not been set yet.
1887 // There's no need otherwise, and it's rather unpractical for testing when
1888 // doing so: we might want to do some actual calcs with the scat_data later
1889 // on. So why throw it away?
1890 const Index N_se = TotalNumberOfElements(scat_data);
1891 if (N_se > 0) {
1892 pnd_field.resize(N_se, np, nlat, nlon);
1893 } else {
1894 pnd_field.resize(1, np, nlat, nlon);
1895
1896 //Resize scat_data and set it to 0:
1897 // Number of scattering elements
1898 scat_data.resize(1);
1899 scat_data[0].resize(1);
1900 scat_data[0][0].ptype = PTYPE_TOTAL_RND;
1901 scat_data[0][0].description = " ";
1902 // Grids which contain full ranges which one wants to calculate
1903 Index nf = f_grid.nelem();
1904 scat_data[0][0].f_grid.resize(nf);
1905 scat_data[0][0].f_grid = f_grid;
1906 Index nT = 1;
1907 scat_data[0][0].T_grid.resize(nT);
1908 scat_data[0][0].T_grid = 270.;
1909 Index nza = 5;
1910 nlinspace(scat_data[0][0].za_grid, 0, 180, nza);
1911 // Resize the data arrays
1912 scat_data[0][0].pha_mat_data.resize(nf, nT, nza, 1, 1, 1, 6);
1913 scat_data[0][0].pha_mat_data = 0.;
1914 scat_data[0][0].ext_mat_data.resize(nf, nT, 1, 1, 1);
1915 scat_data[0][0].ext_mat_data = 0.;
1916 scat_data[0][0].abs_vec_data.resize(nf, nT, 1, 1, 1);
1917 scat_data[0][0].abs_vec_data = 0.;
1918 }
1919
1920 pnd_field = 0.;
1921}
This file contains the definition of Array.
Index TotalNumberOfElements(const Array< Array< base > > &aa)
Determine total number of elements in an ArrayOfArray.
Definition: array.h:343
The global header file for ARTS.
void * data
void GriddedFieldPRegrid(GriddedField3 &out, const Vector &p_grid, const GriddedField3 &in, const Index &interp_order, const Index &zeropadding, const Verbosity &verbosity)
WORKSPACE METHOD: GriddedFieldPRegrid.
void FieldFromGriddedField(Matrix &out, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const GriddedField2 &in, const Verbosity &verbosity)
WORKSPACE METHOD: FieldFromGriddedField.
void scat_dataCheck(const ArrayOfArrayOfSingleScatteringData &scat_data, const String &check_type, const Numeric &sca_mat_threshold, const Verbosity &verbosity)
WORKSPACE METHOD: scat_dataCheck.
void abs_speciesAdd(ArrayOfArrayOfSpeciesTag &abs_species, Index &propmat_clearsky_agenda_checked, Index &abs_xsec_agenda_checked, const ArrayOfString &species, const Verbosity &verbosity)
WORKSPACE METHOD: abs_speciesAdd.
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_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_increasing(const String &x_name, const ArrayOfIndex &x)
chk_if_increasing
Definition: check_input.cc:111
The Agenda class.
Definition: agenda_class.h:44
void set_name(const String &nname)
Set agenda name.
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 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
bool empty() const
Check if variable is empty.
Definition: matpackIV.cc:50
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: matpackVII.cc:48
Index nrows() const
Returns the number of rows.
Definition: matpackVII.cc:57
Index nvitrines() const
Returns the number of vitrines.
Definition: matpackVII.cc:45
Index npages() const
Returns the number of pages.
Definition: matpackVII.cc:54
Index ncols() const
Returns the number of columns.
Definition: matpackVII.cc:60
Index nlibraries() const
Returns the number of libraries.
Definition: matpackVII.cc:42
Index nbooks() const
Returns the number of books.
Definition: matpackVII.cc:51
A constant view of a Vector.
Definition: matpackI.h:489
bool empty() const ARTS_NOEXCEPT
Returns true if variable size is zero.
Definition: matpackI.cc:46
Index nelem() const ARTS_NOEXCEPT
Returns the number of elements.
Definition: matpackI.cc:48
The Matrix class.
Definition: matpackI.h:1225
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1056
The range class.
Definition: matpackI.h:165
The Tensor3 class.
Definition: matpackIII.h:339
The Tensor4 class.
Definition: matpackIV.h:421
void resize(Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackIV.cc:1058
void resize(Index s, Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackV.cc:1741
The Tensor7 class.
Definition: matpackVII.h:2382
void resize(Index l, Index v, Index s, Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackVII.cc:5482
The Vector class.
Definition: matpackI.h:876
void split(Array< my_basic_string< charT > > &aos, const my_basic_string< charT > &delim) const
Split string into substrings.
Definition: mystring.h:213
void chk_scat_species_field(bool &empty_flag, const Tensor3 &scat_species_field, const String &fieldname, const Index &dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid)
Check whether field of a specific scattering species zero everywhere.
Definition: cloudbox.cc:607
void chk_scattering_meta_data(const ScatteringMetaData &scat_meta_single, const String &scat_meta_file, const Verbosity &verbosity)
Check scattering data meta.
Definition: cloudbox.cc:260
void find_cloudlimits(Index &lower, Index &upper, const Tensor3 &scat_species_field, const Index &atmosphere_dim, const Numeric &cloudbox_margin)
Adjust uppermost and lowermost cloudy level for one scat_species_*_*_field.
Definition: cloudbox.cc:666
void chk_pnd_field_raw_only_in_cloudbox(const Index &dim, const ArrayOfGriddedField3 &pnd_field_raw, ConstVectorView p_grid, ConstVectorView lat_grid, ConstVectorView lon_grid, const ArrayOfIndex &cloudbox_limits)
chk_pnd_field_raw_only_in_cloudbox
Definition: cloudbox.cc:140
void parse_partfield_name(String &partfield_name, const String &part_string, const String &delim)
Definition: cloudbox.cc:891
void chk_pnd_raw_data(const ArrayOfGriddedField3 &pnd_field_raw, const String &pnd_field_file, const Index &atmosphere_dim, const Verbosity &verbosity)
Check particle number density files (pnd_field_raw)
Definition: cloudbox.cc:110
void chk_scattering_data(const ArrayOfSingleScatteringData &scat_data, const ArrayOfScatteringMetaData &scat_meta, const Verbosity &)
Check scattering data general.
Definition: cloudbox.cc:241
void chk_pnd_data(const GriddedField3 &pnd_field_raw, const String &pnd_field_file, const Index &atmosphere_dim, const Verbosity &verbosity)
Check particle number density files.
Definition: cloudbox.cc:67
Internal cloudbox functions.
#define DEBUG_ONLY(...)
Definition: debug.h:52
#define ARTS_ASSERT(condition,...)
Definition: debug.h:83
#define ARTS_USER_ERROR(...)
Definition: debug.h:150
#define ARTS_USER_ERROR_IF(condition,...)
Definition: debug.h:134
void find_xml_file(String &filename, const Verbosity &verbosity)
Find an xml file.
Definition: file.cc:414
bool file_exists(const String &filename)
Checks if the given file exists.
Definition: file.cc:303
This file contains basic functions to handle ASCII files.
Implementation of gridded fields.
void gridpos(ArrayOfGridPos &gp, ConstVectorView old_grid, ConstVectorView new_grid, const Numeric &extpolfac)
Set up a grid position Array.
void gridpos_upperend_check(GridPos &gp, const Index &ie)
gridpos_upperend_check
void interpweights(VectorView itw, const GridPos &tc)
Red 1D interpolation weights.
bool is_gridpos_at_index_i(const GridPos &gp, const Index &i, const bool &strict)
is_gridpos_at_index_i
Numeric interp(ConstVectorView itw, ConstVectorView a, const GridPos &tc)
Red 1D Interpolate.
Numeric fractional_gp(const GridPos &gp)
fractional_gp
Header file for interpolation.cc.
Interpolation::Lagrange LagrangeInterpolation
#define abs(x)
#define min(a, b)
#define max(a, b)
Linear algebra functions.
bool is_size(ConstVectorView x, const Index &n)
Verifies that the size of x is l.
Definition: logic.cc:79
bool is_same_within_epsilon(const Numeric &a, const Numeric &b, const Numeric &epsilon)
Check, if two numbers agree within a given epsilon.
Definition: logic.cc:351
Header file for logic.cc.
void cloudbox_fieldCrop(Tensor7 &cloudbox_field, ArrayOfIndex &cloudbox_limits, const Index &atmosphere_dim, const Index &cloudbox_on, const Index &new_limit0, const Index &new_limit1, const Index &new_limit2, const Index &new_limit3, const Index &new_limit4, const Index &new_limit5, const Verbosity &)
WORKSPACE METHOD: cloudbox_fieldCrop.
Definition: m_cloudbox.cc:898
void cloudboxSetManually(Index &cloudbox_on, ArrayOfIndex &cloudbox_limits, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Numeric &p1, const Numeric &p2, const Numeric &lat1, const Numeric &lat2, const Numeric &lon1, const Numeric &lon2, const Verbosity &)
WORKSPACE METHOD: cloudboxSetManually.
Definition: m_cloudbox.cc:332
void ScatSpeciesExtendTemperature(ArrayOfArrayOfSingleScatteringData &scat_data_raw, const ArrayOfString &scat_species, const String &species, const String &scat_species_delim, const Numeric &T_low, const Numeric &T_high, const Verbosity &)
WORKSPACE METHOD: ScatSpeciesExtendTemperature.
Definition: m_cloudbox.cc:1461
void cloudboxSetAutomatically(Index &cloudbox_on, ArrayOfIndex &cloudbox_limits, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor4 &particle_field, const Numeric &cloudbox_margin, const Verbosity &verbosity)
WORKSPACE METHOD: cloudboxSetAutomatically.
Definition: m_cloudbox.cc:113
const Numeric LAT_LON_MIN
const Index GFIELD3_LAT_GRID
void ScatSpeciesPndAndScatAdd(ArrayOfArrayOfSingleScatteringData &scat_data_raw, ArrayOfGriddedField3 &pnd_field_raw, const Index &atmosphere_dim, const ArrayOfString &scat_data_files, const String &pnd_fieldarray_file, const Verbosity &verbosity)
WORKSPACE METHOD: ScatSpeciesPndAndScatAdd.
Definition: m_cloudbox.cc:1076
const Index GFIELD3_LON_GRID
void iyInterpCloudboxField(Matrix &iy, const Tensor7 &cloudbox_field, const Vector &rte_pos, const Vector &rte_los, const Index &jacobian_do, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const Matrix &z_surface, const Index &stokes_dim, const Vector &za_grid, const Vector &aa_grid, const Vector &f_grid, const Index &za_interp_order, const Index &za_restrict, const Index &cos_za_interp, const Numeric &za_extpolfac, const Index &aa_interp_order, const Verbosity &)
WORKSPACE METHOD: iyInterpCloudboxField.
Definition: m_cloudbox.cc:531
void ScatElementsSelect(ArrayOfArrayOfSingleScatteringData &scat_data_raw, ArrayOfArrayOfScatteringMetaData &scat_meta, const ArrayOfString &scat_species, const String &species, const String &sizeparam, const Numeric &sizemin, const Numeric &sizemax, const Numeric &tolerance, const String &delim, const Verbosity &)
WORKSPACE METHOD: ScatElementsSelect.
Definition: m_cloudbox.cc:1361
void ScatSpeciesScatAndMetaRead(ArrayOfArrayOfSingleScatteringData &scat_data_raw, ArrayOfArrayOfScatteringMetaData &scat_meta, const ArrayOfString &scat_data_files, const Verbosity &verbosity)
WORKSPACE METHOD: ScatSpeciesScatAndMetaRead.
Definition: m_cloudbox.cc:1227
void ScatElementsPndAndScatAdd(ArrayOfArrayOfSingleScatteringData &scat_data_raw, ArrayOfGriddedField3 &pnd_field_raw, const Index &atmosphere_dim, const ArrayOfString &scat_data_files, const ArrayOfString &pnd_field_files, const Verbosity &verbosity)
WORKSPACE METHOD: ScatElementsPndAndScatAdd.
Definition: m_cloudbox.cc:1013
void cloudboxSetFullAtm(Index &cloudbox_on, ArrayOfIndex &cloudbox_limits, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Verbosity &)
WORKSPACE METHOD: cloudboxSetFullAtm.
Definition: m_cloudbox.cc:259
void ScatElementsToabs_speciesAdd(ArrayOfArrayOfSingleScatteringData &scat_data_raw, ArrayOfGriddedField3 &vmr_field_raw, ArrayOfArrayOfSpeciesTag &abs_species, Index &propmat_clearsky_agenda_checked, Index &abs_xsec_agenda_checked, const Index &atmosphere_dim, const Vector &f_grid, const ArrayOfString &scat_data_files, const ArrayOfString &pnd_field_files, const Verbosity &verbosity)
WORKSPACE METHOD: ScatElementsToabs_speciesAdd.
Definition: m_cloudbox.cc:1121
void particle_fieldCleanup(Tensor4 &particle_field_out, const Tensor4 &particle_field_in, const Numeric &threshold, const Verbosity &)
WORKSPACE METHOD: particle_fieldCleanup.
Definition: m_cloudbox.cc:972
const Numeric RAD2DEG
void pnd_fieldExpand1D(Tensor4 &pnd_field, const Index &atmosphere_dim, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Index &nzero, const Verbosity &)
WORKSPACE METHOD: pnd_fieldExpand1D.
Definition: m_cloudbox.cc:1803
const Numeric PI
void pnd_fieldZero(Tensor4 &pnd_field, ArrayOfTensor4 &dpnd_field_dx, ArrayOfArrayOfSingleScatteringData &scat_data, const Index &atmosphere_dim, const Vector &f_grid, const ArrayOfIndex &cloudbox_limits, const ArrayOfRetrievalQuantity &jacobian_quantities, const Verbosity &)
WORKSPACE METHOD: pnd_fieldZero.
Definition: m_cloudbox.cc:1854
const Index GFIELD3_P_GRID
void cloudboxSetManuallyAltitude(Index &cloudbox_on, ArrayOfIndex &cloudbox_limits, const Index &atmosphere_dim, const Tensor3 &z_field, const Vector &lat_grid, const Vector &lon_grid, const Numeric &z1, const Numeric &z2, const Numeric &lat1, const Numeric &lat2, const Numeric &lon1, const Numeric &lon2, const Verbosity &)
WORKSPACE METHOD: cloudboxSetManuallyAltitude.
Definition: m_cloudbox.cc:432
void pnd_fieldCalcFrompnd_field_raw(Tensor4 &pnd_field, ArrayOfTensor4 &dpnd_field_dx, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const ArrayOfGriddedField3 &pnd_field_raw, const Index &atmosphere_dim, const ArrayOfIndex &cloudbox_limits, const ArrayOfRetrievalQuantity &jacobian_quantities, const Index &zeropadding, const Verbosity &verbosity)
WORKSPACE METHOD: pnd_fieldCalcFrompnd_field_raw.
Definition: m_cloudbox.cc:1603
const Numeric DEG2RAD
const Numeric DENSITY_OF_ICE
void cloudboxOff(Index &cloudbox_on, Index &ppath_inside_cloudbox_do, ArrayOfIndex &cloudbox_limits, Agenda &iy_cloudbox_agenda, Tensor4 &pnd_field, ArrayOfTensor4 &dpnd_field_dx, ArrayOfString &scat_species, ArrayOfArrayOfSingleScatteringData &scat_data, ArrayOfArrayOfSingleScatteringData &scat_data_raw, Index &scat_data_checked, Matrix &particle_masses, const ArrayOfRetrievalQuantity &jacobian_quantities, const Verbosity &)
WORKSPACE METHOD: cloudboxOff.
Definition: m_cloudbox.cc:82
void ScatSpeciesInit(ArrayOfString &scat_species, ArrayOfArrayOfSingleScatteringData &scat_data_raw, ArrayOfArrayOfScatteringMetaData &scat_meta, Index &scat_data_checked, ArrayOfGriddedField3 &pnd_field_raw, const Verbosity &)
WORKSPACE METHOD: ScatSpeciesInit.
Definition: m_cloudbox.cc:998
void nlinspace(Vector &x, const Numeric start, const Numeric stop, const Index n)
nlinspace
Definition: math_funcs.cc:225
Numeric last(ConstVectorView x)
last
Definition: math_funcs.cc:159
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
Declarations having to do with the four output streams.
#define CREATE_OUT1
Definition: messages.h:205
#define CREATE_OUT3
Definition: messages.h:207
#define CREATE_OUT2
Definition: messages.h:206
#define CREATE_OUT0
Definition: messages.h:204
Internal functions for microphysics calculations (size distributions etc.)
Index nelem(const Lines &l)
Number of lines.
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
void check_lagrange_interpolation(const SortedVectorType &xi, const Index polyorder=1, const std::pair< Numeric, Numeric > x={std::numeric_limits< Numeric >::infinity(), -std::numeric_limits< Numeric >::infinity()}, const Numeric extrapol=0.5, const GridType type=GridType::Standard, const std::pair< Numeric, Numeric > cycle={-180, 180})
Scattering database structure and functions.
@ PTYPE_TOTAL_RND
Definition: optproperties.h:55
This file contains header information for the dealing with command line parameters.
Numeric barometric_heightformula(const Numeric &p, const Numeric &dh)
barometric_heightformula
This file contains declerations of functions of physical character.
#define d
Declaration of functions in rte.cc.
Contains sorting routines.
void p2gridpos(ArrayOfGridPos &gp, ConstVectorView old_pgrid, ConstVectorView new_pgrid, const Numeric &extpolfac)
Calculates grid positions for pressure values.
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.
Header file for special_interp.cc.
Structure to store a grid position.
Definition: interpolation.h:73
Index idx
Definition: interpolation.h:74
This file contains basic functions to handle XML data files.
void xml_read_from_file(const String &filename, T &type, const Verbosity &verbosity)
Reads data from XML file.
Definition: xml_io.h:160