ARTS  2.0.49
m_cloudbox.cc
Go to the documentation of this file.
1 /* Copyright (C) 2002-2008
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 
7  This program is free software; you can redistribute it and/or modify it
8  under the terms of the GNU General Public License as published by the
9  Free Software Foundation; either version 2, or (at your option) any
10  later version.
11 
12  This program is distributed in the hope that it will be useful,
13  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  GNU General Public License for more details.
16 
17  You should have received a copy of the GNU General Public License
18  along with this program; if not, write to the Free Software
19  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
20  USA. */
21 
22 
23 
24 /*===========================================================================
25  === File description
26  ===========================================================================*/
27 
40 /*===========================================================================
41  === External declarations
42  ===========================================================================*/
43 
44 #include <stdexcept>
45 #include <cstdlib>
46 #include <cmath>
47 
48 #include "arts.h"
49 #include "array.h"
50 #include "auto_md.h"
51 #include "check_input.h"
52 #include "xml_io.h"
53 #include "messages.h"
54 #include "gridded_fields.h"
55 #include "logic.h"
56 #include "rte.h"
57 #include "interpolation.h"
58 #include "special_interp.h"
59 #include "cloudbox.h"
60 #include "optproperties.h"
61 #include "math_funcs.h"
62 #include "physics_funcs.h"
63 #include "sorting.h"
64 
65 extern const Index GFIELD3_P_GRID;
66 extern const Index GFIELD3_LAT_GRID;
67 extern const Index GFIELD3_LON_GRID;
68 extern const Numeric DEG2RAD;
69 extern const Numeric PI;
70 
71 
72 /*===========================================================================
73  === The functions (in alphabetical order)
74  ===========================================================================*/
75 
76 /* Workspace method: Doxygen documentation will be auto-generated */
77 void cloudboxOff (// WS Output:
78  Index& cloudbox_on,
79  ArrayOfIndex& cloudbox_limits,
80  Agenda& iy_cloudbox_agenda,
81  const Verbosity&)
82 {
83  cloudbox_on = 0;
84  cloudbox_limits.resize ( 0 );
85  iy_cloudbox_agenda = Agenda();
86  iy_cloudbox_agenda.set_name ( "iy_cloudbox_agenda" );
87 }
88 
89 
90 
91 /* Workspace method: Doxygen documentation will be auto-generated */
92 void cloudboxSetAutomatically (// WS Output:
93  //Workspace& /* ws */,
94  Index& cloudbox_on,
95  ArrayOfIndex& cloudbox_limits,
96  //Agenda& iy_cloudbox_agenda,
97  // WS Input:
98  const Index& atmosphere_dim,
99  const ArrayOfString& part_species,
100  const Vector& p_grid,
101  const Vector& lat_grid,
102  const Vector& lon_grid,
103  const Tensor4& massdensity_field,
104  // Control Parameters
105  const Numeric& cloudbox_margin,
106  const Verbosity& verbosity)
107 {
108  // Variables
109  //const Index np = massdensity_field.npages();
110  Index p1 = massdensity_field.npages()-1;
111  Index p2 = 0;
112  DEBUG_ONLY(
113  Index lat1 = massdensity_field.nrows()-1;
114  Index lat2 = 0;
115  Index lon1 = massdensity_field.ncols()-1;
116  Index lon2 = 0;
117  )
118 
119  Index type_flag=0, i=0, j=0, k=0, l=0;
120  // initialize flag, telling if all selected *massdensity_fields* are
121  // zero(false) or not(true)
122  bool x = false;
123 
124 
125  // Check existing WSV
126  chk_if_in_range ( "atmosphere_dim", atmosphere_dim, 1, 3 );
127  // includes p_grid chk_if_decresasing
128  chk_atm_grids ( atmosphere_dim, p_grid, lat_grid, lon_grid );
129  // Set cloudbox_on
130  cloudbox_on = 1;
131 
132  // Allocate cloudbox_limits
133  cloudbox_limits.resize ( atmosphere_dim*2 );
134 
135  //--------- Start loop over particles ----------------------------------------
136  for ( l=0; l<part_species.nelem(); l++ )
137  {
138  String part_type;
139 
140  //split String and copy to ArrayOfString
141  // split part_species string at "-" and write to ArrayOfString
142  parse_part_type( part_type, part_species[l]);
143 
144  // select hydrometeor type according to user input
145  // Index nhyd, describes the column index in *massdensity_field*
146  if ( part_type == "LWC" ) type_flag = 0;
147  else if ( part_type == "IWC" ) type_flag = 1;
148  else if ( part_type == "Rain" ) type_flag = 2;
149  else if ( part_type == "Snow" ) type_flag = 3;
150 
151  bool y; // empty_flag
152  //y is set to true, if a single value of massdensity_field is unequal zero
154  atmosphere_dim,
155  massdensity_field ( type_flag, joker, joker, joker ),
156  p_grid,
157  lat_grid,
158  lon_grid );
159 
160  //-----------Start setting cloudbox limits----------------------------------
161  if ( y )
162  {
163  //massdensity_field unequal zero -> x is true
164  x = true;
165 
166  if ( atmosphere_dim == 1 )
167  {
168  // Pressure limits
169  ConstVectorView hydro_p = massdensity_field ( type_flag, joker, 0 , 0 );
170 
171 
172  // set lower cloudbox_limit to surface if margin = -1
173  if ( cloudbox_margin == -1 )
174  {
175  cloudbox_limits[0] = 0;
176  i = p1 = 0;
177  }
178  else
179  {
180  // find index of first pressure level where hydromet_field is
181  // unequal 0, starting from the surface
182  for ( i=0; i<hydro_p.nelem(); i++ )
183  {
184  if ( hydro_p[i] != 0.0 )
185  {
186  // check if p1 is the lowest index in all selected
187  // massdensity fields
188  if ( p1 > i )
189  {
190  p1 = i;
191  }
192  break;
193  }
194  }
195 
196  }
197  // find index of highest pressure level, where massdensity_field is
198  // unequal 0, starting from top of the atmosphere
199  for ( j=hydro_p.nelem()-1; j>=i; j-- )
200  {
201  if ( hydro_p[j] != 0.0 )
202  {
203  // check if p2 is the highest index in all selected
204  // massdensity fields
205  if ( p2 < j )
206  {
207  p2 = j;
208  }
209  break;
210  }
211  }
212 
213  }
214  }
215 
216  /* //NOT WORKING YET
217  // Latitude limits
218  else if ( atmosphere_dim == 2 )
219  {
220  MatrixView hydro_lat = hydromet_field ( nhyd, joker, joker, 0 );
221 
222  for ( i=0; i<hydro_lat.nrows(); i++ )
223  {
224  for ( j=0; j<hydro_lat.ncols(); j++ )
225  {
226  if ( hydro_lat[i,j] != 0.0 )
227  {
228 
229  if ( lat1 <= j ) lat1 =j;
230  //cloudbox_limits[2] = lat1;
231  //break;
232  }
233 
234  }
235  if ( p1 <= i ) p1 = i;
236  }
237 
238  for ( k=hydro_lat.nelem()-1; k>=i; k-- )
239  {
240  if ( hydro_lat[k] != 0.0 )
241  {
242  lat2 = k;
243  cloudbox_limits[3] = lat2;
244  break;
245 
246  }
247 
248  }
249  }
250 
251  // Longitude limits
252  if ( atmosphere_dim == 3 )
253  {
254  Tensor3View hydro_lon = hydromet_field ( nhyd, joker, joker, joker );
255 
256  for ( i=0; i<hydro_lon.nelem(); i++ )
257  {
258  if ( hydro_lon[i] != 0.0 )
259  {
260  lon1 = i;
261  cloudbox_limits[4] = lon1;
262  break;
263  }
264 
265  }
266  for ( j=hydro_lon.nelem()-1; j>=i; j-- )
267  {
268  if ( hydro_lon[j] != 0.0 )
269  {
270  lon2 = j;
271  cloudbox_limits[5] = lon2;
272  break;
273 
274  }
275 
276 
277 }*/
278  }
279  // decrease lower cb limit by one to ensure that linear interpolation of
280  // particle number densities is possible.
281  Index p0 = 0; //only for the use of function *max*
282  p1 = max(p1-1, p0);
283 
284  Numeric p_margin1;
285 
286  // alter lower cloudbox_limit by cloudbox_margin, using barometric
287  // height formula
288  p_margin1 = barometric_heightformula ( p_grid[p1], cloudbox_margin );
289  while ( p_grid[k+1] >= p_margin1 && k+1 < p_grid.nelem() ) k++;
290  cloudbox_limits[0]= k;
291 
292  // increase upper cb limit by one to ensure that linear interpolation of
293  // particle number densities is possible.
294  p2 = min(p2+1, massdensity_field.npages()-1);
295  // set upper cloudbox_limit
296  // if cloudbox reaches to the upper most pressure level
297  if ( p2 >= massdensity_field.npages()-1)
298  {
300  out2<<"The cloud reaches to TOA!\n"
301  <<"Check massdensity_field data, if realistic!\n";
302  }
303  cloudbox_limits[1] = p2;
304 
305  //out0<<"\n"<<p2<<"\n"<<p_grid[p2]<<"\n";
306 
307  // check if all selected massdensity fields are zero at each level,
308  // than switch cloudbox off, skipping scattering calculations
309  if ( !x )
310  {
312  //cloudboxOff ( cloudbox_on, cloudbox_limits, iy_cloudbox_agenda );
313  cloudbox_on = 0;
314  out0<<"Cloudbox is switched off!\n";
315 
316  return;
317  }
318 
319 
320  // assert keyword arguments
321 
322  // The pressure in *p1* must be bigger than the pressure in *p2*.
323  assert ( p_grid[p1] > p_grid[p2] );
324  // The pressure in *p1* must be larger than the last value in *p_grid*.
325  assert ( p_grid[p1] > p_grid[p_grid.nelem()-1] );
326  // The pressure in *p2* must be smaller than the first value in *p_grid*."
327  assert ( p_grid[p2] < p_grid[0] );
328 
329  if ( atmosphere_dim >= 2 )
330  {
331  // The latitude in *lat2* must be bigger than the latitude in *lat1*.
332  assert ( lat_grid[lat2] > lat_grid[lat1] );
333  // The latitude in *lat1* must be >= the second value in *lat_grid*.
334  assert ( lat_grid[lat1] >= lat_grid[1] );
335  // The latitude in *lat2* must be <= the next to last value in *lat_grid*.
336  assert ( lat_grid[lat2] <= lat_grid[lat_grid.nelem()-2] );
337  }
338  if ( atmosphere_dim == 3 )
339  {
340  // The longitude in *lon2* must be bigger than the longitude in *lon1*.
341  assert ( lon_grid[lon2] > lon_grid[lon1] );
342  // The longitude in *lon1* must be >= the second value in *lon_grid*.
343  assert ( lon_grid[lon1] >= lon_grid[1] );
344  // The longitude in *lon2* must be <= the next to last value in *lon_grid*.
345  assert ( lon_grid[lon2] <= lon_grid[lon_grid.nelem()-2] );
346  }
347 }
348 
349 /* Workspace method: Doxygen documentation will be auto-generated */
350 void cloudboxSetManually(// WS Output:
351  Index& cloudbox_on,
352  ArrayOfIndex& cloudbox_limits,
353  // WS Input:
354  const Index& atmosphere_dim,
355  const Vector& p_grid,
356  const Vector& lat_grid,
357  const Vector& lon_grid,
358  // Control Parameters:
359  const Numeric& p1,
360  const Numeric& p2,
361  const Numeric& lat1,
362  const Numeric& lat2,
363  const Numeric& lon1,
364  const Numeric& lon2,
365  const Verbosity&)
366 {
367  // Check existing WSV
368  chk_if_in_range( "atmosphere_dim", atmosphere_dim, 1, 3 );
369  chk_atm_grids( atmosphere_dim, p_grid, lat_grid, lon_grid );
370 
371  // Check keyword arguments
372  if( p1 <= p2 )
373  throw runtime_error( "The pressure in *p1* must be bigger than the "
374  "pressure in *p2*." );
375  if( p1 <= p_grid[p_grid.nelem()-1] )
376  throw runtime_error( "The pressure in *p1* must be larger than the "
377  "last value in *p_grid*." );
378  if( p2 >= p_grid[0] )
379  throw runtime_error( "The pressure in *p2* must be smaller than the "
380  "first value in *p_grid*." );
381  if( atmosphere_dim >= 2 )
382  {
383  if( lat2 <= lat1 )
384  throw runtime_error( "The latitude in *lat2* must be bigger than the "
385  "latitude in *lat1*.");
386  if( lat1 < lat_grid[1] )
387  throw runtime_error( "The latitude in *lat1* must be >= the "
388  "second value in *lat_grid*." );
389  if( lat2 > lat_grid[lat_grid.nelem()-2] )
390  throw runtime_error( "The latitude in *lat2* must be <= the "
391  "next to last value in *lat_grid*." );
392  }
393  if( atmosphere_dim == 3 )
394  {
395  if( lon2 <= lon1 )
396  throw runtime_error( "The longitude in *lon2* must be bigger than the "
397  "longitude in *lon1*.");
398  if( lon1 < lon_grid[1] )
399  throw runtime_error( "The longitude in *lon1* must be >= the "
400  "second value in *lon_grid*." );
401  if( lon2 > lon_grid[lon_grid.nelem()-2] )
402  throw runtime_error( "The longitude in *lon2* must be <= the "
403  "next to last value in *lon_grid*." );
404  }
405 
406  // Set cloudbox_on
407  cloudbox_on = 1;
408 
409  // Allocate cloudbox_limits
410  cloudbox_limits.resize( atmosphere_dim*2 );
411 
412  // Pressure limits
413  if( p1 > p_grid[1] )
414  {
415  cloudbox_limits[0] = 0;
416  }
417  else
418  {
419  for( cloudbox_limits[0]=1; p_grid[cloudbox_limits[0]+1]>=p1;
420  cloudbox_limits[0]++ ) {}
421  }
422  if( p2 < p_grid[p_grid.nelem()-2] )
423  {
424  cloudbox_limits[1] = p_grid.nelem() - 1;
425  }
426  else
427  {
428  for( cloudbox_limits[1]=p_grid.nelem()-2;
429  p_grid[cloudbox_limits[1]-1]<=p2; cloudbox_limits[1]-- ) {}
430  }
431 
432  // Latitude limits
433  if( atmosphere_dim >= 2 )
434  {
435  for( cloudbox_limits[2]=1; lat_grid[cloudbox_limits[2]+1]<=lat1;
436  cloudbox_limits[2]++ ) {}
437  for( cloudbox_limits[3]=lat_grid.nelem()-2;
438  lat_grid[cloudbox_limits[3]-1]>=lat2; cloudbox_limits[3]-- ) {}
439  }
440 
441  // Longitude limits
442  if( atmosphere_dim == 3 )
443  {
444  for( cloudbox_limits[4]=1; lon_grid[cloudbox_limits[4]+1]<=lon1;
445  cloudbox_limits[4]++ ) {}
446  for( cloudbox_limits[5]=lon_grid.nelem()-2;
447  lon_grid[cloudbox_limits[5]-1]>=lon2; cloudbox_limits[5]-- ) {}
448  }
449 }
450 
451 
452 /* Workspace method: Doxygen documentation will be auto-generated */
453 void cloudboxSetManuallyAltitude(// WS Output:
454  Index& cloudbox_on,
455  ArrayOfIndex& cloudbox_limits,
456  // WS Input:
457  const Index& atmosphere_dim,
458  const Tensor3& z_field,
459  const Vector& lat_grid,
460  const Vector& lon_grid,
461  // Control Parameters:
462  const Numeric& z1,
463  const Numeric& z2,
464  const Numeric& lat1,
465  const Numeric& lat2,
466  const Numeric& lon1,
467  const Numeric& lon2,
468  const Verbosity&)
469 {
470  // Check existing WSV
471  chk_if_in_range( "atmosphere_dim", atmosphere_dim, 1, 3 );
472 
473  // Check keyword arguments
474  if( z1 >= z2 )
475  throw runtime_error( "The altitude in *z1* must be smaller than the "
476  "altitude in *z2*." );
477  /* These cases are in fact handled
478  if( z1 <= z_field(0, 0, 0) )
479  throw runtime_error( "The altitude in *z1* must be larger than the "
480  "first value in *z_field*." );
481  if( z2 >= z_field(z_field.npages()-1, 0, 0) )
482  throw runtime_error( "The altitude in *z2* must be smaller than the "
483  "last value in *z_field*." );
484  */
485  if( atmosphere_dim == 3 )
486  {
487  if( lat2 <= lat1 )
488  throw runtime_error( "The latitude in *lat2* must be bigger than the "
489  " latitude in *lat1*.");
490  if( lat1 < lat_grid[1] )
491  throw runtime_error( "The latitude in *lat1* must be >= the "
492  "second value in *lat_grid*." );
493  if( lat2 > lat_grid[lat_grid.nelem()-2] )
494  throw runtime_error( "The latitude in *lat2* must be <= the "
495  "next to last value in *lat_grid*." );
496  if( lon2 <= lon1 )
497  throw runtime_error( "The longitude in *lon2* must be bigger than the "
498  "longitude in *lon1*.");
499  if( lon1 < lon_grid[1] )
500  throw runtime_error( "The longitude in *lon1* must be >= the "
501  "second value in *lon_grid*." );
502  if( lon2 > lon_grid[lon_grid.nelem()-2] )
503  throw runtime_error( "The longitude in *lon2* must be <= the "
504  "next to last value in *lon_grid*." );
505  }
506 
507  // Set cloudbox_on
508  cloudbox_on = 1;
509 
510  // Allocate cloudbox_limits
511  cloudbox_limits.resize( atmosphere_dim*2 );
512 
513  // Pressure/altitude limits
514  if( z1 < z_field(1, 0, 0) )
515  {
516  cloudbox_limits[0] = 0;
517  }
518  else
519  {
520  for( cloudbox_limits[0]=1; z_field(cloudbox_limits[0]+1, 0, 0) <= z1;
521  cloudbox_limits[0]++ ) {}
522  }
523  if( z2 > z_field(z_field.npages()-2, 0, 0) )
524  {
525  cloudbox_limits[1] = z_field.npages() - 1;
526  }
527  else
528  {
529  for( cloudbox_limits[1]=z_field.npages()- 2;
530  z_field(cloudbox_limits[1]-1, 0, 0) >= z2; cloudbox_limits[1]-- ) {}
531  }
532 
533  // Latitude limits
534  if( atmosphere_dim >= 2 )
535  {
536  for( cloudbox_limits[2]=1; lat_grid[cloudbox_limits[2]+1]<=lat1;
537  cloudbox_limits[2]++ ) {}
538  for( cloudbox_limits[3]=lat_grid.nelem()-2;
539  lat_grid[cloudbox_limits[3]-1]>=lat2; cloudbox_limits[3]-- ) {}
540  }
541 
542  // Longitude limits
543  if( atmosphere_dim == 3 )
544  {
545  for( cloudbox_limits[4]=1; lon_grid[cloudbox_limits[4]+1]<=lon1;
546  cloudbox_limits[4]++ ) {}
547  for( cloudbox_limits[5]=lon_grid.nelem()-2;
548  lon_grid[cloudbox_limits[5]-1]>=lon2; cloudbox_limits[5]-- ) {}
549  }
550 }
551 
552 
553 
554 /* Workspace method: Doxygen documentation will be auto-generated */
555 void cloudbox_checkedCalc(Index& cloudbox_checked,
556  const Index& basics_checked,
557  const Index& atmosphere_dim,
558  const Vector& p_grid,
559  const Vector& lat_grid,
560  const Vector& lon_grid,
561  const Tensor3& wind_u_field,
562  const Tensor3& wind_v_field,
563  const Tensor3& wind_w_field,
564  const Index& cloudbox_on,
565  const ArrayOfIndex& cloudbox_limits,
566  const Verbosity&)
567 {
568  // Demanded space between cloudbox and lat and lon edges [degrees]
569  const Numeric llmin = 20;
570 
571  if( !basics_checked )
572  throw runtime_error( "The atmosphere and basic control varaibles must be "
573  "flagged to have passed a consistency check (basics_checked=1)." );
574 
575  chk_if_bool( "cloudbox_on", cloudbox_on );
576 
577  if( cloudbox_on )
578  {
579  // Winds, must be empty variables (i.e. no winds allowed)
580  {
581  ostringstream ow;
582  ow << "The scattering methods are not (yet?) handling winds. For this\n"
583  << "reason, the WSVs for wind fields must all be empty with an\n."
584  << "active cloudbox.";
585  if( wind_w_field.npages() > 0 )
586  { throw runtime_error( ow.str() ); }
587  if( wind_v_field.npages() > 0 )
588  { throw runtime_error( ow.str() ); }
589  if( atmosphere_dim > 2 && wind_u_field.npages() > 0 )
590  { throw runtime_error( ow.str() ); }
591  }
592 
593  // Cloudbox limits
594  if( cloudbox_limits.nelem() != atmosphere_dim*2 )
595  {
596  ostringstream os;
597  os << "The array *cloudbox_limits* has incorrect length.\n"
598  << "For atmospheric dim. = " << atmosphere_dim
599  << " the length shall be " << atmosphere_dim*2
600  << " but it is " << cloudbox_limits.nelem() << ".";
601  throw runtime_error( os.str() );
602  }
603  if( cloudbox_limits[1]<=cloudbox_limits[0] || cloudbox_limits[0]<0 ||
604  cloudbox_limits[1]>=p_grid.nelem() )
605  {
606  ostringstream os;
607  os << "Incorrect value(s) for cloud box pressure limit(s) found."
608  << "\nValues are either out of range or upper limit is not "
609  << "greater than lower limit.\nWith present length of "
610  << "*p_grid*, OK values are 0 - " << p_grid.nelem()-1
611  << ".\nThe pressure index limits are set to "
612  << cloudbox_limits[0] << " - " << cloudbox_limits[1] << ".";
613  throw runtime_error( os.str() );
614  }
615  if( atmosphere_dim >= 2 )
616  {
617  const Index n = lat_grid.nelem();
618  if( cloudbox_limits[3]<=cloudbox_limits[2] || cloudbox_limits[2]<1 ||
619  cloudbox_limits[3]>=n-1 )
620  {
621  ostringstream os;
622  os << "Incorrect value(s) for cloud box latitude limit(s) found."
623  << "\nValues are either out of range or upper limit is not "
624  << "greater than lower limit.\nWith present length of "
625  << "*lat_grid*, OK values are 1 - " << n-2
626  << ".\nThe latitude index limits are set to "
627  << cloudbox_limits[2] << " - " << cloudbox_limits[3] << ".";
628  throw runtime_error( os.str() );
629  }
630  if( ( lat_grid[cloudbox_limits[2]] - lat_grid[0] < llmin ) &&
631  ( atmosphere_dim==2 ||
632  ( atmosphere_dim==3 && lat_grid[0]>-90) ) )
633  {
634  ostringstream os;
635  os << "Too small distance between cloudbox and lower end of\n"
636  << "latitude grid. This distance must be " << llmin
637  << " degrees. Cloudbox ends at " << lat_grid[cloudbox_limits[2]]
638  << " and latitude grid starts at " << lat_grid[0] << ".";
639  throw runtime_error( os.str() );
640  }
641  if( ( lat_grid[n-1] - lat_grid[cloudbox_limits[3]] < llmin ) &&
642  ( atmosphere_dim==2 ||
643  (atmosphere_dim==3 && lat_grid[n-1]<90) ) )
644  {
645  ostringstream os;
646  os << "Too small distance between cloudbox and upper end of\n"
647  << "latitude grid. This distance must be " << llmin
648  << "degrees. Cloudbox ends at " << lat_grid[cloudbox_limits[3]]
649  << " and latitude grid ends at " << lat_grid[n-1] << ".";
650  throw runtime_error( os.str() );
651  }
652  }
653  if( atmosphere_dim >= 3 )
654  {
655  const Index n = lon_grid.nelem();
656  if( cloudbox_limits[5]<=cloudbox_limits[4] || cloudbox_limits[4]<1 ||
657  cloudbox_limits[5]>=n-1 )
658  {
659  ostringstream os;
660  os << "Incorrect value(s) for cloud box longitude limit(s) found"
661  << ".\nValues are either out of range or upper limit is not "
662  << "greater than lower limit.\nWith present length of "
663  << "*lon_grid*, OK values are 1 - " << n-2
664  << ".\nThe longitude limits are set to "
665  << cloudbox_limits[4] << " - " << cloudbox_limits[5] << ".";
666  throw runtime_error( os.str() );
667  }
668  if( lon_grid[n-1] - lon_grid[0] < 360 )
669  {
670  const Numeric latmax = max( abs(lat_grid[cloudbox_limits[2]]),
671  abs(lat_grid[cloudbox_limits[3]]) );
672  const Numeric lfac = 1 / cos( DEG2RAD*latmax );
673  if( lon_grid[cloudbox_limits[4]]-lon_grid[0] < llmin/lfac )
674  {
675  ostringstream os;
676  os << "Too small distance between cloudbox and lower end of\n"
677  << "longitude grid. This distance must here be "
678  << llmin/lfac << " degrees.";
679  throw runtime_error( os.str() );
680  }
681  if( lon_grid[n-1] - lon_grid[cloudbox_limits[5]] < llmin/lfac )
682  {
683  ostringstream os;
684  os << "Too small distance between cloudbox and upper end of\n"
685  << "longitude grid. This distance must here be "
686  << llmin/lfac << " degrees.";
687  throw runtime_error( os.str() );
688  }
689  }
690  }
691  }
692 
693  // If here, all OK
694  cloudbox_checked = 1;
695 }
696 
697 
698 
699 /* Workspace method: Doxygen documentation will be auto-generated */
700 void Massdensity_cleanup (//WS Output:
701  Tensor4& massdensity_field,
702  //WS Input:
703  const Numeric& massdensity_threshold,
704  const Verbosity&)
705 {
706  // Check that massdensity_fields contain realistic values
707  //(values smaller than massdensity_threshold will be set to 0)
708  for ( Index i=0; i<massdensity_field.nbooks(); i++ )
709  {
710  for ( Index j=0; j<massdensity_field.npages(); j++ )
711  {
712  for ( Index k=0; k<massdensity_field.nrows(); k++ )
713  {
714  for ( Index l=0; l<massdensity_field.ncols(); l++ )
715  {
716  if ( massdensity_field ( i,j,k,l ) < massdensity_threshold )
717  {
718  massdensity_field ( i,j,k,l ) = 0.0;
719  }
720  }
721  }
722  }
723  }
724 }
725 
726 
727 /* Workspace method: Doxygen documentation will be auto-generated */
728 void ParticleSpeciesInit (ArrayOfString& part_species,
729  const Verbosity&)
730 {
731  part_species.resize ( 0 );
732 }
733 
734 
735 /* Workspace method: Doxygen documentation will be auto-generated */
736 void ParticleSpeciesSet (// WS Generic Output:
737  ArrayOfString& part_species,
738  // Control Parameters:
739  const ArrayOfString& names,
740  const Verbosity& verbosity)
741 {
743 
744  part_species.resize ( names.nelem() );
745  //assign input strings to part_species
746  part_species = names;
747 
748  // Print list of particle settings to the most verbose output stream:
749  out3 << " Defined particle settings: ";
750  for ( Index i=0; i<part_species.nelem(); ++i )
751  {
752  out3 << "\n " << i << ": "<<part_species[i];
753 
754  }
755  out3 << '\n';
756 }
757 
758 /* Workspace method: Doxygen documentation will be auto-generated */
759 void ParticleTypeInit (//WS Output:
760  ArrayOfSingleScatteringData& scat_data_raw,
761  ArrayOfGriddedField3& pnd_field_raw,
762  const Verbosity&)
763 {
764  scat_data_raw.reserve ( 20 );
765  pnd_field_raw.reserve ( 20 );
766 }
767 
768 
769 /* Workspace method: Doxygen documentation will be auto-generated */
770 void ParticleTypeAddAll (//WS Output:
771  ArrayOfSingleScatteringData& scat_data_raw,
772  ArrayOfGriddedField3& pnd_field_raw,
773  // WS Input(needed for checking the datafiles):
774  const Index& atmosphere_dim,
775  const Vector& f_grid,
776  const Vector& p_grid,
777  const Vector& lat_grid,
778  const Vector& lon_grid,
779  const ArrayOfIndex& cloudbox_limits,
780  // Keywords:
781  const String& filename_scat_data,
782  const String& pnd_field_file,
783  const Verbosity& verbosity)
784 {
786 
787  //--- Check input ---------------------------------------------------------
788 
789  // Atmosphere
790  chk_if_in_range ( "atmosphere_dim", atmosphere_dim, 1, 3 );
791  chk_atm_grids ( atmosphere_dim, p_grid, lat_grid, lon_grid );
792 
793  // Cloudbox limits
794  if ( cloudbox_limits.nelem() != 2*atmosphere_dim )
795  throw runtime_error (
796  "*cloudbox_limits* is a vector which contains"
797  "the upper and lower\n"
798  "limit of the cloud for all atmospheric dimensions.\n"
799  "So its length must be 2 x *atmosphere_dim*" );
800  // Frequency grid
801  if ( f_grid.nelem() == 0 )
802  throw runtime_error ( "The frequency grid is empty." );
803  chk_if_increasing ( "f_grid", f_grid );
804 
805 
806  //--- Reading the data ---------------------------------------------------
807  ArrayOfString data_files;
808  xml_read_from_file ( filename_scat_data, data_files, verbosity );
809  scat_data_raw.resize ( data_files.nelem() );
810 
811  for ( Index i = 0; i<data_files.nelem(); i++ )
812  {
813 
814  out2 << " Read single scattering data\n";
815  xml_read_from_file ( data_files[i], scat_data_raw[i], verbosity );
816 
817  chk_single_scattering_data ( scat_data_raw[i],
818  data_files[i], f_grid,
819  verbosity );
820 
821  }
822 
823  out2 << " Read particle number density data \n";
824  xml_read_from_file ( pnd_field_file, pnd_field_raw, verbosity );
825 
826  chk_pnd_raw_data ( pnd_field_raw,
827  pnd_field_file, atmosphere_dim, p_grid, lat_grid,
828  lon_grid, cloudbox_limits, verbosity);
829 }
830 
831 /* Workspace method: Doxygen documentation will be auto-generated */
833  ArrayOfSingleScatteringData& scat_data_raw,
834  ArrayOfScatteringMetaData& scat_data_meta_array,
835  const Vector& f_grid,
836  // Keywords:
837  const String& filename_scat_data,
838  const String& filename_scat_meta_data,
839  const Verbosity& verbosity)
840 {
842 
843  //--- Reading the data ---------------------------------------------------
844  ArrayOfString data_files;
845  ArrayOfString meta_data_files;
846 
847  // single scattering data read to temporary ArrayOfSingleScatteringData
848  xml_read_from_file ( filename_scat_data, data_files, verbosity );
849  scat_data_raw.resize ( data_files.nelem() );
850 
851  for ( Index i = 0; i<data_files.nelem(); i++ )
852  {
853  out3 << " Read single scattering data\n";
854  xml_read_from_file ( data_files[i], scat_data_raw[i], verbosity );
855 
856  chk_single_scattering_data ( scat_data_raw[i],
857  data_files[i], f_grid,
858  verbosity );
859 
860  }
861 
862  // scattering meta data read to temporary ArrayOfScatteringMetaData
863  xml_read_from_file ( filename_scat_meta_data, meta_data_files, verbosity );
864  scat_data_meta_array.resize ( meta_data_files.nelem() );
865 
866  for ( Index i = 0; i<meta_data_files.nelem(); i++ )
867  {
868 
869  out3 << " Read scattering meta data\n";
870  xml_read_from_file ( meta_data_files[i], scat_data_meta_array[i], verbosity );
871 
872  chk_scattering_meta_data ( scat_data_meta_array[i],
873  meta_data_files[i], verbosity );
874 
875  }
876 
877  // check if arrays have same size
878  chk_scattering_data ( scat_data_raw,
879  scat_data_meta_array, verbosity );
880 
881 }
882 
883 
884 /* Workspace method: Doxygen documentation will be auto-generated */
885 void ScatteringParticlesSelect (//WS Output:
886  ArrayOfSingleScatteringData& scat_data_raw,
887  ArrayOfScatteringMetaData& scat_data_meta_array,
888  ArrayOfIndex& scat_data_nelem,
889  // WS Input:
890  const ArrayOfString& part_species,
891  const Verbosity& verbosity)
892 {
895  //--- Adjusting data to user specified input (part_species)-------------------
896 
897  String type;
898  Index intarr_total = 0;
899  ArrayOfIndex intarr;
900 
901  // make temporary copy
902  ArrayOfSingleScatteringData scat_data_raw_tmp = scat_data_raw;
903  ArrayOfScatteringMetaData scat_data_meta_array_tmp = scat_data_meta_array;
904 
905  scat_data_nelem.resize( part_species.nelem() );
906 
907  ArrayOfIndex selected;
908  selected.resize(scat_data_meta_array_tmp.nelem());
909  selected = 0;
910  // loop over array of part_species--------------------------------------------
911  for ( Index k=0; k<part_species.nelem(); k++ )
912  {
913 
914  String part_type;
915  Numeric sizemin;
916  Numeric sizemax;
917 
918  //split part_species string and copy values to parameter
919  parse_part_type( part_type, part_species[k]);
920  // set type according to *part_species* input
921  if ( part_type == "IWC" || part_type== "Snow" ) type = "Ice";
922  else if ( part_type== "LWC" || part_type == "Rain" ) type = "Water";
923 
924  //split part_species string and copy values to parameter
925  parse_part_size(sizemin, sizemax, part_species[k]);
926 
927 
928  // choosing the specified SingleScatteringData and ScatteringMetaData
929  for ( Index j=0; j<scat_data_meta_array_tmp.nelem(); j++ )
930  {
931  // check for particle phase type (e.g. "Ice", "Water",...)
932  if ( scat_data_meta_array_tmp[j].type == type )
933  {
934  // particle radius is calculated from particle volume given in
935  // scattering meta data
936  Numeric r_particle =
937  pow ( 3./4. * scat_data_meta_array_tmp[j].V * 1e18 / PI , 1./3. );
938 
939  // check if particle is in size range
940  // (sizemax < 0 results from wildcard usage and means all sizes on the
941  // upper end)
942  if ( r_particle >= sizemin &&
943  ( sizemax >= r_particle || sizemax < 0. ))
944  {
945  // fill ArrayOfIndex with indices of selected scattering data
946  intarr.push_back ( j );
947  }
948  selected[j] = 1;
949  out3 << "Selecting particle " << j+1 << "/" << scat_data_meta_array_tmp.nelem()
950  << " (" << scat_data_meta_array_tmp[j].type << ")\n";
951  }
952  }
953  // WSV scat_data_nelem gets the number of elements of scattering data
954  // connected to each selection String in *part_species*
955  scat_data_nelem[k] = intarr.nelem() - intarr_total;
956  intarr_total = intarr.nelem();
957  }
958  // check if array is empty
959  if ( !intarr.nelem() )
960  {
961  ostringstream os;
962  os << "The selection in " << part_species <<
963  " is NOT choosing any of the given Scattering Data.\n"
964  "--> Does the selection in *part_species* fit any of the "
965  "Single Scattering Data input? \n";
966  throw runtime_error ( os.str() );
967  }
968  // check if we ignored any smd
969  for ( Index j = 0; j<selected.nelem(); j++)
970  {
971  if (selected[j]==0)
972  {
973  out1 << "WARNING! Ignored SMD[" << j << "] (" << scat_data_meta_array_tmp[j].type << ")!\n";
974  }
975  }
976 
977 
978  // resize WSVs to size of intarr
979  scat_data_raw.resize ( intarr.nelem() );
980  scat_data_meta_array.resize ( intarr.nelem() );
981 
982  for ( Index j=0; j<intarr.nelem(); j++ )
983  {
984  //append to WSV Arrays
985  scat_data_meta_array[j] = scat_data_meta_array_tmp[intarr[j]] ;
986  scat_data_raw[j] = scat_data_raw_tmp[intarr[j]] ;
987  }
988 
989 
990 }
991 
992 
993 /* Workspace method: Doxygen documentation will be auto-generated */
994 void ParticleTypeAdd( //WS Output:
995  ArrayOfSingleScatteringData& scat_data_raw,
996  ArrayOfGriddedField3& pnd_field_raw,
997  // WS Input (needed for checking the datafiles):
998  const Index& atmosphere_dim,
999  const Vector& f_grid,
1000  const Vector& p_grid,
1001  const Vector& lat_grid,
1002  const Vector& lon_grid,
1003  const ArrayOfIndex& cloudbox_limits,
1004  // Keywords:
1005  const String& scat_data_file,
1006  const String& pnd_field_file,
1007  const Verbosity& verbosity)
1008 {
1009  CREATE_OUT2
1010 
1011  //--- Check input ---------------------------------------------------------
1012 
1013  // Atmosphere
1014  chk_if_in_range( "atmosphere_dim", atmosphere_dim, 1, 3 );
1015  chk_atm_grids( atmosphere_dim, p_grid, lat_grid, lon_grid );
1016 
1017  // Cloudbox limits
1018  if ( cloudbox_limits.nelem() != 2*atmosphere_dim )
1019  throw runtime_error(
1020  "*cloudbox_limits* is a vector which contains"
1021  "the upper and lower\n"
1022  "limit of the cloud for all atmospheric dimensions.\n"
1023  "So its length must be 2 x *atmosphere_dim*" );
1024  // Frequency grid
1025  if( f_grid.nelem() == 0 )
1026  throw runtime_error( "The frequency grid is empty." );
1027  chk_if_increasing( "f_grid", f_grid );
1028 
1029 
1030  //--- Reading the data ---------------------------------------------------
1031 
1032  // Append *scat_data_raw* and *pnd_field_raw* with empty Arrays of Tensors.
1033  SingleScatteringData scat_data;
1034  scat_data_raw.push_back(scat_data);
1035 
1036  GriddedField3 pnd_field_data;
1037  pnd_field_raw.push_back(pnd_field_data);
1038 
1039  out2 << " Read single scattering data\n";
1040  xml_read_from_file(scat_data_file, scat_data_raw[scat_data_raw.nelem()-1],
1041  verbosity);
1042 
1043  chk_single_scattering_data(scat_data_raw[scat_data_raw.nelem()-1],
1044  scat_data_file, f_grid, verbosity);
1045 
1046  out2 << " Read particle number density field\n";
1047  if (pnd_field_file.nelem() < 1)
1048  {
1049  CREATE_OUT1
1050  out1 << "Warning: No pnd_field_file specified. Ignored. \n";
1051  }
1052  else
1053  {
1054  xml_read_from_file(pnd_field_file, pnd_field_raw[pnd_field_raw.nelem()-1],
1055  verbosity);
1056 
1057  chk_pnd_data(pnd_field_raw[pnd_field_raw.nelem()-1],
1058  pnd_field_file, atmosphere_dim, p_grid, lat_grid,
1059  lon_grid, cloudbox_limits, verbosity);
1060  }
1061 }
1062 
1063 
1064 /* Workspace method: Doxygen documentation will be auto-generated */
1065 void pnd_fieldCalc(//WS Output:
1066  Tensor4& pnd_field,
1067  //WS Input
1068  const Vector& p_grid,
1069  const Vector& lat_grid,
1070  const Vector& lon_grid,
1071  const ArrayOfGriddedField3& pnd_field_raw,
1072  const Index& atmosphere_dim,
1073  const ArrayOfIndex& cloudbox_limits,
1074  const Verbosity&)
1075 {
1076  // Basic checks of input variables
1077  //
1078  // Particle number density data
1079  //
1080  if (pnd_field_raw.nelem() == 0)
1081  throw runtime_error(
1082  "No particle number density data given. Please\n"
1083  "use WSMs *ParticleTypeInit* and \n"
1084  "*ParticleTypeAdd(All)* for reading cloud particle\n"
1085  "data.\n"
1086  );
1087 
1088  chk_atm_grids( atmosphere_dim, p_grid, lat_grid, lon_grid );
1089  if ( cloudbox_limits.nelem()!= 2*atmosphere_dim)
1090  throw runtime_error(
1091  "*cloudbox_limits* is a vector which contains the"
1092  "upper and lower limit of the cloud for all "
1093  "atmospheric dimensions. So its dimension must"
1094  "be 2 x *atmosphere_dim*");
1095 
1096  // Check that pnd_field_raw has at least 2 grid-points in each dimension.
1097  // Otherwise, interpolation further down will fail with assertion.
1098 
1099  for (Index d = 0; d < atmosphere_dim; d++)
1100  {
1101  for (Index i = 0; i < pnd_field_raw.nelem(); i++)
1102  {
1103  if (pnd_field_raw[i].get_grid_size(d) < 2)
1104  {
1105  ostringstream os;
1106  os << "Error in pnd_field_raw data. ";
1107  os << "Dimension " << d << " (name: \"";
1108  os << pnd_field_raw[i].get_grid_name(d);
1109  os << "\") has only ";
1110  os << pnd_field_raw[i].get_grid_size(d);
1111  os << " element";
1112  os << ((pnd_field_raw[i].get_grid_size(d)==1) ? "" : "s");
1113  os << ". Must be at least 2.";
1114  throw runtime_error(os.str());
1115  }
1116  }
1117  }
1118  const Index Np_cloud = cloudbox_limits[1]-cloudbox_limits[0]+1;
1119 
1120  ConstVectorView p_grid_cloud = p_grid[Range(cloudbox_limits[0], Np_cloud)];
1121 
1122  // Check that no scatterers exist outside the cloudbox
1123  chk_pnd_field_raw_only_in_cloudbox(atmosphere_dim, pnd_field_raw,
1124  p_grid, lat_grid, lon_grid,
1125  cloudbox_limits);
1126 
1127  //==========================================================================
1128  if ( atmosphere_dim == 1)
1129  {
1130  //Resize variables
1131  pnd_field.resize(pnd_field_raw.nelem(), Np_cloud, 1, 1 );
1132 
1133  // Gridpositions:
1134  ArrayOfGridPos gp_p(Np_cloud);
1135 
1136  // Interpolate pnd_field.
1137  // Loop over the particle types:
1138  for (Index i = 0; i < pnd_field_raw.nelem(); ++ i)
1139  {
1140  // Calculate grid positions:
1141  p2gridpos( gp_p, pnd_field_raw[i].get_numeric_grid(GFIELD3_P_GRID),
1142  p_grid_cloud );
1143 
1144  // Interpolation weights:
1145  Matrix itw(Np_cloud, 2);
1146 
1147  interpweights( itw, gp_p);
1148  // Interpolate:
1149  interp( pnd_field(i,joker,0,0), itw,
1150  pnd_field_raw[i].data(joker,0,0), gp_p );
1151  }
1152  }
1153 
1154  else if(atmosphere_dim == 2)
1155  {
1156  const Index Nlat_cloud = cloudbox_limits[3]-cloudbox_limits[2]+1;
1157 
1158  ConstVectorView lat_grid_cloud =
1159  lat_grid[Range(cloudbox_limits[2],Nlat_cloud)];
1160 
1161  //Resize variables
1162  pnd_field.resize( pnd_field_raw.nelem(), Np_cloud, Nlat_cloud, 1 );
1163 
1164  // Gridpositions:
1165  ArrayOfGridPos gp_p(Np_cloud);
1166  ArrayOfGridPos gp_lat(Nlat_cloud);
1167 
1168  // Interpolate pnd_field.
1169  // Loop over the particle types:
1170  for (Index i = 0; i < pnd_field_raw.nelem(); ++ i)
1171  {
1172  // Calculate grid positions:
1173  p2gridpos( gp_p, pnd_field_raw[i].get_numeric_grid(GFIELD3_P_GRID),
1174  p_grid_cloud);
1175  gridpos( gp_lat, pnd_field_raw[i].get_numeric_grid(GFIELD3_LAT_GRID),
1176  lat_grid_cloud);
1177 
1178  // Interpolation weights:
1179  Tensor3 itw( Np_cloud, Nlat_cloud, 4 );
1180  interpweights( itw, gp_p, gp_lat );
1181 
1182  // Interpolate:
1183  interp( pnd_field(i,joker,joker,0), itw,
1184  pnd_field_raw[i].data(joker,joker,0), gp_p, gp_lat );
1185  }
1186  }
1187  else
1188  {
1189  const Index Nlat_cloud = cloudbox_limits[3]-cloudbox_limits[2]+1;
1190  const Index Nlon_cloud = cloudbox_limits[5]-cloudbox_limits[4]+1;
1191 
1192  ConstVectorView lat_grid_cloud =
1193  lat_grid[Range(cloudbox_limits[2],Nlat_cloud)];
1194  ConstVectorView lon_grid_cloud =
1195  lon_grid[Range(cloudbox_limits[4],Nlon_cloud)];
1196 
1197  //Resize variables
1198  pnd_field.resize( pnd_field_raw.nelem(), Np_cloud, Nlat_cloud,
1199  Nlon_cloud );
1200 
1201  // Gridpositions:
1202  ArrayOfGridPos gp_p(Np_cloud);
1203  ArrayOfGridPos gp_lat(Nlat_cloud);
1204  ArrayOfGridPos gp_lon(Nlon_cloud);
1205 
1206  // Interpolate pnd_field.
1207  // Loop over the particle types:
1208  for (Index i = 0; i < pnd_field_raw.nelem(); ++ i)
1209  {
1210  // Calculate grid positions:
1211  p2gridpos( gp_p, pnd_field_raw[i].get_numeric_grid(GFIELD3_P_GRID),
1212  p_grid_cloud);
1213  gridpos( gp_lat, pnd_field_raw[i].get_numeric_grid(GFIELD3_LAT_GRID),
1214  lat_grid_cloud);
1215  gridpos( gp_lon, pnd_field_raw[i].get_numeric_grid(GFIELD3_LON_GRID),
1216  lon_grid_cloud);
1217 
1218  // Interpolation weights:
1219  Tensor4 itw( Np_cloud, Nlat_cloud, Nlon_cloud, 8 );
1220  interpweights( itw, gp_p, gp_lat, gp_lon );
1221 
1222  // Interpolate:
1223  interp( pnd_field(i,joker,joker,joker), itw,
1224  pnd_field_raw[i].data, gp_p, gp_lat, gp_lon );
1225  }
1226  }
1227 }
1228 
1229 
1230 /* Workspace method: Doxygen documentation will be auto-generated */
1231 void pnd_fieldExpand1D(Tensor4& pnd_field,
1232  const Index& atmosphere_dim,
1233  const Index& cloudbox_checked,
1234  const Index& cloudbox_on,
1235  const ArrayOfIndex& cloudbox_limits,
1236  const Index& nzero,
1237  const Verbosity&)
1238 {
1239  if( !cloudbox_checked )
1240  throw runtime_error( "The cloudbox must be flagged to have passed a "
1241  "consistency check (cloudbox_checked=1)." );
1242 
1243  if( atmosphere_dim == 1 )
1244  { throw runtime_error( "No use in calling this method for 1D." ); }
1245  if( !cloudbox_on )
1246  { throw runtime_error(
1247  "No use in calling this method with cloudbox off." ); }
1248 
1249  if( nzero < 1 )
1250  { throw runtime_error( "The argument *nzero must be > 0." ); }
1251 
1252  // Sizes
1253  const Index npart = pnd_field.nbooks();
1254  const Index np = cloudbox_limits[1] - cloudbox_limits[0] + 1;
1255  const Index nlat = cloudbox_limits[3] - cloudbox_limits[2] + 1;
1256  Index nlon = 1;
1257  if( atmosphere_dim == 3 )
1258  { nlon = cloudbox_limits[5] - cloudbox_limits[4] + 1; }
1259 
1260  if( pnd_field.npages() != np || pnd_field.nrows() != 1 ||
1261  pnd_field.ncols() != 1 )
1262  { throw runtime_error( "The input *pnd_field* is either not 1D or does not "
1263  "match pressure size of cloudbox." );}
1264 
1265  // Temporary container
1266  Tensor4 pnd_temp = pnd_field;
1267 
1268  // Resize and fill
1269  pnd_field.resize( npart, np, nlat, nlon );
1270  pnd_field = 0;
1271  //
1272  for( Index ilon=nzero; ilon<nlon-nzero; ilon++ )
1273  {
1274  for( Index ilat=nzero; ilat<nlat-nzero; ilat++ )
1275  {
1276  for( Index ip=0; ip<np; ip++ )
1277  {
1278  for( Index is=0; is<npart; is++ )
1279  { pnd_field(is,ip,ilat,ilon) = pnd_temp(is,ip,0,0); }
1280  }
1281  }
1282  }
1283 }
1284 
1285 
1286 /* Workspace method: Doxygen documentation will be auto-generated */
1287 void pnd_fieldZero(//WS Output:
1288  Tensor4& pnd_field,
1289  ArrayOfSingleScatteringData& scat_data_raw,
1290  //WS Input:
1291  const Vector& p_grid,
1292  const Vector& lat_grid,
1293  const Vector& lon_grid,
1294  const Verbosity&)
1295 {
1296  // 3D atmosphere
1297  if (lat_grid.nelem()>1)
1298  {
1299  //Resize pnd_field and set it to 0:
1300  pnd_field.resize(1, p_grid.nelem(), lat_grid.nelem(), lon_grid.nelem());
1301  pnd_field = 0.;
1302  }
1303  else // 1D atmosphere
1304  {
1305  //Resize pnd_field and set it to 0:
1306  pnd_field.resize(1, p_grid.nelem(), 1, 1);
1307  pnd_field = 0.;
1308  }
1309 
1310  //Resize scat_data_raw and set it to 0:
1311  // Number iof particle types
1312  scat_data_raw.resize(1);
1313  scat_data_raw[0].ptype = PARTICLE_TYPE_MACROS_ISO;
1314  scat_data_raw[0].description = " ";
1315  // Grids which contain full ranges which one wants to calculate
1316  nlinspace(scat_data_raw[0].f_grid, 1e9, 3.848043e+13 , 5);
1317  nlinspace(scat_data_raw[0].T_grid, 0, 400, 5);
1318  nlinspace(scat_data_raw[0].za_grid, 0, 180, 5);
1319  nlinspace(scat_data_raw[0].aa_grid, 0, 360, 5);
1320  // Resize the data arrays
1321  scat_data_raw[0].pha_mat_data.resize(5,5,5,1,1,1,6);
1322  scat_data_raw[0].pha_mat_data = 0.;
1323  scat_data_raw[0].ext_mat_data.resize(5,5,1,1,1);
1324  scat_data_raw[0].ext_mat_data = 0.;
1325  scat_data_raw[0].abs_vec_data.resize(5,5,1,1,1);
1326  scat_data_raw[0].abs_vec_data = 0.;
1327 }
1328 
1329 
1330 /* Workspace method: Doxygen documentation will be auto-generated */
1331 void pnd_fieldSetup (//WS Output:
1332  Tensor4& pnd_field,
1333  //WS Input:
1334  const Index& atmosphere_dim,
1335  const Index& cloudbox_on,
1336  const ArrayOfIndex& cloudbox_limits,
1337  const Tensor4& massdensity_field,
1338  const Tensor3& t_field,
1339  const ArrayOfScatteringMetaData& scat_data_meta_array,
1340  const ArrayOfString& part_species,
1341  const ArrayOfIndex& scat_data_nelem,
1342  const Verbosity& verbosity)
1343 {
1344  // Cloudbox on/off?
1345  if ( !cloudbox_on )
1346  {
1347  /* Must initialise pnd_field anyway; but empty */
1348  pnd_field.resize(0, 0, 0, 0);
1349  return;
1350  }
1351 
1352  // ------- set pnd_field boundaries to cloudbox boundaries -------------------
1353  //initialize pnd_field boundaries
1354  Index p_cbstart = 0;
1355  Index p_cbend = 1;
1356  Index lat_cbstart = 0;
1357  Index lat_cbend = 1;
1358  Index lon_cbstart = 0;
1359  Index lon_cbend = 1;
1360  //pressure
1361  p_cbstart = cloudbox_limits[0];
1362  p_cbend = cloudbox_limits[1]+1;
1363 
1364  //latitude
1365  if ( atmosphere_dim >= 2 )
1366  {
1367  lat_cbstart = cloudbox_limits[2];
1368  lat_cbend = cloudbox_limits[3]+1;
1369  }
1370  //longitude
1371  if ( atmosphere_dim == 3 )
1372  {
1373  lon_cbstart = cloudbox_limits[4];
1374  lon_cbend = cloudbox_limits[5]+1;
1375 
1376  }
1377 
1378  /* Do some checks. Not foolproof, but catches at least some. */
1379 
1380  if ((p_cbend > massdensity_field.npages()) ||
1381  (p_cbend > t_field.npages()) ||
1382  (lat_cbend > massdensity_field.nrows()) ||
1383  (lat_cbend > t_field.nrows()) ||
1384  (lon_cbend > massdensity_field.ncols()) ||
1385  (lon_cbend > t_field.ncols()))
1386  {
1387  ostringstream os;
1388  os << "Cloudbox out of bounds compared to fields. "
1389  << "Upper limits: (p, lat, lon): "
1390  << "(" << p_cbend << ", " << lat_cbend << ", " << lon_cbend << "). "
1391  << "*massdensity_field*: "
1392  << "(" << massdensity_field.npages() << ", "
1393  << massdensity_field.nrows() << ", "
1394  << massdensity_field.ncols() << "). "
1395  << "*t_field*: "
1396  << "(" << t_field.npages() << ", "
1397  << t_field.nrows() << ", "
1398  << t_field.ncols() << ").";
1399  throw runtime_error(os.str());
1400  }
1401 
1402  //resize pnd_field to required atmospheric dimension and scatt particles
1403  pnd_field.resize ( scat_data_meta_array.nelem(), p_cbend-p_cbstart,
1404  lat_cbend-lat_cbstart, lon_cbend-lon_cbstart );
1405  Index scat_data_start = 0;
1406  ArrayOfIndex intarr;
1407 
1408  //-------- Start pnd_field calculations---------------------------------------
1409 
1410  // loop over nelem of part_species
1411  for ( Index k=0; k<part_species.nelem(); k++ )
1412  {
1413 
1414  String psd_param;
1415 
1416  //split String and copy to ArrayOfString
1417  parse_psd_param( psd_param, part_species[k]);
1418 
1419  // initialize control parameters
1420  Vector vol_unsorted ( scat_data_nelem[k], 0.0 );
1421  Vector d_max_unsorted (scat_data_nelem[k], 0.0);
1422  Vector vol ( scat_data_nelem[k], 0.0 );
1423  Vector dm ( scat_data_nelem[k], 0.0 );
1424  Vector r ( scat_data_nelem[k], 0.0 );
1425  Vector rho ( scat_data_nelem[k], 0.0 );
1426  Vector pnd ( scat_data_nelem[k], 0.0 );
1427  //Vector pnd2 ( scat_data_nelem[k], 0.0 ); //temporary
1428  Vector dN ( scat_data_nelem[k], 0.0 );
1429  //Vector dN2 ( scat_data_nelem[k], 0.0 ); //temporary
1430  //Vector dlwc ( scat_data_nelem[k], 0.0 ); //temporary
1431 
1432 
1433  //---- start pnd_field calculations for MH97 -------------------------------
1434  if ( psd_param == "MH97" )
1435  {
1436 
1437  for ( Index i=0; i < scat_data_nelem[k]; i++ )
1438  {
1439  //m^3
1440  vol_unsorted[i] = ( scat_data_meta_array[i+scat_data_start].V );
1441  }
1442  get_sorted_indexes(intarr, vol_unsorted);
1443  //cout<<"intarr\t"<<intarr<<endl;
1444 
1445  // NOTE: the order of scattering particle profiles in *massdensity_field*
1446  // is HARD WIRED!
1447  // extract IWC_field and convert from kg/m^3 to g/m^3
1448  Tensor3 IWC_field = massdensity_field ( 1, joker, joker, joker );
1449  IWC_field*=1000; //IWC [g/m^3]
1450  //out0<<"\n"<<IWC_field<<"\n";
1451 
1452  // extract scattering meta data
1453  for ( Index i=0; i< scat_data_nelem[k]; i++ )
1454  {
1455  vol[i] = ( scat_data_meta_array[intarr[i]+scat_data_start].V ); //m^3
1456  // calculate melted diameter from volume [m]
1457  dm[i] = pow (
1458  ( (6*scat_data_meta_array[intarr[i]+scat_data_start].V) /PI ),
1459  ( 1./3. ) );
1460  // get density from meta data [g/m^3]
1461  rho[i] = scat_data_meta_array[intarr[i]+scat_data_start].density * 1000;
1462 
1463  //check for correct particle phase
1464  if ( scat_data_meta_array[intarr[i]+scat_data_start].type != "Ice" )
1465  {
1466  throw runtime_error ( "The particle phase is unequal 'Ice'.\n"
1467  "MH97 can only be applied to ice particles.\n"
1468  "Check ScatteringMetaData!" );
1469  }
1470  }
1471 
1472 
1473  // itertation over all atm. levels
1474  for ( Index p=p_cbstart; p<p_cbend; p++ )
1475  {
1476  for ( Index lat=lat_cbstart; lat<lat_cbend; lat++ )
1477  {
1478  for ( Index lon=lon_cbstart; lon<lon_cbend; lon++ )
1479  {
1480  // iteration over all given size bins
1481  for ( Index i=0; i<dm.nelem(); i++ )
1482  {
1483  // calculate particle size distribution with MH97
1484  // [# m^-3 m^-1]
1485  dN[i] = IWCtopnd_MH97 ( IWC_field ( p, lat, lon ), dm[i],
1486  t_field ( p, lat, lon ), rho[i] );
1487  //dN2[i] = dN[i] * vol[i] * rho[i];
1488  }
1489  //out0<<"level: "<<p<<"\n"<<dN<<"\n";
1490 
1491  // scale pnds by bin width
1492  if (dm.nelem() > 1)
1493  {
1494  scale_pnd( pnd, dm, dN );
1495  } else
1496  {
1497  pnd = dN;
1498  }
1499 
1500 
1501  //out0<<"level: "<<p<<"\n"<<pnd<<"\n"<<"IWC: "<<IWC_field ( p, lat, lon )<<"\n"<<"T: "<<t_field ( p, lat, lon )<<"\n";
1502  // calculate error of pnd sum and real XWC
1503  chk_pndsum ( pnd, IWC_field ( p,lat,lon ), vol, rho, p, lat, lon, verbosity );
1504  //chk_pndsum2 (pnd2, IWC_field ( p,lat,lon ));
1505 
1506  // writing pnd vector to wsv pnd_field
1507  for ( Index i = 0; i< scat_data_nelem[k]; i++ )
1508  {
1509  pnd_field ( intarr[i]+scat_data_start, p-p_cbstart,
1510  lat-lat_cbstart, lon-lon_cbstart ) = pnd[i];
1511  }
1512 
1513  }
1514  }
1515  }
1516 
1517  }
1518 
1519  //---- start pnd_field calculations for H11 ----------------------------
1520  else if ( psd_param == "H11" )
1521  {
1522  String part_type;
1523  Tensor3 X_field;
1524 
1525  for ( Index i=0; i < scat_data_nelem[k]; i++ )
1526  {
1527  //m
1528  d_max_unsorted[i] = ( scat_data_meta_array[i+scat_data_start].d_max );
1529  }
1530  get_sorted_indexes(intarr, d_max_unsorted);
1531 
1532  //get particle type to decide if H11 gets apllied on 'IWC' profile or 'Snow' profile
1533  parse_part_type( part_type, part_species[k]);
1534 
1535  if ( part_type == "IWC" )
1536  {
1537  // NOTE: the order of scattering particle profiles in *massdensity_field*
1538  // is HARD WIRED!
1539  // extract IWC and convert from kg/m^3 to g/m^3
1540  X_field = massdensity_field ( 1, joker, joker, joker );
1541  X_field *= 1000; //IWC [g/m^3]
1542  }
1543  else if ( part_type == "Snow" )
1544  {
1545  // NOTE: the order of scattering particle profiles in *massdensity_field*
1546  // is HARD WIRED!
1547  // extract Snow rate and convert from kg/(m2*s) to g/(m2*s)
1548  X_field = massdensity_field ( 3, joker, joker, joker );
1549  X_field *= 1000; //Snow [g/(m2*s)]
1550  }
1551  // extract scattering meta data
1552  for ( Index i=0; i< scat_data_nelem[k]; i++ )
1553  {
1554  vol[i]= scat_data_meta_array[intarr[i]+scat_data_start].V; //[m^3]
1555 
1556  // get maximum diameter from meta data [m]
1557  dm[i] = scat_data_meta_array[intarr[i]+scat_data_start].d_max;
1558 
1559  // get density from meta data [g/m^3]
1560  rho[i] = scat_data_meta_array[intarr[i]+scat_data_start].density * 1000;
1561 
1562  //check for correct particle phase
1563  if ( scat_data_meta_array[intarr[i]+scat_data_start].type != "Ice" )
1564  {
1565  throw runtime_error (
1566  "The particle phase is unequal 'Ice'.\n"
1567  "H11 can only be applied to ice/snow particles.\n"
1568  "Check ScatteringMetaData!" );
1569  }
1570  }
1571  // itertation over all atm. levels
1572  for ( Index p=p_cbstart; p<p_cbend; p++ )
1573  {
1574  for ( Index lat=lat_cbstart; lat<lat_cbend; lat++ )
1575  {
1576  for ( Index lon=lon_cbstart; lon<lon_cbend; lon++ )
1577  {
1578  // iteration over all given size bins
1579  for ( Index i=0; i<dm.nelem(); i++ ) //loop over number of particles
1580  {
1581  // calculate particle size distribution for H11
1582  // [# m^-3 m^-1]
1583  dN[i] = psd_H11 ( X_field ( p, lat, lon), dm[i], t_field ( p, lat, lon ) );
1584  }
1585  // scale pnds by scale width
1586  if (dm.nelem() > 1)
1587  {
1588  scale_pnd( pnd, dm, dN ); //[# m^-3]
1589  } else
1590  {
1591  pnd = dN;
1592  }
1593 
1594  //cout<<"\nlevel: "<<p<<"\n"<<"X: "<<X_field ( p, lat, lon )<<"\n"<<"T: "<<t_field ( p, lat, lon )<<"\n"<< dN<<"\n"<<pnd<<"\n";
1595  // scale H11 distribution (which is independent of Ice or
1596  // Snow massdensity) to current massdensity.
1597  // Output pnd: still in [# m^-3]
1598  scale_H11 ( pnd, X_field ( p,lat,lon ), vol, rho );
1599 
1600  // calculate error of pnd sum and real XWC
1601  chk_pndsum ( pnd, X_field ( p,lat,lon ), vol, rho, p, lat, lon, verbosity );
1602 
1603  //cout<<pnd<<"\n";
1604 
1605  // writing pnd vector to wsv pnd_field
1606  for ( Index i =0; i< scat_data_nelem[k]; i++ )
1607  {
1608  pnd_field ( intarr[i]+scat_data_start, p-p_cbstart,
1609  lat-lat_cbstart, lon-lon_cbstart ) = pnd[i];
1610  //dlwc[q] = pnd2[q]*vol[q]*rho[q];
1611  }
1612  }
1613  }
1614  }
1615 
1616  }
1617 
1618  // ---- start pnd_field calculations for liquid ----------------------------
1619  else if ( psd_param == "liquid" )
1620  {
1621  for ( Index i=0; i < scat_data_nelem[k]; i++ )
1622  {
1623  //m^3
1624  vol_unsorted[i] = ( scat_data_meta_array[i+scat_data_start].V );
1625  }
1626  get_sorted_indexes(intarr, vol_unsorted);
1627  //cout<<"intarr\t"<<intarr<<endl;
1628 
1629  // NOTE: the order of scattering particle profiles in *massdensity_field*
1630  // is HARD WIRED!
1631  // extract LWC_field and convert from kg/m^3 to g/m^3
1632  Tensor3 LWC_field = massdensity_field ( 0, joker, joker, joker );
1633  LWC_field *= 1000; //LWC [g/m^3]
1634 
1635  // extract scattering meta data
1636  for ( Index i=0; i< scat_data_nelem[k]; i++ )
1637  {
1638  vol[i]= scat_data_meta_array[intarr[i]+scat_data_start].V; //m^3
1639  // calculate diameter from volume [m]
1640  dm[i] = pow (
1641  ( 6*scat_data_meta_array[intarr[i]+scat_data_start].V/PI ),
1642  ( 1./3. ) );
1643  // diameter to radius
1644  r[i] = dm[i]/2; // [m]
1645  // get density from meta data [g/m^3]
1646  rho[i] = scat_data_meta_array[intarr[i]+scat_data_start].density * 1000;
1647 
1648  //check for correct particle phase
1649  if ( scat_data_meta_array[intarr[i]+scat_data_start].type != "Water" )
1650  {
1651  throw runtime_error (
1652  "The particle phase is unequal 'Water'.\n"
1653  "All particles must be of liquid phase to apply this PSD.\n"
1654  "Check ScatteringMetaData!" );
1655  }
1656  }
1657  //cout<<"\nr:\t"<<r<<endl;
1658  // itertation over all atm. levels
1659  for ( Index p=p_cbstart; p<p_cbend; p++ )
1660  {
1661  for ( Index lat=lat_cbstart; lat<lat_cbend; lat++ )
1662  {
1663  for ( Index lon=lon_cbstart; lon<lon_cbend; lon++ )
1664  {
1665  // iteration over all given size bins
1666  for ( Index i=0; i<r.nelem(); i++ ) //loop over number of particles
1667  {
1668  // calculate particle size distribution for liquid
1669  // [# m^-3 m^-1]
1670  dN[i] = LWCtopnd ( LWC_field ( p,lat,lon ), r[i] );
1671  //dN2[i] = LWCtopnd2 ( r[i] ); // [# m^-3 m^-1]
1672  //dN2[i] = dN[i] * vol[i] * rho[i];
1673  }
1674  //dlwc *= LWC_field(p, lat, lon)/dlwc.sum();
1675  //out0<<"\n"<<dN<<"\n"<<dN2<<"\n";
1676 
1677  // scale pnds by scale width
1678  if (r.nelem() > 1)
1679  {
1680  scale_pnd( pnd, r, dN ); //[# m^-3]
1681  } else
1682  {
1683  pnd = dN;
1684  }
1685  //scale_pnd( pnd2, r, dN2 );
1686  //trapezoid_integrate ( pnd2, r, dN2 );//[# m^-3]
1687  //out0<<"\n"<<"HIER!"<<"\n"<<pnd<<"\n"<<pnd2<<"\n";
1688 
1689  // calculate error of pnd sum and real XWC
1690  chk_pndsum ( pnd, LWC_field ( p,lat,lon ), vol, rho, p, lat, lon, verbosity );
1691  //chk_pndsum2 (pnd2, LWC_field ( p,lat,lon ));
1692  //chk_pndsum (pnd, testiwc[p], vol, rho);
1693 
1694 
1695  // writing pnd vector to wsv pnd_field
1696  for ( Index i =0; i< scat_data_nelem[k]; i++ )
1697  {
1698  pnd_field ( intarr[i]+scat_data_start, p-p_cbstart,
1699  lat-lat_cbstart, lon-lon_cbstart ) = pnd[i];
1700  //dlwc[q] = pnd2[q]*vol[q]*rho[q];
1701  }
1702 
1703  //out0<<"\n"<<LWC_field(p, lat, lon)<<"\n"<< dlwc.sum()<<"\n";
1704 
1705  //pnd2 *= ( LWC_field ( p, lat, lon ) /dlwc.sum() );
1706 
1707  //out0<<"\n"<<pnd<<"\n"<<pnd2<<"\n";
1708 
1709  }
1710  }
1711  }
1712  }
1713 
1714  // alter starting index of current scattering data array to starting index
1715  // of next iteration step
1716  scat_data_start = scat_data_start + scat_data_nelem[k];
1717 
1718  }
1719 }
1720 
1721 
Matrix
The Matrix class.
Definition: matpackI.h:767
gridded_fields.h
Implementation of gridded fields.
pnd_fieldSetup
void pnd_fieldSetup(Tensor4 &pnd_field, const Index &atmosphere_dim, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Tensor4 &massdensity_field, const Tensor3 &t_field, const ArrayOfScatteringMetaData &scat_data_meta_array, const ArrayOfString &part_species, const ArrayOfIndex &scat_data_nelem, const Verbosity &verbosity)
WORKSPACE METHOD: pnd_fieldSetup.
Definition: m_cloudbox.cc:1331
Tensor4::resize
void resize(Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackIV.cc:1404
chk_pndsum
void chk_pndsum(Vector &pnd, const Numeric xwc, const Vector &vol, const Vector &density, const Index &p, const Index &lat, const Index &lon, const Verbosity &verbosity)
Definition: cloudbox.cc:1080
psd_H11
Numeric psd_H11(const Numeric xwc, const Numeric d, const Numeric t)
Definition: cloudbox.cc:927
scale_pnd
void scale_pnd(Vector &w, const Vector &x, const Vector &y)
Definition: cloudbox.cc:1040
ScatteringParticleTypeAndMetaRead
void ScatteringParticleTypeAndMetaRead(ArrayOfSingleScatteringData &scat_data_raw, ArrayOfScatteringMetaData &scat_data_meta_array, const Vector &f_grid, const String &filename_scat_data, const String &filename_scat_meta_data, const Verbosity &verbosity)
WORKSPACE METHOD: ScatteringParticleTypeAndMetaRead.
Definition: m_cloudbox.cc:832
auto_md.h
chk_if_increasing
void chk_if_increasing(const String &x_name, const ArrayOfIndex &x)
chk_if_increasing
Definition: check_input.cc:126
gridpos
void gridpos(ArrayOfGridPos &gp, ConstVectorView old_grid, ConstVectorView new_grid, const Numeric &extpolfac)
Set up a grid position Array.
Definition: interpolation.cc:167
Tensor3
The Tensor3 class.
Definition: matpackIII.h:340
chk_scattering_meta_data
void chk_scattering_meta_data(const ScatteringMetaData &scat_data_meta, const String &scat_data_meta_file, const Verbosity &verbosity)
Check scattering data meta files.
Definition: cloudbox.cc:454
ParticleTypeInit
void ParticleTypeInit(ArrayOfSingleScatteringData &scat_data_raw, ArrayOfGriddedField3 &pnd_field_raw, const Verbosity &)
WORKSPACE METHOD: ParticleTypeInit.
Definition: m_cloudbox.cc:759
cloudboxOff
void cloudboxOff(Index &cloudbox_on, ArrayOfIndex &cloudbox_limits, Agenda &iy_cloudbox_agenda, const Verbosity &)
WORKSPACE METHOD: cloudboxOff.
Definition: m_cloudbox.cc:77
chk_if_bool
void chk_if_bool(const String &x_name, const Index &x)
chk_if_bool
Definition: check_input.cc:67
ParticleSpeciesSet
void ParticleSpeciesSet(ArrayOfString &part_species, const ArrayOfString &names, const Verbosity &verbosity)
WORKSPACE METHOD: ParticleSpeciesSet.
Definition: m_cloudbox.cc:736
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:350
p2gridpos
void p2gridpos(ArrayOfGridPos &gp, ConstVectorView old_pgrid, ConstVectorView new_pgrid, const Numeric &extpolfac)
p2gridpos
Definition: special_interp.cc:810
joker
const Joker joker
interpolation.h
Header file for interpolation.cc.
GFIELD3_LAT_GRID
const Index GFIELD3_LAT_GRID
interp
Numeric interp(ConstVectorView itw, ConstVectorView a, const GridPos &tc)
Red 1D Interpolate.
Definition: interpolation.cc:1000
ScatteringParticlesSelect
void ScatteringParticlesSelect(ArrayOfSingleScatteringData &scat_data_raw, ArrayOfScatteringMetaData &scat_data_meta_array, ArrayOfIndex &scat_data_nelem, const ArrayOfString &part_species, const Verbosity &verbosity)
WORKSPACE METHOD: ScatteringParticlesSelect.
Definition: m_cloudbox.cc:885
PARTICLE_TYPE_MACROS_ISO
@ PARTICLE_TYPE_MACROS_ISO
Definition: optproperties.h:57
pnd_fieldExpand1D
void pnd_fieldExpand1D(Tensor4 &pnd_field, const Index &atmosphere_dim, const Index &cloudbox_checked, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Index &nzero, const Verbosity &)
WORKSPACE METHOD: pnd_fieldExpand1D.
Definition: m_cloudbox.cc:1231
DEG2RAD
const Numeric DEG2RAD
CREATE_OUT2
#define CREATE_OUT2
Definition: messages.h:207
Tensor4
The Tensor4 class.
Definition: matpackIV.h:375
array.h
This file contains the definition of Array.
pnd_fieldZero
void pnd_fieldZero(Tensor4 &pnd_field, ArrayOfSingleScatteringData &scat_data_raw, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Verbosity &)
WORKSPACE METHOD: pnd_fieldZero.
Definition: m_cloudbox.cc:1287
CREATE_OUT1
#define CREATE_OUT1
Definition: messages.h:206
Agenda
The Agenda class.
Definition: agenda_class.h:44
cloudboxSetAutomatically
void cloudboxSetAutomatically(Index &cloudbox_on, ArrayOfIndex &cloudbox_limits, const Index &atmosphere_dim, const ArrayOfString &part_species, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor4 &massdensity_field, const Numeric &cloudbox_margin, const Verbosity &verbosity)
WORKSPACE METHOD: cloudboxSetAutomatically.
Definition: m_cloudbox.cc:92
GriddedField3
Definition: gridded_fields.h:274
SingleScatteringData
Structure which describes the single scattering properties of a.
Definition: optproperties.h:74
cloudbox_checkedCalc
void cloudbox_checkedCalc(Index &cloudbox_checked, const Index &basics_checked, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &wind_u_field, const Tensor3 &wind_v_field, const Tensor3 &wind_w_field, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Verbosity &)
WORKSPACE METHOD: cloudbox_checkedCalc.
Definition: m_cloudbox.cc:555
GFIELD3_LON_GRID
const Index GFIELD3_LON_GRID
ConstTensor3View::npages
Index npages() const
Returns the number of pages.
Definition: matpackIII.h:151
chk_scattering_data
void chk_scattering_data(const ArrayOfSingleScatteringData &scat_data_raw, const ArrayOfScatteringMetaData &scat_data_meta_array, const Verbosity &)
Check scattering data general.
Definition: cloudbox.cc:427
DEBUG_ONLY
#define DEBUG_ONLY(...)
Definition: arts.h:146
parse_part_size
void parse_part_size(Numeric &sizemin, Numeric &sizemax, const String &part_string)
Definition: cloudbox.cc:1270
Array< Index >
get_sorted_indexes
void get_sorted_indexes(ArrayOfIndex &sorted, const T &data)
get_sorted_indexes
Definition: sorting.h:79
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:850
CREATE_OUT3
#define CREATE_OUT3
Definition: messages.h:208
CREATE_OUT0
#define CREATE_OUT0
Definition: messages.h:205
PI
const Numeric PI
scale_H11
void scale_H11(Vector &pnd, const Numeric xwc, const Vector &density, const Vector &vol)
Definition: cloudbox.cc:1158
messages.h
Declarations having to do with the four output streams.
LWCtopnd
Numeric LWCtopnd(const Numeric lwc, const Numeric r)
Definition: cloudbox.cc:983
my_basic_string
The implementation for String, the ARTS string class.
Definition: mystring.h:62
abs
#define abs(x)
Definition: continua.cc:13094
ConstVectorView::nelem
Index nelem() const
Returns the number of elements.
Definition: matpackI.cc:175
physics_funcs.h
ConstTensor4View::ncols
Index ncols() const
Returns the number of columns.
Definition: matpackIV.cc:78
optproperties.h
Scattering database structure and functions.
Numeric
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
Verbosity
Definition: messages.h:50
parse_psd_param
void parse_psd_param(String &psd_param, const String &part_string)
Definition: cloudbox.cc:1240
ConstTensor4View::npages
Index npages() const
Returns the number of pages.
Definition: matpackIV.cc:66
chk_pnd_data
void chk_pnd_data(const GriddedField3 &pnd_field_raw, const String &pnd_field_file, const Index &atmosphere_dim, ConstVectorView p_grid, ConstVectorView lat_grid, ConstVectorView lon_grid, const ArrayOfIndex &cloudbox_limits, const Verbosity &verbosity)
Check particle number density files.
Definition: cloudbox.cc:283
ConstTensor4View::nbooks
Index nbooks() const
Returns the number of books.
Definition: matpackIV.cc:60
math_funcs.h
ParticleTypeAdd
void ParticleTypeAdd(ArrayOfSingleScatteringData &scat_data_raw, ArrayOfGriddedField3 &pnd_field_raw, const Index &atmosphere_dim, const Vector &f_grid, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const ArrayOfIndex &cloudbox_limits, const String &scat_data_file, const String &pnd_field_file, const Verbosity &verbosity)
WORKSPACE METHOD: ParticleTypeAdd.
Definition: m_cloudbox.cc:994
IWCtopnd_MH97
Numeric IWCtopnd_MH97(const Numeric iwc, const Numeric dm, const Numeric t, const Numeric density)
Definition: cloudbox.cc:828
ConstTensor3View::nrows
Index nrows() const
Returns the number of rows.
Definition: matpackIII.h:154
ParticleTypeAddAll
void ParticleTypeAddAll(ArrayOfSingleScatteringData &scat_data_raw, ArrayOfGriddedField3 &pnd_field_raw, const Index &atmosphere_dim, const Vector &f_grid, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const ArrayOfIndex &cloudbox_limits, const String &filename_scat_data, const String &pnd_field_file, const Verbosity &verbosity)
WORKSPACE METHOD: ParticleTypeAddAll.
Definition: m_cloudbox.cc:770
max
#define max(a, b)
Definition: continua.cc:13097
my_basic_string::nelem
Index nelem() const
Number of elements.
Definition: mystring.h:242
nlinspace
void nlinspace(Vector &x, const Numeric start, const Numeric stop, const Index n)
nlinspace
Definition: math_funcs.cc:261
cloudbox.h
Internal cloudbox functions.
ConstTensor4View::nrows
Index nrows() const
Returns the number of rows.
Definition: matpackIV.cc:72
Range
The range class.
Definition: matpackI.h:148
pnd_fieldCalc
void pnd_fieldCalc(Tensor4 &pnd_field, 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 Verbosity &)
WORKSPACE METHOD: pnd_fieldCalc.
Definition: m_cloudbox.cc:1065
logic.h
Header file for logic.cc.
ConstTensor3View::ncols
Index ncols() const
Returns the number of columns.
Definition: matpackIII.h:157
cloudboxSetManuallyAltitude
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:453
chk_pnd_field_raw_only_in_cloudbox
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: check_input.cc:1211
chk_atm_grids
void chk_atm_grids(const Index &dim, ConstVectorView p_grid, ConstVectorView lat_grid, ConstVectorView lon_grid)
chk_atm_grids
Definition: check_input.cc:603
min
#define min(a, b)
Definition: continua.cc:13096
rte.h
Declaration of functions in rte.cc.
chk_single_scattering_data
void chk_single_scattering_data(const SingleScatteringData &scat_data_raw, const String &scat_data_file, ConstVectorView f_grid, const Verbosity &verbosity)
Check single scattering data files.
Definition: cloudbox.cc:486
GFIELD3_P_GRID
const Index GFIELD3_P_GRID
Agenda::set_name
void set_name(const String &nname)
Set agenda name.
Definition: agenda_class.cc:602
special_interp.h
Header file for special_interp.cc.
chk_pnd_raw_data
void chk_pnd_raw_data(const ArrayOfGriddedField3 &pnd_field_raw, const String &pnd_field_file, const Index &atmosphere_dim, ConstVectorView p_grid, ConstVectorView lat_grid, ConstVectorView lon_grid, const ArrayOfIndex &cloudbox_limits, const Verbosity &verbosity)
Check particle number density files (pnd_field_raw)
Definition: cloudbox.cc:394
Index
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
chk_if_in_range
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:95
Massdensity_cleanup
void Massdensity_cleanup(Tensor4 &massdensity_field, const Numeric &massdensity_threshold, const Verbosity &)
WORKSPACE METHOD: Massdensity_cleanup.
Definition: m_cloudbox.cc:700
check_input.h
ParticleSpeciesInit
void ParticleSpeciesInit(ArrayOfString &part_species, const Verbosity &)
WORKSPACE METHOD: ParticleSpeciesInit.
Definition: m_cloudbox.cc:728
Vector
The Vector class.
Definition: matpackI.h:555
ConstVectorView
A constant view of a Vector.
Definition: matpackI.h:300
sorting.h
Contains sorting routines.
Array::nelem
Index nelem() const
Number of elements.
Definition: array.h:172
interpweights
void interpweights(VectorView itw, const GridPos &tc)
Red 1D interpolation weights.
Definition: interpolation.cc:756
chk_massdensity_field
void chk_massdensity_field(bool &empty_flag, const Index &dim, const Tensor3 &massdensity, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid)
Check whether hydromet grid size is equal to atmospheric grid size and if hydromet profile is zero (n...
Definition: cloudbox.cc:62
parse_part_type
void parse_part_type(String &part_type, const String &part_string)
Definition: cloudbox.cc:1205
arts.h
The global header file for ARTS.
xml_io.h
This file contains basic functions to handle XML data files.
barometric_heightformula
Numeric barometric_heightformula(const Numeric &p, const Numeric &dh)
Definition: cloudbox.cc:786