ARTS 2.5.9 (git: 825fa5f2)
m_batch.cc
Go to the documentation of this file.
1/* Copyright (C) 2004-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 "gridded_fields.h"
41using namespace std;
42
43#include "arts.h"
44#include "arts_omp.h"
45#include "auto_md.h"
46#include "math_funcs.h"
47#include "physics_funcs.h"
48#include "rte.h"
49#include "xml_io.h"
50
51inline constexpr Numeric PI=Constant::pi;
55
56/*===========================================================================
57 === The functions (in alphabetical order)
58 ===========================================================================*/
59
60/* Workspace method: Doxygen documentation will be auto-generated
61
62 2008-07-21 Stefan Buehler */
64 // WS Input:
65 const Agenda& forloop_agenda,
66 // Control Parameters:
67 const Index& start,
68 const Index& stop,
69 const Index& step,
70 const Verbosity& verbosity) {
72
73 for (Index i = start; i <= stop; i += step) {
74 out1 << " Executing for loop body, index: " << i << "\n";
75 forloop_agendaExecute(ws, i, forloop_agenda);
76 }
77}
78
79/* Workspace method: Doxygen documentation will be auto-generated */
81 // WS Output:
82 ArrayOfVector& ybatch,
83 ArrayOfArrayOfVector& ybatch_aux,
84 ArrayOfMatrix& ybatch_jacobians,
85 // WS Input:
86 const Index& ybatch_start,
87 const Index& ybatch_n,
88 const Agenda& ybatch_calc_agenda,
89 // Control Parameters:
90 const Index& robust,
91 const Verbosity& verbosity) {
93
94 Index first_ybatch_index = 0;
95
96 ArrayOfString fail_msg;
97 bool do_abort = false;
98
99 // We allow a start index ybatch_start that is different from 0. We
100 // will calculate ybatch_n jobs starting at the start
101 // index. Internally, we count from zero, which is the right
102 // counting for the output array ybatch. When we call
103 // ybatch_calc_agenda, we add ybatch_start to the internal index
104 // count.
105
106 // We create a counter, so that we can generate nice output about
107 // how many jobs are already done. (All parallel threads later will
108 // increment this, so that we really get an accurate total count!)
109 Index job_counter = 0;
110
111 // Resize the output arrays:
112 ybatch.resize(ybatch_n);
113 ybatch_aux.resize(ybatch_n);
114 ybatch_jacobians.resize(ybatch_n);
115 for (Index i = 0; i < ybatch_n; i++) {
116 ybatch[i].resize(0);
117 ybatch_aux[i].resize(0);
118 ybatch_jacobians[i].resize(0, 0);
119 }
120
121 // Go through the batch:
122
123 if (ybatch_n) {
125
126#pragma omp parallel for schedule(dynamic) if (!arts_omp_in_parallel() && \
127 ybatch_n > 1) firstprivate(wss)
128 for (Index ybatch_index = first_ybatch_index; ybatch_index < ybatch_n;
129 ybatch_index++) {
130 Index l_job_counter; // Thread-local copy of job counter.
131
132 if (do_abort) continue;
133#pragma omp critical(ybatchCalc_job_counter)
134 { l_job_counter = ++job_counter; }
135
136 {
137 ostringstream os;
138 os << " Job " << l_job_counter << " of " << ybatch_n << ", Index "
139 << ybatch_start + ybatch_index << ", Thread-Id "
140 << arts_omp_get_thread_num() << "\n";
141 out2 << os.str();
142 }
143
144 try {
145 Vector y;
146 ArrayOfVector y_aux;
147 Matrix jacobian;
148
150 y,
151 y_aux,
152 jacobian,
153 ybatch_start + ybatch_index,
154 ybatch_calc_agenda);
155
156 if (y.nelem()) {
157#pragma omp critical(ybatchCalc_assign_y)
158 ybatch[ybatch_index] = y;
159#pragma omp critical(ybatchCalc_assign_y_aux)
160 ybatch_aux[ybatch_index] = y_aux;
161
162 // Dimensions of Jacobian:
163 const Index Knr = jacobian.nrows();
164 const Index Knc = jacobian.ncols();
165
166 if (Knr != 0 || Knc != 0) {
167 if (Knr != y.nelem()) {
168 ostringstream os;
169 os << "First dimension of Jacobian must have same length as the measurement *y*.\n"
170 << "Length of *y*: " << y.nelem() << "\n"
171 << "Dimensions of *jacobian*: (" << Knr << ", " << Knc
172 << ")\n";
173 // A mismatch of the Jacobian dimension is a fatal error
174 // and should result in program termination. By setting abort
175 // to true, this will result in a runtime error in the catch
176 // block even if robust == 1
177#pragma omp critical(ybatchCalc_setabort)
178 do_abort = true;
179
180 throw runtime_error(os.str());
181 }
182
183 ybatch_jacobians[ybatch_index] = jacobian;
184
185 // After creation, all individual Jacobi matrices in the array will be
186 // empty (size zero). No need for explicit initialization.
187 }
188 }
189 } catch (const std::exception& e) {
190 if (robust && !do_abort) {
191 ostringstream os;
192 os << "WARNING! Job at ybatch_index " << ybatch_start + ybatch_index
193 << " failed.\n"
194 << "y Vector in output variable ybatch will be empty for this job.\n"
195 << "The runtime error produced was:\n"
196 << e.what() << "\n";
197 out0 << os.str();
198 } else {
199 // The user wants the batch job to fail if one of the
200 // jobs goes wrong.
201#pragma omp critical(ybatchCalc_setabort)
202 do_abort = true;
203
204 ostringstream os;
205 os << " Job at ybatch_index " << ybatch_start + ybatch_index
206 << " failed. Aborting...\n";
207 out1 << os.str();
208 }
209 ostringstream os;
210 os << "Run-time error at ybatch_index " << ybatch_start + ybatch_index
211 << ": \n"
212 << e.what();
213#pragma omp critical(ybatchCalc_push_fail_msg)
214 fail_msg.push_back(os.str());
215 }
216 }
217 }
218
219 if (fail_msg.nelem()) {
220 ostringstream os;
221
222 if (!do_abort) os << "\nError messages from failed batch cases:\n";
223 for (ArrayOfString::const_iterator it = fail_msg.begin();
224 it != fail_msg.end();
225 it++)
226 os << *it << '\n';
227
228 if (do_abort)
229 throw runtime_error(os.str());
230 else
231 out0 << os.str();
232 }
233}
234
235/* Workspace method: Doxygen documentation will be auto-generated */
237 //Output
238 ArrayOfVector& ybatch,
239 //Input
240 const ArrayOfArrayOfSpeciesTag& abs_species,
241 const Agenda& met_profile_calc_agenda,
242 const Vector& f_grid,
243 const Matrix& met_amsu_data,
244 const Matrix& sensor_pos,
245 const Vector& refellipsoid,
246 const Vector& lat_grid,
247 const Vector& lon_grid,
248 const Index& atmosphere_dim,
249 const ArrayOfArrayOfSingleScatteringData& scat_data,
250 //Keyword
251 const Index& nelem_p_grid,
252 const String& met_profile_path,
253 const String& met_profile_pnd_path,
254 const Verbosity& verbosity) {
255 GriddedField3 t_field_raw;
256 GriddedField3 z_field_raw;
257 ArrayOfGriddedField3 vmr_field_raw;
258 ArrayOfGriddedField3 pnd_field_raw;
259 Vector p_grid;
260 Matrix sensor_los;
261 Index cloudbox_on;
262 ArrayOfIndex cloudbox_limits;
263 Matrix z_surface;
264 Vector y;
265 Index no_profiles = met_amsu_data.nrows();
266
267 // *vmr_field_raw* is an ArrayOfArrayOfTensor3 where the first array
268 //holds the gaseous species.
269 //Resize *vmr_field_raw* according to the number of gaseous species
270 //elements
271 vmr_field_raw.resize(abs_species.nelem());
272
273 //pnd_field_raw is an ArrayOfArrayOfTensor3 where the first array
274 //holds the scattering elements.
275 // Number of scattering elements:
276 const Index N_se = TotalNumberOfElements(scat_data);
277
278 pnd_field_raw.resize(N_se);
279
280 // The satellite zenith angle is read in from the amsu data
281 // and converted to arts sensor_los
282 ConstVectorView sat_za_from_data = met_amsu_data(Range(joker), 3);
283
284 sensor_los.resize(1, 1);
285
286 // The lat and lon are extracted to get the proper file names of
287 // profiles
288 ConstVectorView lat = met_amsu_data(Range(joker), 0);
289 ConstVectorView lon = met_amsu_data(Range(joker), 1);
290
291 z_surface.resize(1, 1);
292
293 // The spectra .
294 y.resize(f_grid.nelem());
295
296 // The batch spectra.
297 ybatch.resize(no_profiles);
298
299 // Loop over the number of profiles.
300 for (Index i = 0; i < no_profiles; ++i) {
301 ostringstream lat_os, lon_os;
302
303 Index lat_prec = 3;
304 if (lat[i] < 0) lat_prec--;
305 if (abs(lat[i]) >= 10) {
306 lat_prec--;
307 if (abs(lat[i]) >= 100) lat_prec--;
308 }
309
310 lat_os.setf(ios::showpoint | ios::fixed);
311 lat_os << setprecision((int)lat_prec) << lat[i];
312
313 Index lon_prec = 4;
314 if (lon[i] < 0) lon_prec--;
315 if (abs(lon[i]) >= 10) {
316 lon_prec--;
317 if (abs(lon[i]) >= 100) lon_prec--;
318 }
319 lon_os.setf(ios::showpoint | ios::fixed);
320 lon_os << setprecision((int)lon_prec) << lon[i];
321
322 sensor_los(0, 0) =
323 180.0 - (asin(refellipsoid[0] * sin(sat_za_from_data[i] * DEG2RAD) /
324 sensor_pos(0, 0))) *
325 RAD2DEG;
326
327 //Reads the t_field_raw from file
328 xml_read_from_file(met_profile_path + "profile.lat_" + lat_os.str() +
329 ".lon_" + lon_os.str() + ".t.xml",
330 t_field_raw,
331 verbosity);
332
333 //Reads the z_field_raw from file
334 xml_read_from_file(met_profile_path + "profile.lat_" + lat_os.str() +
335 ".lon_" + lon_os.str() + ".z.xml",
336 z_field_raw,
337 verbosity);
338
339 //Reads the humidity from file - it is only an ArrayofTensor3
340 xml_read_from_file(met_profile_path + "profile.lat_" + lat_os.str() +
341 ".lon_" + lon_os.str() + ".H2O.xml",
342 vmr_field_raw[0],
343 verbosity);
344
345 //Reads the pnd_field_raw for one scattering element
346 //xml_read_from_file("/rinax/storage/users/rekha/uk_data/profiles/new_obs/newest_forecastfields/reff100/profile.lat_"+lat_os.str()+".lon_"+lon_os.str() + ".pnd100.xml", pnd_field_raw[0]);
347
348 //xml_read_from_file(met_profile_pnd_path +"reff100_newformat/profile.lat_"+lat_os.str()+".lon_"+lon_os.str() + ".pnd100.xml", pnd_field_raw[0]);
349
350 xml_read_from_file(met_profile_pnd_path + "lwc_reff15/profile.lat_" +
351 lat_os.str() + ".lon_" + lon_os.str() + ".pnd15.xml",
352 pnd_field_raw[0],
353 verbosity);
354 //Write the profile number into a file.
355 // xml_write_to_file("profile_number.xml", i);
356
357 // Set z_surface from lowest level of z_field
358 z_surface(0, 0) = z_field_raw.data(0, 0, 0);
359
360 /* The vmr_field_raw is an ArrayofArrayofTensor3 where the outer
361 array is for species.
362
363 The oxygen and nitrogen VMRs are set to constant values of 0.209
364 and 0.782, respectively and are used along with humidity field
365 to generate *vmr_field_raw*.*/
366
367 /*The second element of the species. The first 3 Tensors in the
368 array are the same . They are pressure grid, latitude grid and
369 longitude grid. The third tensor which is the vmr is set to a
370 constant value of 0.782, corresponding to N2.*/
371
372 vmr_field_raw[1].resize(vmr_field_raw[0]);
373 vmr_field_raw[1].copy_grids(vmr_field_raw[0]);
374 vmr_field_raw[1] = 0.782; //vmr of N2
375
376 /*the third element of the species. the first 3 Tensors in the
377 array are the same . They are pressure grid, latitude grid and
378 longitude grid. The third tensor which is the vmr is set to a
379 constant value of 0.209, corresponding to O2.*/
380 vmr_field_raw[2].resize(vmr_field_raw[0]);
381 vmr_field_raw[2].copy_grids(vmr_field_raw[0]);
382 vmr_field_raw[2] = 0.209; //vmr of O2
383
384 const Vector& tfr_p_grid =
385 t_field_raw.get_numeric_grid(GFIELD3_P_GRID);
386 // N_p is the number of elements in the pressure grid
387 Index N_p = tfr_p_grid.nelem();
388
389 //Making a p_grid with the first and last element taken from the profile.
391 p_grid, nelem_p_grid, tfr_p_grid[0], tfr_p_grid[N_p - 1], verbosity);
392
393 /*To set the cloudbox limits, the lower and upper cloudbox limits
394 are to be set. The lower cloudbox limit is set to the lowest
395 pressure level. The upper level is the highest level where the
396 ice water content is non-zero.*/
397 Numeric cl_grid_min, cl_grid_max;
398
399 //Lower limit = lowest pressure level of the original grid.
400 //Could it be the interpolated p_grid? FIXME STR
401 cl_grid_min = tfr_p_grid[0];
402
403 // A counter for non-zero ice content
404 Index level_counter = 0;
405
406 // Loop over all pressure levels
407 for (Index ip = 0; ip < N_p; ++ip) {
408 //Checking for non-zero ice content. 0.001 is a threshold for
409 //ice water content.
410 // if((pnd_field_raw[0].data()(ip, 0, 0) > 0.001) || (pnd_field_raw[1](ip, 0, 0) > 0.001))
411 if (pnd_field_raw[0].data(ip, 0, 0) > 0.001) {
412 ++level_counter;
413 //if non-zero ice content is found, it is set to upper
414 //cloudbox limit. Moreover, we take one level higher
415 // than the upper limit because we want the upper limit
416 //to have 0 pnd.
417 cl_grid_max = tfr_p_grid[ip + 1];
418 }
419 }
420
421 //cloudbox limits have dimensions 2*atmosphere_dim
422 cloudbox_limits.resize(atmosphere_dim * 2);
423
424 //if there is no cloud in the considered profile, still we
425 //need to set the upper limit. I here set the first level
426 //for the upper cloudbox limit.
427 if (level_counter == 0) {
428 cl_grid_max = p_grid[1];
429 }
430
431 //Cloudbox is set.
432 cloudboxSetManually(cloudbox_on,
433 cloudbox_limits,
434 atmosphere_dim,
435 p_grid,
436 lat_grid,
437 lon_grid,
438 cl_grid_min,
439 cl_grid_max,
440 0,
441 0,
442 0,
443 0,
444 verbosity);
445
446 /*executing the met_profile_calc_agenda
447 Agenda communication variables are
448 Output of met_profile_calc_agenda : y
449 Input to met_profile_calc_agenda : t_field_raw,
450 z_field_raw, vmr_field_raw, pnd_field_raw, p_grid,
451 sensor_los, cloudbox_on, cloudbox_limits, z_surface, */
452
454 y,
455 t_field_raw,
456 vmr_field_raw,
457 z_field_raw,
458 pnd_field_raw,
459 p_grid,
460 sensor_los,
461 cloudbox_on,
462 cloudbox_limits,
463 z_surface,
464 met_profile_calc_agenda);
465
466 //putting in the spectra *y* for each profile, thus assigning y
467 //to the ith row of ybatch
468 ybatch[i] = y;
469
470 } // closing the loop over profile basenames
471}
472
473/* Workspace method: Doxygen documentation will be auto-generated */
475 //Output
476 ArrayOfVector& ybatch,
477 //Input
478 const ArrayOfArrayOfSpeciesTag& abs_species,
479 const Agenda& met_profile_calc_agenda,
480 const Vector& f_grid,
481 const Matrix& met_amsu_data,
482 const Matrix& sensor_pos,
483 const Vector& refellipsoid,
484 //Keyword
485 const Index& nelem_p_grid,
486 const String& met_profile_path,
487 const Verbosity& verbosity) {
488 GriddedField3 t_field_raw;
489 GriddedField3 z_field_raw;
490 ArrayOfGriddedField3 vmr_field_raw;
491 ArrayOfGriddedField3 pnd_field_raw;
492 Vector p_grid;
493 Matrix sensor_los;
494 Index cloudbox_on = 0;
495 ArrayOfIndex cloudbox_limits;
496 Matrix z_surface;
497 Vector y;
498 Index no_profiles = met_amsu_data.nrows();
499 //Index no_profiles = met_profile_basenames.nelem();
500 // The humidity data is stored as an ArrayOfTensor3 whereas
501 // vmr_field_raw is an ArrayOfArrayOfTensor3
502 GriddedField3 vmr_field_raw_h2o;
503
504 vmr_field_raw.resize(abs_species.nelem());
505
506 y.resize(f_grid.nelem());
507 ybatch.resize(no_profiles);
508
509 Vector sat_za_from_profile;
510 sat_za_from_profile = met_amsu_data(Range(joker), 3);
511 Numeric sat_za;
512
513 sensor_los.resize(1, 1);
514
515 Vector lat, lon;
516 lat = met_amsu_data(Range(joker), 0);
517 lon = met_amsu_data(Range(joker), 1);
518
519// Vector oro_height;
520// oro_height = met_amsu_data(Range(joker), 5);
521
522 z_surface.resize(1, 1);
523 for (Index i = 0; i < no_profiles; ++i) {
524 ostringstream lat_os, lon_os;
525
526 Index lat_prec = 3;
527 if (lat[i] < 0) lat_prec--;
528 if (abs(lat[i]) >= 10) {
529 lat_prec--;
530 if (abs(lat[i]) >= 100) lat_prec--;
531 }
532
533 lat_os.setf(ios::showpoint | ios::fixed);
534 lat_os << setprecision((int)lat_prec) << lat[i];
535
536 Index lon_prec = 4;
537 if (lon[i] < 0) lon_prec--;
538 if (abs(lon[i]) >= 10) {
539 lon_prec--;
540 if (abs(lon[i]) >= 100) lon_prec--;
541 }
542 lon_os.setf(ios::showpoint | ios::fixed);
543 lon_os << setprecision((int)lon_prec) << lon[i];
544 cout << lat_os.str() << endl;
545 cout << lon_os.str() << endl;
546
547 sat_za = sat_za_from_profile[i];
548
549 sensor_los(Range(joker), 0) =
550 180.0 -
551 (asin(refellipsoid[0] * sin(sat_za * PI / 180.) / sensor_pos(0, 0))) *
552 180. / PI;
553 cout << "sensor_los" << sat_za_from_profile[i] << endl;
554 cout << "sensor_los" << sat_za << endl;
555 cout << "sensor_los" << sensor_los << endl;
556 //Reads the t_field_raw from file
557
558 xml_read_from_file(met_profile_path + "profile.lat_" + lat_os.str() +
559 ".lon_" + lon_os.str() + ".t.xml",
560 t_field_raw,
561 verbosity);
562 //Reads the z_field_raw from file
563 xml_read_from_file(met_profile_path + "profile.lat_" + lat_os.str() +
564 ".lon_" + lon_os.str() + ".z.xml",
565 z_field_raw,
566 verbosity);
567
568 //Reads the humidity from file - it is only an ArrayofTensor3
569 // The vmr_field_raw is an ArrayofArrayofTensor3 where the outer
570 // array is for species
571 xml_read_from_file(met_profile_path + "profile.lat_" + lat_os.str() +
572 ".lon_" + lon_os.str() + ".H2O.xml",
573 vmr_field_raw_h2o,
574 verbosity);
575 //xml_read_from_file("/home/home01/rekha/uk/profiles/sat_vmr/profile.lat_"+lat_os.str()//+".lon_"+lon_os.str() + ".H2O_es.xml",
576 // vmr_field_raw_h2o, verbosity);
577
578 cout
579 << "--------------------------------------------------------------------------"
580 << endl;
581 cout << "The file"
582 << met_profile_path + "profile.lat_" + lat_os.str() + ".lon_" +
583 lon_os.str()
584 << "is executed now" << endl;
585 cout
586 << "--------------------------------------------------------------------------"
587 << endl;
588 xml_write_to_file("profile_number.xml", i, FILE_TYPE_ASCII, 0, verbosity);
589 // the first element of the species is water vapour.
590
591 // N_p is the number of elements in the pressure grid
592 //z_surface(0,0) = oro_height[i]+ 0.01;
593 z_surface(0, 0) = z_field_raw.data(0, 0, 0);
594 cout << "z_surface" << z_surface << endl;
595 const Vector& tfr_p_grid =
596 t_field_raw.get_numeric_grid(GFIELD3_P_GRID);
597 Index N_p = tfr_p_grid.nelem();
598
599 vmr_field_raw[0] = vmr_field_raw_h2o;
600
601 // the second element of the species. the first 3 Tensors in the
602 //array are the same . They are pressure grid, latitude grid and
603 // longitude grid. The third tensor which is the vmr is set to a
604 // constant value of 0.782.
605 vmr_field_raw[1].resize(vmr_field_raw[0]);
606 vmr_field_raw[1].copy_grids(vmr_field_raw[0]);
607 vmr_field_raw[1].data(joker, joker, joker) = 0.782;
608
609 // the second element of the species. the first 3 Tensors in the
610 //array are the same . They are pressure grid, latitude grid and
611 // longitude grid. The third tensor which is the vmr is set to a
612 // constant value of 0.209.
613 vmr_field_raw[2].resize(vmr_field_raw[0]);
614 vmr_field_raw[2].copy_grids(vmr_field_raw[0]);
615 vmr_field_raw[2].data(joker, joker, joker) = 0.209;
616
617 //xml_write_to_file(met_profile_basenames[i]+ ".N2.xml", vmr_field_raw[1]);
618 //xml_write_to_file(met_profile_basenames[i]+ ".O2.xml", vmr_field_raw[2]);
619
620 //Making a p_grid with the first and last element taken from the profile.
621 // this is because of the extrapolation problem.
622
624 p_grid, nelem_p_grid, tfr_p_grid[0], tfr_p_grid[N_p - 1], verbosity);
625 cout << "t_field_raw[0](0,0,0)" << tfr_p_grid[0] << endl;
626 cout << "t_field_raw[0](N_p -1,0,0)" << tfr_p_grid[N_p - 1] << endl;
627 xml_write_to_file("p_grid.xml", p_grid, FILE_TYPE_ASCII, 0, verbosity);
628
629 // executing the met_profile_calc_agenda
631 y,
632 t_field_raw,
633 vmr_field_raw,
634 z_field_raw,
635 pnd_field_raw,
636 p_grid,
637 sensor_los,
638 cloudbox_on,
639 cloudbox_limits,
640 z_surface,
641 met_profile_calc_agenda);
642
643 //putting in the spectra *y* for each profile
644 ybatch[i] = y;
645
646 } // closing the loop over profile basenames
647}
648
649/* Workspace method: Doxygen documentation will be auto-generated */
651 ArrayOfTensor7& dobatch_cloudbox_field,
652 ArrayOfTensor5& dobatch_radiance_field,
653 ArrayOfTensor4& dobatch_irradiance_field,
654 ArrayOfTensor5& dobatch_spectral_irradiance_field,
655 const Index& ybatch_start,
656 const Index& ybatch_n,
657 const Agenda& dobatch_calc_agenda,
658 const Index& robust,
659 const Verbosity& verbosity) {
661
662 Index first_ybatch_index = 0;
663
664 ArrayOfString fail_msg;
665 bool do_abort = false;
666
667 // We allow a start index ybatch_start that is different from 0. We
668 // will calculate ybatch_n jobs starting at the start
669 // index. Internally, we count from zero, which is the right
670 // counting for the output array ybatch. When we call
671 // ybatch_calc_agenda, we add ybatch_start to the internal index
672 // count.
673
674 // We create a counter, so that we can generate nice output about
675 // how many jobs are already done. (All parallel threads later will
676 // increment this, so that we really get an accurate total count!)
677 Index job_counter = 0;
678
679 // Resize the output arrays:
680 dobatch_cloudbox_field.resize(ybatch_n);
681 dobatch_radiance_field.resize(ybatch_n);
682 dobatch_irradiance_field.resize(ybatch_n);
683 dobatch_spectral_irradiance_field.resize(ybatch_n);
684
685 // Go through the batch:
686
687 if (ybatch_n) {
689
690#pragma omp parallel for schedule(dynamic) if (!arts_omp_in_parallel() && \
691 ybatch_n > 1) firstprivate(wss)
692 for (Index ybatch_index = first_ybatch_index; ybatch_index < ybatch_n;
693 ybatch_index++) {
694 Index l_job_counter; // Thread-local copy of job counter.
695
696 if (do_abort) continue;
697#pragma omp critical(dobatchCalc_job_counter)
698 { l_job_counter = ++job_counter; }
699
700 {
701 ostringstream os;
702 os << " Job " << l_job_counter << " of " << ybatch_n << ", Index "
703 << ybatch_start + ybatch_index << ", Thread-Id "
704 << arts_omp_get_thread_num() << "\n";
705 out2 << os.str();
706 }
707
708 try {
709 Tensor7 cloudbox_field;
710 Tensor5 radiance_field;
711 Tensor4 irradiance_field;
712 Tensor5 spectral_irradiance_field;
713
715 cloudbox_field,
716 radiance_field,
717 irradiance_field,
718 spectral_irradiance_field,
719 ybatch_start + ybatch_index,
720 dobatch_calc_agenda);
721
722#pragma omp critical(dobatchCalc_assign_cloudbox_field)
723 dobatch_cloudbox_field[ybatch_index] = cloudbox_field;
724#pragma omp critical(dobatchCalc_assign_radiance_field)
725 dobatch_radiance_field[ybatch_index] = radiance_field;
726#pragma omp critical(dobatchCalc_assign_irradiance_field)
727 dobatch_irradiance_field[ybatch_index] = irradiance_field;
728#pragma omp critical(dobatchCalc_assign_spectral_irradiance_field)
729 dobatch_spectral_irradiance_field[ybatch_index] =
730 spectral_irradiance_field;
731
732 } catch (const std::exception& e) {
733 if (robust && !do_abort) {
734 ostringstream os;
735 os << "WARNING! Job at ybatch_index " << ybatch_start + ybatch_index
736 << " failed.\n"
737 << "element in output variables will be empty for this job.\n"
738 << "The runtime error produced was:\n"
739 << e.what() << "\n";
740 out0 << os.str();
741 } else {
742 // The user wants the batch job to fail if one of the
743 // jobs goes wrong.
744#pragma omp critical(dobatchCalc_setabort)
745 do_abort = true;
746
747 ostringstream os;
748 os << " Job at ybatch_index " << ybatch_start + ybatch_index
749 << " failed. Aborting...\n";
750 out1 << os.str();
751 }
752 ostringstream os;
753 os << "Run-time error at ybatch_index " << ybatch_start + ybatch_index
754 << ": \n"
755 << e.what();
756#pragma omp critical(dobatchCalc_push_fail_msg)
757 fail_msg.push_back(os.str());
758 }
759 }
760 }
761
762 if (fail_msg.nelem()) {
763 ostringstream os;
764
765 if (!do_abort) os << "\nError messages from failed batch cases:\n";
766 for (ArrayOfString::const_iterator it = fail_msg.begin();
767 it != fail_msg.end();
768 it++)
769 os << *it << '\n';
770
771 if (do_abort)
772 throw runtime_error(os.str());
773 else
774 out0 << os.str();
775 }
776}
Index TotalNumberOfElements(const Array< Array< base > > &aa)
Determine total number of elements in an ArrayOfArray.
Definition: array.h:241
The global header file for ARTS.
int arts_omp_get_thread_num()
Wrapper for omp_get_thread_num.
Definition: arts_omp.cc:76
Header file for helper functions for OpenMP.
void forloop_agendaExecute(Workspace &ws, const Index forloop_index, const Agenda &input_agenda)
Definition: auto_md.cc:24792
void dobatch_calc_agendaExecute(Workspace &ws, Tensor7 &spectral_radiance_field, Tensor5 &radiance_field, Tensor4 &irradiance_field, Tensor5 &spectral_irradiance_field, const Index ybatch_index, const Agenda &input_agenda)
Definition: auto_md.cc:24620
void ybatch_calc_agendaExecute(Workspace &ws, Vector &y, ArrayOfVector &y_aux, Matrix &jacobian, const Index ybatch_index, const Agenda &input_agenda)
Definition: auto_md.cc:26076
void met_profile_calc_agendaExecute(Workspace &ws, Vector &y, const GriddedField3 &t_field_raw, const ArrayOfGriddedField3 &vmr_field_raw, const GriddedField3 &z_field_raw, const ArrayOfGriddedField3 &pnd_field_raw, const Vector &p_grid, const Matrix &sensor_los, const Index cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Matrix &z_surface, const Agenda &input_agenda)
Definition: auto_md.cc:25441
The Agenda class.
Definition: agenda_class.h:69
Index nelem() const ARTS_NOEXCEPT
Definition: array.h:92
Index nrows() const noexcept
Definition: matpackI.h:1079
Index ncols() const noexcept
Definition: matpackI.h:1080
A constant view of a Vector.
Definition: matpackI.h:521
Index nelem() const noexcept
Returns the number of elements.
Definition: matpackI.h:547
const Vector & get_numeric_grid(Index i) const
Get a numeric grid.
The Matrix class.
Definition: matpackI.h:1285
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1011
The range class.
Definition: matpackI.h:160
The Tensor4 class.
Definition: matpackIV.h:435
The Tensor5 class.
Definition: matpackV.h:524
The Tensor7 class.
Definition: matpackVII.h:2407
The Vector class.
Definition: matpackI.h:910
void resize(Index n)
Resize function.
Definition: matpackI.cc:390
Workspace class.
Definition: workspace_ng.h:53
Implementation of gridded fields.
void VectorNLogSpace(Vector &x, const Index &n, const Numeric &start, const Numeric &stop, const Verbosity &verbosity)
WORKSPACE METHOD: VectorNLogSpace.
constexpr Numeric DEG2RAD
Definition: m_batch.cc:52
constexpr Numeric RAD2DEG
Definition: m_batch.cc:53
void ForLoop(Workspace &ws, const Agenda &forloop_agenda, const Index &start, const Index &stop, const Index &step, const Verbosity &verbosity)
WORKSPACE METHOD: ForLoop.
Definition: m_batch.cc:63
void ybatchMetProfilesClear(Workspace &ws, ArrayOfVector &ybatch, const ArrayOfArrayOfSpeciesTag &abs_species, const Agenda &met_profile_calc_agenda, const Vector &f_grid, const Matrix &met_amsu_data, const Matrix &sensor_pos, const Vector &refellipsoid, const Index &nelem_p_grid, const String &met_profile_path, const Verbosity &verbosity)
WORKSPACE METHOD: ybatchMetProfilesClear.
Definition: m_batch.cc:474
void ybatchCalc(Workspace &ws, ArrayOfVector &ybatch, ArrayOfArrayOfVector &ybatch_aux, ArrayOfMatrix &ybatch_jacobians, const Index &ybatch_start, const Index &ybatch_n, const Agenda &ybatch_calc_agenda, const Index &robust, const Verbosity &verbosity)
WORKSPACE METHOD: ybatchCalc.
Definition: m_batch.cc:80
void ybatchMetProfiles(Workspace &ws, ArrayOfVector &ybatch, const ArrayOfArrayOfSpeciesTag &abs_species, const Agenda &met_profile_calc_agenda, const Vector &f_grid, const Matrix &met_amsu_data, const Matrix &sensor_pos, const Vector &refellipsoid, const Vector &lat_grid, const Vector &lon_grid, const Index &atmosphere_dim, const ArrayOfArrayOfSingleScatteringData &scat_data, const Index &nelem_p_grid, const String &met_profile_path, const String &met_profile_pnd_path, const Verbosity &verbosity)
WORKSPACE METHOD: ybatchMetProfiles.
Definition: m_batch.cc:236
constexpr Numeric PI
Definition: m_batch.cc:51
void DOBatchCalc(Workspace &ws, ArrayOfTensor7 &dobatch_cloudbox_field, ArrayOfTensor5 &dobatch_radiance_field, ArrayOfTensor4 &dobatch_irradiance_field, ArrayOfTensor5 &dobatch_spectral_irradiance_field, const Index &ybatch_start, const Index &ybatch_n, const Agenda &dobatch_calc_agenda, const Index &robust, const Verbosity &verbosity)
WORKSPACE METHOD: DOBatchCalc.
Definition: m_batch.cc:650
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:347
void abs(Sparse &A, const Sparse &B)
Absolute value of sparse matrix elements.
Definition: matpackII.cc:394
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
const Joker joker
#define CREATE_OUT1
Definition: messages.h:204
#define CREATE_OUTS
Definition: messages.h:208
constexpr Numeric pi
The following mathematical constants are generated in python Decimal package by the code:
constexpr auto deg2rad(auto x) noexcept
Converts degrees to radians.
constexpr auto rad2deg(auto x) noexcept
Converts radians to degrees.
constexpr Index GFIELD3_P_GRID
Global constant, Index of the pressure grid in GriddedField3.
VectorView std(VectorView std, const Vector &y, const ArrayOfVector &ys, const Index start, const Index end_tmp)
Compute the standard deviation of the ranged ys.
Definition: raw.cc:205
This file contains declerations of functions of physical character.
Declaration of functions in rte.cc.
This file contains basic functions to handle XML data files.
void xml_write_to_file(const String &filename, const T &type, const FileType ftype, const Index no_clobber, const Verbosity &verbosity)
Write data to XML file.
Definition: xml_io.h:172
void xml_read_from_file(const String &filename, T &type, const Verbosity &verbosity)
Reads data from XML file.
Definition: xml_io.h:151
@ FILE_TYPE_ASCII
Definition: xml_io_base.h:43