OB.DAAC Logo
NASA Logo
Ocean Color Science Software

ocssw V2022
TmProcess.cpp
Go to the documentation of this file.
1 /*******************************************************************************
2  *
3  * NAME: TmProcess.cpp
4  *
5  * DESCRIPTION: .
6  *
7  * Created on: November 3, 2016
8  * Author: Sam Anderson
9  *
10  * Modified:
11  *
12  *******************************************************************************/
13 
14 #include <new> // nothrow
15 #include <algorithm> // sort
16 #include <iostream> // cout
17 #include <functional> // bind
18 #include <vector>
19 
20 #include <libgen.h>
21 
22 #include <TmProcess.h>
23 #include <TmParamsReader.h>
24 
25 using namespace std;
26 
27 extern "C" {
28 int calcrand_(double* RAT, double* LAM, double* MRR, double* MRI, double* EPS,
29  int* NP, double* DDELT, int* NDGS, int* NPNAX, double* AXMAX, double* B,
30  double* GAM, int* NDISTR, int* NKMAX, int* NPNA, int* NCOEFF,
31  double* REFF, double* VEFF, double* CEXTIN, double* CSCAT, double* WALB,
32  double* ASYMM, double* ALPHA, double* BETA, double* F, double* SZ);
33 }
34 
35 const double pace_wavelengths[NUM_PACE_BANDS] = {
36  310.0, 315.0, 320.0, 325.0, 330.0, 335.0, 340.0, 345.0,
37  350.0, 355.0, 360.0, 365.0, 370.0, 375.0, 380.0, 385.0,
38  390.0, 395.0, 400.0, 405.0, 410.0, 415.0, 420.0, 425.0,
39  430.0, 435.0, 440.0, 445.0, 450.0, 455.0, 460.0, 465.0,
40  470.0, 475.0, 480.0, 485.0, 490.0, 495.0, 500.0, 505.0,
41  510.0, 515.0, 520.0, 525.0, 530.0, 535.0, 540.0, 545.0,
42  550.0, 555.0, 560.0, 565.0, 570.0, 575.0, 580.0, 585.0,
43  590.0, 595.0, 600.0, 605.0, 610.0, 615.0, 620.0, 625.0,
44  630.0, 635.0, 640.0, 645.0, 650.0, 655.0, 660.0, 665.0,
45  670.0, 675.0, 680.0, 685.0, 690.0, 695.0, 700.0, 705.0,
46  710.0, 715.0, 720.0, 725.0, 730.0, 735.0, 740.0, 745.0,
47  750.0, 755.0, 760.0, 765.0, 770.0, 775.0, 780.0, 785.0,
48  790.0, 795.0, 800.0, 805.0, 810.0, 815.0, 820.0, 825.0,
49  830.0, 835.0, 840.0, 845.0, 850.0, 855.0, 860.0, 865.0,
50  870.0, 875.0, 880.0, 885.0, 890.0, 895.0, 940.0, 1040.0,
51  1250.0, 1250.0, 1378.0, 1615.0, 1615.0, 2130.0, 2260.0
52  };
53 
54 const double dtdb_wavelengths[NUM_DTDB_BANDS] = {
55  410.0, 490.0, 550.0, 670.0, 865.0, 1040.0,
56  1250.0, 1378.0, 1615.0, 1615.0, 2130.0, 2260.0
57  };
58 
59 /**************************************************************************
60  * NAME: TmProcess()
61  *
62  * DESCRIPTION: Class Constructor
63  *
64  *************************************************************************/
65 
67 }
68 
69 /**************************************************************************
70  * NAME: ~TmProcess()
71  *
72  * DESCRIPTION: Class Destructor
73  *
74  *************************************************************************/
75 
77  delete tin_;
78  delete tout_;
79 }
80 
81 /**************************************************************************
82  * NAME: initialize()
83  *
84  * DESCRIPTION: Virtual function initializes data and object classes for
85  * process operations.
86  *
87  *************************************************************************/
88 
90  int status = TM_SUCCESS;
91 
92  tin_ = new (tmatrix_in);
93 
94  string str = tm_get_option(INPUT_RAT);
95  rad_type_ = str;
96  if (str.empty() || (str == "EQUAL_VOLUME")) {
97  tin_->rat = EQUAL_VOLUME;
98  } else if (str == "EQUAL_AREA") {
99  tin_->rat = EQUAL_AREA;
100  }
102  shape_ = str;
103  if (str.empty() || (str == "SPHEROID")) {
104  tin_->np = SPHEROID;
105  } else if (str == "CYLINDER") {
106  tin_->np = CYLINDER;
107  } else if (str == "CHEBYSHEV") {
108  tin_->np = CHEBYSHEV;
109  }
111  distribution_ = str;
112  if (str.empty() || (str == "MODIFIED_GAMMA")) {
113  tin_->ndistr = MODIFIED_GAMMA;
114  } else if (str == "LOGNORMAL") {
115  tin_->ndistr = LOGNORMAL;
116  } else if (str == "POWERLAW") {
117  tin_->ndistr = POWERLAW;
118  } else if (str == "GAMMA") {
119  tin_->ndistr = GAMMA;
120  } else if (str == "MODIFIED_POWERLAW") {
121  tin_->ndistr = MODIFIED_POWERLAW;
122  }
123  tin_->lam = tm_get_option_doubles(INPUT_LAM, &num_bands_);
124  tin_->eps_max = tm_get_option_double(INPUT_EPS_MAX);
125  tin_->eps_min = tm_get_option_double(INPUT_EPS_MIN);
126  tin_->eps_num = tm_get_option_int(INPUT_EPS_NUM);
127  tin_->mrr_max = tm_get_option_double(INPUT_MRR_MAX);
128  tin_->mrr_min = tm_get_option_double(INPUT_MRR_MIN);
129  tin_->mrr_num = tm_get_option_int(INPUT_MRR_NUM);
130  tin_->mri_max = tm_get_option_double(INPUT_MRI_MAX);
131  tin_->mri_min = tm_get_option_double(INPUT_MRI_MIN);
132  tin_->mri_num = tm_get_option_int(INPUT_MRI_NUM);
133  tin_->b_max = tm_get_option_double(INPUT_B_MAX);
134  tin_->b_min = tm_get_option_double(INPUT_B_MIN);
135  tin_->b_num = tm_get_option_int(INPUT_B_NUM);
136  tin_->ddelt = tm_get_option_double(INPUT_DDELT);
137  tin_->axmax = tm_get_option_double(INPUT_AXMAX);
138  tin_->gam = tm_get_option_double(INPUT_GAM);
139  tin_->ndgs = tm_get_option_int(INPUT_NDGS);
140  tin_->npnax = tm_get_option_int(INPUT_NPNAX);
141  tin_->nkmax = tm_get_option_int(INPUT_NKMAX);
142  tin_->npna = tm_get_option_int(INPUT_NPNA);
143  tin_->ncoeff = tm_get_option_int(INPUT_NCOEFF);
144 
145  int nwl = num_bands_;
146  int neps = tin_->eps_num;
147  int nmrr = tin_->mrr_num;
148  int nmri = tin_->mri_num;
149  int nb = tin_->b_num;
150  int npax = tin_->npnax;
151  int nang = tin_->npna;
152  tout_ = new (tmatrix_out);
153  tout_->wl.resize(boost::extents[nwl]);
154  tout_->eps.resize(boost::extents[neps]);
155  tout_->mrr.resize(boost::extents[nmrr]);
156  tout_->mri.resize(boost::extents[nmri]);
157  tout_->b.resize(boost::extents[nb]);
158  tout_->a.resize(boost::extents[npax]);
159  tout_->angles.resize(boost::extents[nang]);
160  tout_->reff.resize(boost::extents[1][neps][nmrr][nmri][nb][npax]);
161  tout_->veff.resize(boost::extents[1][neps][nmrr][nmri][nb][npax]);
162  tout_->cext.resize(boost::extents[1][neps][nmrr][nmri][nb][npax]);
163  tout_->cscat.resize(boost::extents[1][neps][nmrr][nmri][nb][npax]);
164  tout_->walb.resize(boost::extents[1][neps][nmrr][nmri][nb][npax]);
165  tout_->asymm.resize(boost::extents[1][neps][nmrr][nmri][nb][npax]);
166  tout_->alpha.resize(boost::extents[neps][nmrr][nmri][nb][npax][NPL][ALPHA_COEFF]);
167  tout_->beta.resize(boost::extents[neps][nmrr][nmri][nb][npax][NPL][BETA_COEFF]);
168  tout_->f.resize(
169  boost::extents[1][neps][nmrr][nmri][nb][npax][nang][STOKES6]);
170 
171  for (int iB = 0; iB < num_bands_; iB++) {
172  tout_->wl[iB] = tin_->lam[iB] * MILLI;
173  }
174  if (tin_->eps_num < 2) {
175  tout_->eps[0] = tin_->eps_max;
176  tout_->eps[0] = (tout_->eps[0]==1.0) ? tout_->eps[0]+1.0e-6 : tout_->eps[0];
177  } else {
178  for (int i = 0; i < tin_->eps_num; i++) {
179  tout_->eps[i] = tin_->eps_min +
180  i*(tin_->eps_max - tin_->eps_min)/(tin_->eps_num-1);
181  tout_->eps[i] = (tout_->eps[i]==1.0) ? tout_->eps[i]+1.0e-6 : tout_->eps[i];
182  }
183  }
184  if (tin_->mrr_num < 2) {
185  tout_->mrr[0] = tin_->mrr_max;
186  } else {
187  for (int i = 0; i < tin_->mrr_num; i++) {
188  tout_->mrr[i] = tin_->mrr_min +
189  i*(tin_->mrr_max - tin_->mrr_min)/(tin_->mrr_num-1);
190  }
191  }
192  if (tin_->mri_num < 2) {
193  tout_->mri[0] = tin_->mri_max;
194  } else {
195  for (int i = 0; i < tin_->mri_num; i++) {
196  tout_->mri[i] = tin_->mri_min +
197  i*(tin_->mri_max - tin_->mri_min)/(tin_->mri_num-1);
198  }
199  }
200  if (tin_->b_num < 2) {
201  tout_->b[0] = tin_->b_max;
202  } else {
203  for (int i = 0; i < tin_->b_num; i++) {
204  tout_->b[i] = tin_->b_min +
205  i*(tin_->b_max - tin_->b_min)/(tin_->b_num-1);
206  }
207  }
208  if (tin_->npna < 2) {
209  tout_->angles[0] = 0.0;
210  } else {
211  for (int i = 0; i < tin_->npna; i++) {
212  tout_->angles[i] = i*180.0/(tin_->npna-1);
213  }
214  }
215 
216  status = create_output();
217 
218  return status;
219 }
220 
221 /**************************************************************************
222  * NAME: compute()
223  *
224  * DESCRIPTION: Virtual function executes process algorithm.
225  *
226  *************************************************************************/
227 
229  int status = TM_SUCCESS;
230 
231  for (int iB = 0; iB < NUM_PACE_BANDS; iB++) {
232 
233  double lam = dtdb_wavelengths[iB] * MILLI;
234 
235  calcrand_(&tin_->rat, &lam,
236  &tin_->mrr_max, &tin_->mri_max,
237  &tin_->eps_max, &tin_->np,
238  &tin_->ddelt, &tin_->ndgs,
239  &tin_->npnax, &tin_->axmax,
240  &tin_->b_max, &tin_->gam,
241  &tin_->ndistr, &tin_->nkmax,
242  &tin_->npna, &tin_->ncoeff,
243  &tout_->reff[iB][0][0][0][0][0],
244  &tout_->veff[iB][0][0][0][0][0],
245  &tout_->cext[iB][0][0][0][0][0],
246  &tout_->cscat[iB][0][0][0][0][0],
247  &tout_->walb[iB][0][0][0][0][0],
248  &tout_->asymm[iB][0][0][0][0][0],
249  &tout_->alpha[0][0][0][0][0][0][0],
250  &tout_->beta[0][0][0][0][0][0][0],
251  &tout_->f[iB][0][0][0][0][0][0][0],
252  &tout_->a[0]);
253 
254  cout << "Wavelength = " << lam << endl;
255  }
256 
257  return status;
258 }
259 
260 /**************************************************************************
261  * NAME: compute()
262  *
263  * DESCRIPTION: Virtual function executes process algorithm.
264  *
265  *************************************************************************/
266 
268  int status = TM_SUCCESS;
269 
270  for (int iB = 0; iB < num_bands_; iB++) {
271  for (int iE = 0; iE < tin_->eps_num; iE++) {
272  for (int iR = 0; iR < tin_->mrr_num; iR++) {
273  for (int iI = 0; iI < tin_->mri_num; iI++) {
274  for (int iV = 0; iV < tin_->b_num; iV++) {
275 
276  double ddelt = tin_->ddelt;
277  int status = calcrand_(&tin_->rat, &tout_->wl[iB],
278  &tout_->mrr[iR], &tout_->mri[iI],
279  &tout_->eps[iE], &tin_->np,
280  &ddelt, &tin_->ndgs,
281  &tin_->npnax, &tin_->axmax,
282  &tout_->b[iV], &tin_->gam,
283  &tin_->ndistr, &tin_->nkmax,
284  &tin_->npna, &tin_->ncoeff,
285  &tout_->reff[0][iE][iR][iI][iV][0],
286  &tout_->veff[0][iE][iR][iI][iV][0],
287  &tout_->cext[0][iE][iR][iI][iV][0],
288  &tout_->cscat[0][iE][iR][iI][iV][0],
289  &tout_->walb[0][iE][iR][iI][iV][0],
290  &tout_->asymm[0][iE][iR][iI][iV][0],
291  &tout_->alpha[iE][iR][iI][iV][0][0][0],
292  &tout_->beta[iE][iR][iI][iV][0][0][0],
293  &tout_->f[0][iE][iR][iI][iV][0][0][0],
294  &tout_->a[0]);
295 
296  if (status != TM_SUCCESS) {
297  for( int iP=0; iP<tin_->npnax; iP++) {
298  tout_->reff[0][iE][iR][iI][iV][iP] = fillvalue;
299  tout_->veff[0][iE][iR][iI][iV][iP] = fillvalue;
300  tout_->cext[0][iE][iR][iI][iV][iP] = fillvalue;
301  tout_->cscat[0][iE][iR][iI][iV][iP] = fillvalue;
302  tout_->walb[0][iE][iR][iI][iV][iP] = fillvalue;
303  tout_->asymm[0][iE][iR][iI][iV][iP] = fillvalue;
304  for( int iA=0; iA<tin_->npna; iA++) {
305  for( int iS=0; iS<STOKES6; iS++) {
306  tout_->f[0][iE][iR][iI][iV][iP][iA][iS] = fillvalue;
307  }
308  }
309  }
310  cout << "\n Run failed for wl = " << tout_->wl[iB] <<
311  " axial ratio = " << tout_->eps[iE] <<
312  " mrr = " << tout_->mrr[iR] <<
313  " mri = " << tout_->mri[iI] <<
314  " var = " << tout_->b[iV] <<
315  " max_size = " << tin_->axmax << "\n" << endl;
316  } else {
317  cout << "\n Run complete for wl = " << tout_->wl[iB] <<
318  " axial ratio = " << tout_->eps[iE] <<
319  " mrr = " << tout_->mrr[iR] <<
320  " mri = " << tout_->mri[iI] <<
321  " var = " << tout_->b[iV] <<
322  " max_size = " << tin_->axmax << "\n" << endl;
323  }
324  }
325  }
326  }
327  }
328  status = write_wavelength(iB);
329  }
330 
331  return status;
332 }
333 
334 /**************************************************************************
335  * NAME: create_output()
336  *
337  * DESCRIPTION: Create NetCDF4 file containing 4x4 scattering matrices for
338  * configured wavelength and scattering angle
339  *
340  *************************************************************************/
341 
343  int status = TM_SUCCESS;
344 
345  output_filepath_ = tm_get_option(OUTPUT_NC4);
346  if (output_filepath_.empty()) {
347  cerr << "\nTmProcess:: Failure locating output file path.\n" << endl;
348  return TM_FAIL;
349  }
350  size_t pos = output_filepath_.rfind("/");
351  prod_name_ = output_filepath_.substr(pos + 1);
352  NcFile* nc_output;
353  try {
354  nc_output = new NcFile(output_filepath_, NcFile::replace);
355  } catch (NcException& e) {
356  e.what();
357  cerr
358  << "\nTmProcess:: Failure creating product file: "
359  + output_filepath_ + "\n" << endl;
360  return TM_FAIL;
361  }
362 
363  write_global_attributes(nc_output);
364 
365  NcDim wave_dim = nc_output->addDim("dim_wavelength", num_bands_);
366  NcDim eps_dim = nc_output->addDim("dim_axial_ratio", tin_->eps_num);
367  NcDim mrr_dim = nc_output->addDim("dim_mrr", tin_->mrr_num);
368  NcDim mri_dim = nc_output->addDim("dim_mri", tin_->mri_num);
369  NcDim var_dim = nc_output->addDim("dim_variance", tin_->b_num);
370  NcDim size_dim = nc_output->addDim("dim_radius", tin_->npnax);
371  NcDim scat_dim = nc_output->addDim("dim_angles", tin_->npna);
372  NcDim stokes_dim = nc_output->addDim("dim_stokes", STOKES6);
373  vector<NcDim> scat_dims;
374  scat_dims.push_back(wave_dim);
375  scat_dims.push_back(eps_dim);
376  scat_dims.push_back(mrr_dim);
377  scat_dims.push_back(mri_dim);
378  scat_dims.push_back(var_dim);
379  scat_dims.push_back(size_dim);
380  scat_dims.push_back(scat_dim);
381  scat_dims.push_back(stokes_dim);
382  vector<NcDim> size_dims;
383  size_dims.push_back(wave_dim);
384  size_dims.push_back(eps_dim);
385  size_dims.push_back(mrr_dim);
386  size_dims.push_back(mri_dim);
387  size_dims.push_back(var_dim);
388  size_dims.push_back(size_dim);
389 
390  NcVar var;
391 
392  var = nc_output->addVar("pts_wavelength", ncDouble, wave_dim);
393  var.putVar(&tout_->wl[0]);
394  var = nc_output->addVar("pts_axial_ratio", ncDouble, eps_dim);
395  var.putVar(&tout_->eps[0]);
396  var = nc_output->addVar("pts_mrr", ncDouble, mrr_dim);
397  var.putVar(&tout_->mrr[0]);
398  var = nc_output->addVar("pts_mri", ncDouble, mri_dim);
399  var.putVar(&tout_->mri[0]);
400  var = nc_output->addVar("pts_variance", ncDouble, var_dim);
401  var.putVar(&tout_->b[0]);
402  var = nc_output->addVar("pts_radius", ncDouble, size_dim);
403  var = nc_output->addVar("pts_angles", ncDouble, scat_dim);
404  var.putVar(&tout_->angles[0]);
405  var = nc_output->addVar("cext", ncDouble, size_dims);
406  var = nc_output->addVar("cscat", ncDouble, size_dims);
407  var = nc_output->addVar("walb", ncDouble, size_dims);
408  var = nc_output->addVar("asymmetry", ncDouble, size_dims);
409  var = nc_output->addVar("reff", ncDouble, size_dims);
410  var = nc_output->addVar("veff", ncDouble, size_dims);
411  var = nc_output->addVar("F", ncDouble, scat_dims);
412  var.putAtt("long_name", "Single scattering matrix parameters (F11,F22,F33,F44,F12,F34)");
413  var.putAtt("dimensions", "In order, wavelength, axial_ratio, ior_real, ior_imag, size_spread, size_radius, scattering_angles, scattering matrix parameters");
414  var.setFill(true, -999.9);
415 
416  delete nc_output;
417  return status;
418 }
419 
420 /**************************************************************************
421  * NAME: write_wavelength()
422  *
423  * DESCRIPTION: Create NetCDF4 LUT containing 4x4 scattering matrices for
424  * configured wavelength and scattering angle
425  *
426  *************************************************************************/
427 
428 int TmProcess::write_wavelength( const int wl) {
429  int status = TM_SUCCESS;
430 
431  output_filepath_ = tm_get_option(OUTPUT_NC4);
432  if (output_filepath_.empty()) {
433  cerr << "\nTmProcess:: Failure locating output file path.\n" << endl;
434  return TM_FAIL;
435  }
436  size_t pos = output_filepath_.rfind("/");
437  prod_name_ = output_filepath_.substr(pos + 1);
438  NcFile* nc_output;
439  try {
440  nc_output = new NcFile(output_filepath_, NcFile::write);
441  } catch (NcException& e) {
442  e.what();
443  cerr
444  << "\nTmProcess:: Failure creating product file: "
445  + output_filepath_ + "\n" << endl;
446  return TM_FAIL;
447  }
448 
449  vector<size_t> start;
450  start.push_back(wl);
451  start.push_back(0);
452  start.push_back(0);
453  start.push_back(0);
454  start.push_back(0);
455  start.push_back(0);
456  vector<size_t> count;
457  count.push_back(1);
458  count.push_back(tin_->eps_num);
459  count.push_back(tin_->mrr_num);
460  count.push_back(tin_->mri_num);
461  count.push_back(tin_->b_num);
462  count.push_back(tin_->npnax);
463 
464  NcVar var = nc_output->getVar("pts_radius");
465  var.putVar(&tout_->a[0]);
466 
467  var = nc_output->getVar("cext");
468  var.putVar(start, count, &tout_->cext[0][0][0][0][0][0]);
469 
470  var = nc_output->getVar("cscat");
471  var.putVar(start, count, &tout_->cscat[0][0][0][0][0][0]);
472 
473  var = nc_output->getVar("walb");
474  var.putVar(start, count, &tout_->walb[0][0][0][0][0][0]);
475 
476  var = nc_output->getVar("asymmetry");
477  var.putVar(start, count, &tout_->asymm[0][0][0][0][0][0]);
478 
479  var = nc_output->getVar("reff");
480  var.putVar(start, count, &tout_->reff[0][0][0][0][0][0]);
481 
482  var = nc_output->getVar("veff");
483  var.putVar(start, count, &tout_->veff[0][0][0][0][0][0]);
484 
485  start.push_back(0);
486  start.push_back(0);
487  count.push_back(tin_->npna);
488  count.push_back(STOKES6);
489  var = nc_output->getVar("F");
490  var.putVar(start, count, &tout_->f[0][0][0][0][0][0][0][0]);
491 
492  delete nc_output;
493  return status;
494 }
495 
496 /**************************************************************************
497  * write_global_attributes()
498  *
499  * Write global attributes to specified netCDF file ID
500  *
501  **************************************************************************/
502 
503 int TmProcess::write_global_attributes(NcFile* nc_output) {
504 // global attributes
505  title_ = "Particle Scattering Characteristics";
506  processing_version_ = "v0.0.0";
507  Conventions_ = "CF-1.6";
508  institution_ =
509  "NASA Goddard Space Flight Center, Ocean Biology Group";
510  license_ =
511  "http://science.nasa.gov/earth-science/earth-science-data/data-information-policy/";
512  naming_authority_ = "gov.nasa.gsfc";
513 
514  time_t rawtime;
515  struct tm * timeinfo;
516  char buffer[80];
517  time(&rawtime);
518  timeinfo = localtime(&rawtime);
519  strftime(buffer, 80, "%Y-%m-%dT%I:%M:%SZ", timeinfo);
520  date_created_ = (string) buffer;
521  nc_output->putAtt("processing_version", processing_version_);
522 
523  nc_output->putAtt("Conventions", Conventions_);
524  nc_output->putAtt("institution", institution_);
525  nc_output->putAtt("license", license_);
526  nc_output->putAtt("naming_authority", naming_authority_);
527  nc_output->putAtt("date_created", date_created_);
528  nc_output->putAtt("history", history_);
529  nc_output->putAtt("source", source_files_);
530  string versionid_ = basename((char*) tm_get_option("VersionID").c_str());
531  if (!versionid_.empty()) {
532  nc_output->putAtt("VersionId", versionid_);
533  }
534 
535  nc_output->putAtt("title", title_);
536  nc_output->putAtt("product_name", prod_name_);
537  nc_output->putAtt("particle_shape", shape_);
538  nc_output->putAtt("size_distribution", distribution_);
539  nc_output->putAtt("radius_type", rad_type_);
540 
541  return TM_SUCCESS;
542 }
543 
const std::string INPUT_MRR_NUM
const int TM_FAIL
Definition: TmConstants.h:48
const std::string INPUT_NPNA
int status
Definition: l1_czcs_hdf.c:32
const std::string INPUT_NDISTR
const double MILLI
const double pace_wavelengths[NUM_PACE_BANDS]
Definition: TmProcess.cpp:35
double tm_get_option_double(const std::string &name)
const std::string INPUT_MRI_NUM
const std::string INPUT_NDGS
const int TM_SUCCESS
Definition: TmConstants.h:47
float32 * pos
Definition: l1_czcs_hdf.c:35
float tm[MODELMAX]
#define NP
Definition: l1_mos_hdf.c:25
const std::string INPUT_GAM
@ string
#define EPS
Definition: fuzzy_func_v3.c:26
virtual int compute_pace()
Definition: TmProcess.cpp:228
const double F
const std::string INPUT_MRR_MIN
double * tm_get_option_doubles(const std::string &name, int *count)
int write_wavelength(const int wl)
Definition: TmProcess.cpp:428
const std::string INPUT_EPS_MAX
virtual ~TmProcess()
Definition: TmProcess.cpp:76
int create_output()
Definition: TmProcess.cpp:342
int calcrand_(double *RAT, double *LAM, double *MRR, double *MRI, double *EPS, int *NP, double *DDELT, int *NDGS, int *NPNAX, double *AXMAX, double *B, double *GAM, int *NDISTR, int *NKMAX, int *NPNA, int *NCOEFF, double *REFF, double *VEFF, double *CEXTIN, double *CSCAT, double *WALB, double *ASYMM, double *ALPHA, double *BETA, double *F, double *SZ)
virtual int compute_dtdb()
Definition: TmProcess.cpp:267
const std::string INPUT_B_NUM
const std::string INPUT_DDELT
const std::string INPUT_MRR_MAX
const std::string INPUT_NP
const std::string INPUT_LAM
int write_global_attributes(NcFile *nc_output)
Definition: TmProcess.cpp:503
const std::string INPUT_EPS_MIN
virtual int initialize()
Definition: TmProcess.cpp:89
const string OUTPUT_NC4
Definition: DDOptions.cpp:49
const std::string INPUT_RAT
#define basename(s)
Definition: l0chunk_modis.c:29
const char * str
Definition: l1c_msi.cpp:35
int tm_get_option_int(const std::string &name)
const double dtdb_wavelengths[NUM_DTDB_BANDS]
Definition: TmProcess.cpp:54
const std::string INPUT_MRI_MIN
this program makes no use of any feature of the SDP Toolkit that could generate such a then geolocation is calculated at that and then aggregated up to Resolved feature request Bug by adding three new int8 SDSs for each high resolution offsets between the high resolution geolocation and a bi linear interpolation extrapolation of the positions This can be used to reconstruct the high resolution geolocation Resolved Bug by delaying cumulation of gflags until after validation of derived products Resolved Bug by setting Latitude and Longitude to the correct fill resolving to support Near Real Time because they may be unnecessary if use of entrained ephemeris and attitude data is turned resolving bug report Corrected to filter out Aqua attitude records with missing status helping resolve bug MOD_PR03 will still correctly write scan and pixel data that does not depend upon the start time
Definition: HISTORY.txt:248
const std::string INPUT_MRI_MAX
const std::string INPUT_NPNAX
const std::string INPUT_B_MAX
const std::string INPUT_B_MIN
int i
Definition: decode_rs.h:71
const std::string INPUT_NCOEFF
std::string tm_get_option(const std::string &name)
const std::string INPUT_AXMAX
const std::string INPUT_NKMAX
int32_t nb
Definition: atrem_corl1.h:132
const std::string INPUT_EPS_NUM
int count
Definition: decode_rs.h:79