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