A sphere with three ocean waves in differing shades of blue. Next to the sphere, there is the bolded text "Ocean Color". Under that, there are two acronyms, separated by a vertical pipe: OB.DAAC (Ocean Biology Distributed Active Archive Center) and OBPG (Ocean Biology Processing Group). Lorem ipsum dolor sit amet
NASA Logo
Ocean Color Science Software

ocssw V2022
l1_modis.c
Go to the documentation of this file.
1 /* ============================================================================ */
2 /* module l1_hmodis_hdf.c - functions to read MODIS HIRES L1B for MSL12 */
3 /* */
4 /* Written By: B. Franz, NASA/SIMBIOS, January 2003. */
5 /* Modified By: J. Gales, NASA/OBPG, August 2005. */
6 /* Conversion for 16-band MODIS: B. Franz, J. Gales, 2006. */
7 /* Completely restructured by G. Fireman in 2014. */
8 /* */
9 /* ============================================================================ */
10 
11 #define _XOPEN_SOURCE_EXTENDED /* for strdup() */
12 #include "libnav.h"
13 #include "hdf4utils.h"
14 #include "l1.h"
15 #include <libgen.h>
16 #include <hdf.h>
17 #include <mfhdf.h>
18 #undef _XOPEN_SOURCE_EXTENDED
19 
20 #define SCAN_TIME_INTERVAL 1.4771
21 
22 /***********************************************************************/
23 #define TRYMEM(file,line,memstat) { \
24  if (memstat == NULL) { \
25  fprintf(stderr, \
26  "-E- %s line %d: Memory allocation error.\n", \
27  file, line); \
28  exit(1); } \
29  }
30 /***********************************************************************/
31 
32 /* SDSs of interest in each type of input file */
33 typedef struct {
34  int32_t index;
35  int32_t scandim;
36  char *name;
37 } sdslist;
38 
39 /* Positioning info to preload from MODIS Geolocation file */
40 enum REF_SDS {
41  GEO_TAISEC,
42  GEO_MSIDE,
43  GEO_ANGLES,
44  GEO_MNORM,
45  REF_NUM_SDS
46 };
47 static const sdslist SDSLIST_REF[] = {
48  { GEO_TAISEC, 0, "EV start time"},
49  { GEO_MSIDE, 0, "Mirror side"},
50  { GEO_ANGLES, 0, "attitude_angles"},
51  { GEO_MNORM, 0, "T_inst2ECR"},
52  { 0, 0, NULL}
53 };
54 
55 /* MODIS Geolocation swath SDSs */
56 enum GEO_SDS {
57  GEO_LON,
58  GEO_LAT,
59  GEO_HGT,
60  GEO_SOLZ,
61  GEO_SOLA,
62  GEO_SENZ,
63  GEO_SENA,
64  GEO_NUM_SDS
65 };
66 static const sdslist SDSLIST_GEO[] = {
67  { GEO_LON, 0, "Longitude"}, /* float32 */
68  { GEO_LAT, 0, "Latitude"}, /* float32 */
69  { GEO_HGT, 0, "Height"}, /* int16 */
70  { GEO_SOLZ, 0, "SolarZenith"}, /* int16 */
71  { GEO_SOLA, 0, "SolarAzimuth"}, /* int16 */
72  { GEO_SENZ, 0, "SensorZenith"}, /* int16 */
73  { GEO_SENA, 0, "SensorAzimuth"}, /* int16 */
74  { 0, 0, NULL}
75 };
76 
77 enum GEO_COEFFS {
78  GEO_SCALE,
79  GEO_NUM_COEFFS
80 };
81 static const char* geo_coeff[] = {/* index with enum GEO_COEFFS */
82  "scale_factor"
83 };
84 
85 /* MODIS L1B Band groupings */
86 enum L1B_SDS {
87  RSB_250,
88  RSB_500,
89  RSB_1KM,
90  CIR_1KM,
91  TEB_1KM,
92  L1B_NUM_SDS
93 };
94 static sdslist SDSLIST_1KM[] = {// 1KM file contents
95  { RSB_250, 1, "EV_250_Aggr1km_RefSB"},
96  { RSB_500, 1, "EV_500_Aggr1km_RefSB"},
97  { RSB_1KM, 1, "EV_1KM_RefSB"},
98  { CIR_1KM, 0, "EV_Band26"}, // 1st dim=scan, not band
99  { TEB_1KM, 1, "EV_1KM_Emissive"},
100  { 0, 0, NULL}
101 };
102 static sdslist SDSLIST_HKM[] = {// HKM file contents
103  { RSB_250, 1, "EV_250_Aggr500_RefSB"},
104  { RSB_500, 1, "EV_500_RefSB"},
105  { 0, 0, NULL}
106 };
107 static sdslist SDSLIST_QKM[] = {// QKM file
108  { RSB_250, 1, "EV_250_RefSB"},
109  { 0, 0, NULL}
110 };
111 
112 enum L1B_COEFFS {
113  REFL_SCALE,
114  REFL_OFFSET,
115  RAD_SCALE,
116  RAD_OFFSET,
117  L1B_NUM_COEFFS
118 };
119 static const char* l1b_coeff[] = {/* index with enum L1B_COEFFS */
120  "reflectance_scales",
121  "reflectance_offsets",
122  "radiance_scales",
123  "radiance_offsets"
124 };
125 
126 /***********************************************************************/
127 
128 /* Structure holding info for each MODIS input file */
129 typedef struct {
130  int32_t id;
131  char file[FILENAME_MAX];
132  char shortname[FILENAME_MAX];
133  double starttime;
134  int16_t resolution;
135  int32_t nlines;
136  int32_t npixls;
137  int32_t n_sds;
138  sds_struct *sds;
139 } modis_file;
140 
141 /***********************************************************************/
142 
143 /* Structure holding info needed for each swath-based SDS */
144 typedef struct {
145  sds_struct sds;
146  int16_t resolution;
147  int32_t nbands;
148  int32_t nscans;
149  int32_t ndets;
150  int32_t nframes;
151  int32_t nsf;
152  int32_t scandim;
153 } modis_sds;
154 
155 /*
156  Notes:
157  All GEO arrays are dimensioned [nscans*ndets, nframes]
158  Most L1B arrays are dimensioned [nbands, nscans*ndets, nframes*nsf]
159  where:
160  nscans and nframes are at 1KM resolution
161  nsf = #subframes = 1000/resolution
162  ndets = 10*nsubframes
163 
164  Band 26 drops the 1st dimension.
165  Geolocation is inherently 1KM resolution: nbands==1, ndets==10, nsf==1.
166 
167  ndets and scandim are needed for reading data one scan at a time.
168 
169  There is some redundancy between resolution, ndets, and nsf,
170  but it's helpful to have them all defined.
171  */
172 
173 /***********************************************************************/
174 /* Global variables */
175 
176 static modis_file file_geo, file_1km, file_hkm, file_qkm;
177 static modis_sds ref[REF_NUM_SDS], geo[GEO_NUM_SDS], l1b[L1B_NUM_SDS];
178 
179 /***********************************************************************/
180 
181 /* Substitute provided string for LAC/1KM part of filename */
182 int modpath_1km(const char *oldpath, const char* newchars, char* newpath) {
183  int32_t status;
184  char *ptr = NULL;
185  char *tmpbuf = NULL;
186  char *filename = NULL;
187 
188  /* directory name */
189  tmpbuf = strdup(oldpath);
190  strcpy(newpath, dirname(tmpbuf));
191  strcat(newpath, "/");
192 
193  /* filename up to substring */
194  tmpbuf = strdup(oldpath);
195  filename = basename(tmpbuf);
196  status = ((ptr = strstr(filename, "_MODIS.")) == NULL);
197  if (status) {
198  status = (((ptr = strstr(filename, "LAC")) == NULL) &&
199  ((ptr = strstr(filename, "1KM")) == NULL));
200 
201  if (status) {
202  fprintf(stderr,
203  "Input file %s does not follow the standard name convention;"
204  " can't find \"%s\" file.\n",
205  oldpath, newchars);
206  newpath[0] = '\0';
207  return status;
208  }
209 
210  strncat(newpath, filename, ptr - filename);
211 
212  /* add new substring and the rest of the filename */
213  strcat(newpath, newchars);
214  strcat(newpath, ptr + strlen(newchars));
215  } else {
216  strncat(newpath, filename, ptr - filename + 6);
217 
218  /* add new substring and the rest of the filename */
219  strcat(newpath, "_");
220  strcat(newpath, newchars);
221  strcat(newpath, ptr + strlen(newchars) + 3);
222  }
223  free(tmpbuf);
224 
225  return status;
226 }
227 
228 /***********************************************************************/
229 
230 int open_modis_l1bfile(modis_file *l1bfile) {
231  int32_t status = FAIL;
232  int32_t i;
233  char result[FILENAME_MAX];
234  sdslist *file_content;
235 
236  /* Open input file */
237  status = ((l1bfile->id = SDstart(l1bfile->file, DFACC_RDONLY)) == FAIL);
238  if (status) {
239  fopen_warn(l1bfile->file, __FILE__, __LINE__);
240  return status;
241  }
242 
243  /*----- Read Global Attributes to determine file type -----*/
244  get_hdfeos_meta(l1bfile->id, "CoreMetadata.0", "SHORTNAME", result);
245  if ((strcmp(result, "MYD021KM") == 0) ||
246  (strcmp(result, "MOD021KM") == 0)) { // 1KM
247  l1bfile->resolution = 1000;
248  file_content = SDSLIST_1KM;
249  } else if ((strcmp(result, "MYD02HKM") == 0) ||
250  (strcmp(result, "MOD02HKM") == 0)) { // HKM
251  l1bfile->resolution = 500;
252  file_content = SDSLIST_HKM;
253  } else if ((strcmp(result, "MYD02QKM") == 0) ||
254  (strcmp(result, "MOD02QKM") == 0)) { // QKM
255  l1bfile->resolution = 250;
256  file_content = SDSLIST_QKM;
257  } else {
258  fprintf(stderr,
259  "Input file %s has type %s; "
260  "does not contain MODIS calibrated radiances.\n",
261  l1bfile->file, result);
262  return FAIL;
263  }
264  strcpy(l1bfile->shortname, result);
265 
266  /*----- Data set info -----*/
267  l1bfile->n_sds = 0;
268  while (file_content[l1bfile->n_sds].name)
269  l1bfile->n_sds++;
270  TRYMEM(__FILE__, __LINE__,
271  (l1bfile->sds = calloc(l1bfile->n_sds, sizeof (sds_struct))));
272  for (i = 0; i < l1bfile->n_sds; i++)
273  init_sds_byname(l1bfile->id, file_content[i].name, &(l1bfile->sds[i]));
274 
275  /*----- Native dimensions -----*/
276  l1bfile->nlines = l1bfile->sds[RSB_250].dimlen[1]; //(bands,NSCANS*NDETS,nframes*nsf)
277  l1bfile->npixls = l1bfile->sds[RSB_250].dimlen[2]; //(bands,nscans*ndets,NFRAMES*NSF)
278 
279  /* Note: Every MODIS L1B file contains a 250m-band SDS, aggregated for
280  coarser resolutions, so we can use [RSB_250] to reference native
281  dimensions for every file. Later we'll use [RSB_250] to index the
282  highest-resolution SDS, to specify output dimensions. */
283 
284  return SUCCESS;
285 }
286 
287 int init_l1b(const char filename[FILENAME_MAX], int32_t *max_resolution) {
288  int32_t i, j;
289  int32_t status = FAIL;
290 
291  /* Start fresh */
292  memset(&file_1km, 0, sizeof (modis_file));
293  memset(&file_hkm, 0, sizeof (modis_file));
294  memset(&file_qkm, 0, sizeof (modis_file));
295 
296  /* Initialize L1B 1KM file info */
297  strcpy(file_1km.file, filename);
298  status = open_modis_l1bfile(&file_1km);
299  if (status) {
300  fprintf(stderr,
301  "Error reading %s; please specify a valid Level 1B file.\n",
302  filename);
303  return status;
304  }
305 
306  /* Input L1B file must have 1KM resolution */
307  status = (file_1km.resolution != 1000);
308  if (status) {
309  fprintf(stderr,
310  "Input file %s has %dm resolution; "
311  "please specify a 1km-resolution Level 1B file.\n",
312  filename, file_1km.resolution);
313  return status;
314  }
315 
316  /* Open other L1B files as needed for specified resolution */
317  switch (*max_resolution) {
318 
319  case 250: // QKM
320  status = modpath_1km(filename, "QKM", file_qkm.file);
321  if (status)
322  return status;
323  status = open_modis_l1bfile(&file_qkm);
324  if (status) {
325  fprintf(stderr, "File not found: %s\n", file_qkm.file);
326  fprintf(stderr,
327  "Processing at %im resolution requires a QKM file to be "
328  "present in the same directory as the 1KM L1B file.\n",
329  *max_resolution);
330  return status;
331  }
332  /* no break: also load HKM */
333 
334  case 500: // HKM
335  status = modpath_1km(filename, "HKM", file_hkm.file);
336  if (status)
337  return status;
338  status = open_modis_l1bfile(&file_hkm);
339  if (status) {
340  fprintf(stderr, "File not found: %s\n", file_hkm.file);
341  fprintf(stderr,
342  "Processing at %im resolution requires a HKM file to be "
343  "present in the same directory as the 1KM L1B file.\n",
344  *max_resolution);
345  return status;
346  }
347  break;
348 
349  case -1: // default to 1KM
350  *max_resolution = 1000;
351  /* no break: proceed as for 1KM */
352 
353  case 1000: // 1KM
354  /* file_1km has already been initialized */
355  break;
356 
357  default:
358  fprintf(stderr, "Invalid resolution %i; ", *max_resolution);
359  fprintf(stderr, "must be 1000, 500 or 250 m.\n");
360  status = FAIL;
361  exit(status);
362  break;
363  }
364 
365  /* Populate global info for each L1B SDS */
366  for (i = 0; i < L1B_NUM_SDS; i++) {
367  l1b[i].sds = file_1km.sds[i]; // initialize to 1KM inputs
368  l1b[i].resolution = file_1km.resolution;
369  l1b[i].nscans = file_1km.nlines / 10;
370  l1b[i].nframes = file_1km.npixls;
371  l1b[i].scandim = SDSLIST_1KM[i].scandim;
372  }
373  if (*max_resolution < 1000) {
374  l1b[RSB_500].sds = file_hkm.sds[RSB_500]; // HKM or QKM
375  l1b[RSB_500].resolution = file_hkm.resolution;
376 
377  l1b[RSB_250].sds = file_hkm.sds[RSB_250]; // HKM only
378  l1b[RSB_250].resolution = file_hkm.resolution;
379  }
380  if (*max_resolution < 500) {
381  l1b[RSB_250].sds = file_qkm.sds[RSB_250]; // QKM only
382  l1b[RSB_250].resolution = file_qkm.resolution;
383  }
384 
385  /* Derived values */
386  for (i = 0; i < L1B_NUM_SDS; i++) {
387  l1b[i].nbands =
388  (l1b[i].scandim == 0) ? 1 // Band 26 only
389  : l1b[i].sds.dimlen[0];
390  l1b[i].nsf = 1000 / l1b[i].resolution;
391  l1b[i].ndets = 10 * l1b[i].nsf;
392 
393  /* Load scaling coefficient attribute(s) */
394  TRYMEM(__FILE__, __LINE__,
395  (l1b[i].sds.atts = calloc(L1B_NUM_COEFFS, sizeof (att_struct))));
396  l1b[i].sds.natts = L1B_NUM_COEFFS;
397  for (j = 0; j < L1B_NUM_COEFFS; j++)
398  load_att_byname(l1b[i].sds.id, l1b_coeff[j], &l1b[i].sds.atts[j]);
399  }
400 
401  return status;
402 }
403 
404 /***********************************************************************/
405 
406 int open_modis_geofile(modis_file *geofile) {
407  int32_t status = FAIL;
408  int32_t i;
409  char result[FILENAME_MAX];
410 
411  /* Open input file */
412  status = ((geofile->id = SDstart(geofile->file, DFACC_RDONLY)) == FAIL);
413  if (status) {
414  fopen_warn(geofile->file, __FILE__, __LINE__);
415  return status;
416  }
417 
418  /*----- Read Global Attributes to determine file type -----*/
419  get_hdfeos_meta(geofile->id, "CoreMetadata.0", "SHORTNAME", result);
420  if ((strcmp(result, "MYD03") == 0) ||
421  (strcmp(result, "MOD03") == 0)) { // GEO
422  } else {
423  fprintf(stderr, "Input file %s has type %s; "
424  "does not contain MODIS geolocation.\n",
425  geofile->file, result);
426  return FAIL;
427  }
428  strcpy(geofile->shortname, result);
429 
430  /*----- Data set info -----*/
431  geofile->n_sds = GEO_NUM_SDS;
432  TRYMEM(__FILE__, __LINE__,
433  (geofile->sds = calloc(GEO_NUM_SDS, sizeof (sds_struct))));
434  for (i = 0; i < GEO_NUM_SDS; i++)
435  init_sds_byname(geofile->id, SDSLIST_GEO[i].name, &(geofile->sds[i]));
436 
437  /*----- Native dimensions -----*/
438  geofile->nlines = geofile->sds[GEO_LAT].dimlen[0];
439  geofile->npixls = geofile->sds[GEO_LAT].dimlen[1];
440 
441  return SUCCESS;
442 }
443 
444 int init_geo(const char filename[FILENAME_MAX]) {
445  int32_t i, j;
446  int32_t status = FAIL;
447 
448  /* Start fresh */
449  memset(&file_geo, 0, sizeof (modis_file));
450 
451  /* Initialize GEO file info */
452  strcpy(file_geo.file, filename);
453  status = open_modis_geofile(&file_geo);
454  if (status) {
455  fprintf(stderr,
456  "Error reading %s; please specify a valid geolocation file.\n",
457  filename);
458  return status;
459  }
460 
461  /* Populate global info for each GEO SDS */
462  for (i = 0; i < GEO_NUM_SDS; i++) {
463  geo[i].sds = file_geo.sds[i];
464  geo[i].resolution = 1000;
465  geo[i].nbands = 1;
466  geo[i].nscans = file_geo.nlines / 10;
467  geo[i].ndets = 10;
468  geo[i].nframes = file_geo.npixls;
469  geo[i].nsf = 1;
470  geo[i].scandim = SDSLIST_GEO[i].scandim;
471 
472  /* Load scaling coefficient attribute(s) */
473  TRYMEM(__FILE__, __LINE__,
474  (geo[i].sds.atts = calloc(GEO_NUM_COEFFS, sizeof (att_struct))));
475  geo[i].sds.natts = GEO_NUM_COEFFS;
476  for (j = 0; j < GEO_NUM_COEFFS; j++)
477  load_att_byname(geo[i].sds.id, geo_coeff[j], &geo[i].sds.atts[j]);
478  }
479 
480  /*----- Read preloaded variables -----*/
481  for (i = 0; i < REF_NUM_SDS; i++) {
482  init_sds_byname(file_geo.id, SDSLIST_REF[i].name, &(ref[i].sds));
483  readall_sds(&(ref[i].sds));
484  }
485 
486  return status;
487 }
488 
489 int read_sds_1scan(sds_struct *sds,
490  int32_t iscan,
491  int32_t scandim,
492  int32_t per_scan) {
493  int32_t i;
494  int32_t status = FAIL;
495  int32_t nvals = 1;
496  int32_t *start, *edges;
497 
498  /* don't bother if SDS handle is invalid */
499  if (sds->id != FAIL) {
500 
501  /* set array indices */
502  TRYMEM(__FILE__, __LINE__,
503  (start = malloc(sds->ndims * sizeof (int32_t))));
504  TRYMEM(__FILE__, __LINE__,
505  (edges = malloc(sds->ndims * sizeof (int32_t))));
506  for (i = 0; i < sds->ndims; i++) {
507  if (i == scandim) {
508  start[i] = iscan * per_scan;
509  edges[i] = per_scan;
510  } else {
511  start[i] = 0;
512  edges[i] = sds->dimlen[i];
513  }
514  nvals *= edges[i];
515  }
516 
517  /* allocate one scan as needed (first call) */
518  if (sds->data == NULL)
519  TRYMEM(__FILE__, __LINE__,
520  (sds->data = malloc(nvals * DFKNTsize(sds->ntype))));
521 
522  /* read data */
523  status = SDreaddata(sds->id, start, NULL, edges, sds->data);
524  free(start);
525  free(edges);
526  }
527  if (status != FAIL)
528  return nvals;
529  else
530  return status;
531 }
532 
533 /***********************************************************************/
534 /* Everything below this line is specific to OCSSW implementation. */
535 /***********************************************************************/
536 #define MIN_BAD_SI 65500
537 #define MAX_ATTERR 0.017453293 /* radians, ~= 1 degree */
538 #define LOOPI for (i = 0; i < nvals; i++)
539 
540 /* Global variables from l1_hmodis_hdf.c */
541 static int32_t spixl = 0; /* first valid pixel of extracted granule */
542 static float *Fobar;
543 
544 /* RSMAS corrections to Terra IR radiances (20,22,23,27,28,29,31,32) */
545 /* latband v6 sst coefficients were estimated from Bt's from which these */
546 /* corrections have been applied, so we have to keep using them. */
547 static float radoff[][10] = {
548  {-0.000972, 0.014200, 0.000022, 0.000238, -0.000003,
549  0.003760, 0.002337, 0.000084, 0.008848, 0.005050},
550  {-0.002300, 0.000600, 0.000236, -0.000280, -0.000003,
551  0.002798, 0.003496, 0.018035, 0.002942, 0.001787},
552  {-0.003833, 0.003657, 0.002118, 0.000798, -0.000003,
553  0.000208, 0.000399, 0.000553, 0.000258, 0.021128},
554  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
555  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
556  {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
557  {-0.000423, -0.000242, -0.000330, -0.000065, -0.000001,
558  -0.000006, 0.000064, 0.000398, 0.000362, 0.000322},
559  {-0.000433, -0.000246, -0.000222, -0.000148, -0.000001,
560  0.000068, 0.000221, 0.000221, 0.000437, 0.000303}
561 };
562 static float mfact[] = {0.0, 0.0, 0.0005, 0.0, 0.0, 0.0, 0.0028, 0.00278};
563 
567 typedef struct {
568  int16_t sdsindex;
569  int16_t bandindex;
570  int16_t flagindex;
571  int16_t wavelength;
572  char *bandname;
573 } bandinfo;
574 
575 /* Bands used by l2gen, in expected order */
576 static const bandinfo BANDLIST[] = {
577  /*** Reflective Solar Bands ***/
578  { RSB_1KM, 0, 180, 412, "8"}, // 0
579  { RSB_1KM, 1, 190, 443, "9"}, // 1
580  { RSB_500, 0, 80, 469, "3"}, // 2
581  { RSB_1KM, 2, 200, 488, "10"}, // 3
582  { RSB_1KM, 3, 210, 531, "11"}, // 4
583  { RSB_1KM, 4, 220, 547, "12"}, // 5
584  { RSB_500, 1, 100, 555, "4"}, // 6
585  { RSB_250, 0, 0, 645, "1"}, // 7
586  { RSB_1KM, 5, 230, 667, "13lo"}, // 8
587  { RSB_1KM, 7, 250, 678, "14lo"}, // 9
588  { RSB_1KM, 9, 270, 748, "15"}, // 10
589  { RSB_250, 1, 40, 859, "2"}, // 11
590  { RSB_1KM, 10, 280, 869, "16"}, // 12
591  { RSB_500, 2, 120, 1240, "5"}, // 13
592  { RSB_500, 3, 140, 1640, "6"}, // 14
593  { RSB_500, 4, 160, 2130, "7"}, // 15
594  /*** Thermal Emissive Bands ***/
595  { TEB_1KM, 0, 330, 3750, "20"}, // 16
596  { TEB_1KM, 2, 350, 3959, "22"}, // 17
597  { TEB_1KM, 3, 360, 4050, "23"}, // 18
598  { TEB_1KM, 6, 390, 6715, "27"}, // 19
599  { TEB_1KM, 7, 400, 7325, "28"}, // 20
600  { TEB_1KM, 8, 410, 8550, "29"}, // 21
601  { TEB_1KM, 10, 430, 11000, "31"}, // 22
602  { TEB_1KM, 11, 440, 12000, "32"}, // 23
603  /*** Cirrus Band ***/
604  { CIR_1KM, 0, 320, 1375, "26"} // 24
605 };
606 
607 /***********************************************************************/
608 
614 int openl1_modis(filehandle *l1file) {
615  char geometa[EOSMETALEN] = "";
616  char result[EOSMETALEN] = "";
617  char datetime[32] = "";
618  int32_t resolution;
619  int32_t intval;
620  uint32_t gflags[8];
621 
622  /*----- Initialize L1B and GEO info -----*/
623  if (init_l1b(l1file->name, &l1_input->resolution))
624  exit(EXIT_FAILURE);
625  if (init_geo(l1file->geofile))
626  exit(EXIT_FAILURE);
627 
628  /* Make sure this is not an ocean-band subsetted L1A file */
629  int i;
630  for (i = 0; i < file_1km.n_sds; i++) {
631  if (strcmp(file_1km.sds[i].name, "EV_250_Aggr1km_RefSB") == 0) {
632  if (SDfindattr(file_1km.sds[i].id, "Rescaled Ocean R") != FAIL) {
633  fprintf(stderr, "\n"
634  "This L1B file contains only the ocean band subset.\n"
635  "Processing is no longer supported; exiting.\n");
636  exit(EXIT_FAILURE);
637  }
638  break;
639  }
640  }
641 
642  /*----- Populate filehandle structure -----*/
643  l1file->ndets = l1b[RSB_250].ndets;
644  l1file->nscan = l1b[RSB_250].nscans * l1b[RSB_250].ndets;
645  l1file->npix = l1b[RSB_250].nframes * l1b[RSB_250].nsf;
646 
647  /* Capture processing resolution */
648  resolution = l1b[RSB_250].resolution;
649  if (resolution == 1000)
650  strcpy(l1file->spatialResolution, "1 km");
651  else
652  sprintf(l1file->spatialResolution, "%d m", resolution);
653 
654  if (want_verbose) {
655  printf("Processing at %d meter resolution.\n", resolution);
656  printf(" 1000-meter file: %s\n", file_1km.file);
657  if (resolution < 1000) {
658  printf(" 500-meter file: %s\n", file_hkm.file);
659  if (resolution < 500)
660  printf(" 250-meter file: %s\n", file_qkm.file);
661  }
662  printf("\n");
663  }
664 
665  /*----- Read Level 1B metadata -----*/
666  get_hdfeos_meta(file_1km.id, "CoreMetadata.0", "INPUTPOINTER", result);
667  strcpy(l1_input->calfile, result);
668 
669  /*----- Read Geolocation metadata -----*/
670  if (read_att(file_geo.id, "CoreMetadata.0", geometa)) {
671  fprintf(stderr,
672  "-E- %s line %d: Error reading CoreMetadata.0 from %s.\n",
673  __FILE__, __LINE__, file_geo.file);
674  exit(EXIT_FAILURE);
675  }
676 
677  /* Orbit info */
678  parse_odl(geometa, "ORBITNUMBER", result);
679  if (strlen(result) > 0)
680  l1file->orbit_number = (int32_t) atoi(result);
681  parse_odl(geometa, "EQUATORCROSSINGLONGITUDE", result);
682  if (strlen(result) > 0)
683  l1file->orbit_node_lon = atof(result);
684 
685  parse_odl(geometa, "EQUATORCROSSINGDATE", result);
686  strcpy(datetime, result);
687  strcat(datetime, "T");
688  parse_odl(geometa, "EQUATORCROSSINGTIME", result);
689  strcat(datetime, result);
690  l1file->node_crossing_time = isodate2unix(datetime);
691 
692  /* Granule start time */
693  parse_odl(geometa, "RANGEBEGINNINGDATE", result);
694  strcpy(datetime, result);
695  strcat(datetime, "T");
696  parse_odl(geometa, "RANGEBEGINNINGTIME", result);
697  strcat(datetime, result);
698  file_geo.starttime = isodate2unix(datetime);
699 
700  /* Terrain height */
701  if (read_att(file_geo.id, "Cumulated gflags", gflags))
702  fprintf(stderr, "-W- %s line %d: Error reading gflags from %s.\n",
703  __FILE__, __LINE__, file_geo.file);
704  else
705  l1file->terrain_corrected = (((int32_t*) gflags)[5] == 0);
706 
707  /*** Adjust dimensions for extracted granules ***/
708  /*
709  When Level 1 MODIS granules are extracted, the full scan
710  width is retained but pixels outside the ROI are set to fill
711  value. We must keep track of pixel offset and count in
712  order to load the ROI correctly. Line offset and count are
713  handled in the main() loop.
714  */
715  if (read_att(file_1km.id, "Extract Pixel Count", &intval) == 0
716  && intval > 0) {
717  int32_t sline, nline, npixl;
718  npixl = intval * l1b[RSB_250].nsf;
719  read_att(file_1km.id, "Extract Pixel Offset", &intval);
720  spixl = intval * l1b[RSB_250].nsf;
721  read_att(file_1km.id, "Extract Line Offset", &intval);
722  sline = intval * l1b[RSB_250].ndets;
723  if (want_verbose) {
724  nline = l1file->nscan;
725  fprintf(stdout,
726  "File was generated from L1A extract of size %d x %d\n",
727  npixl, nline);
728  fprintf(stdout, " Pixels %d - %d\n", spixl + 1, spixl + npixl);
729  fprintf(stdout, " Lines %d - %d\n", sline + 1, sline + nline);
730  }
731  l1file->npix = npixl;
732  }
733 
734  /* Read reflectance to radiance conversion factors */
735  rdsensorinfo(l1file->sensorID, l1_input->evalmask, "Fobar", (void **) &Fobar);
736 
737  return SUCCESS;
738 }
739 
740 /***********************************************************************/
741 
744 static struct {
745  int32_t iscan;
746  int32_t nvals;
747  int32_t npix;
748  int32_t nbands;
750  /* read from input files */
751  double taisec;
752  int32_t mside;
753  double angles[3];
754  double mnorm[3];
755  float *lon;
756  float *lat;
757  float *hgt;
758  float *solz;
759  float *sola;
760  float *senz;
761  float *sena;
762  float *allbands;
763 } scan;
764 
770 int alloc_scan(const int32_t nvals) {
771 
772  /* set stack variables */
773  memset(&scan, 0, sizeof (scan));
774  scan.mside = -1;
775  scan.nvals = nvals;
776  scan.npix = nvals / l1b[RSB_250].ndets;
777  scan.nbands = sizeof (BANDLIST) / sizeof (BANDLIST[0]); // works on stack vars only
778 
779  /* allocate memory for dynamic variables */
780  TRYMEM(__FILE__, __LINE__, (scan.lon = malloc(nvals * sizeof (*scan.lon))));
781  TRYMEM(__FILE__, __LINE__, (scan.lat = malloc(nvals * sizeof (*scan.lat))));
782  TRYMEM(__FILE__, __LINE__, (scan.hgt = malloc(nvals * sizeof (*scan.hgt))));
783  TRYMEM(__FILE__, __LINE__, (scan.solz = malloc(nvals * sizeof (*scan.solz))));
784  TRYMEM(__FILE__, __LINE__, (scan.sola = malloc(nvals * sizeof (*scan.sola))));
785  TRYMEM(__FILE__, __LINE__, (scan.senz = malloc(nvals * sizeof (*scan.senz))));
786  TRYMEM(__FILE__, __LINE__, (scan.sena = malloc(nvals * sizeof (*scan.sena))));
787  TRYMEM(__FILE__, __LINE__,
788  (scan.allbands = malloc(scan.nbands * nvals * sizeof (*scan.allbands))));
789 
790  return SUCCESS;
791 }
792 
793 /***********************************************************************/
794 
795 typedef struct {
796  int16_t nrad;
797  int16_t nsf;
798  double *SI_val;
799  double *corr;
800 } sfcorr_table;
801 
802 void read_sfcorrtables(sfcorr_table *sfcorr, int32_t sensorID) {
803  char *dataroot;
804  char file[FILENAME_MAX] = "";
805  char sensorName[32];
806  FILE *fp = NULL;
807  char line[81];
808  int32_t iband, isf, nsf, irad, nrad;
809  size_t isds, ical;
810  att_struct *atts;
811  double scale, offset;
812  double radiance, SI, radiance2, SI2, rad_corr[4];
813 
814  /* calibration directory path */
815  if ((dataroot = getenv("OCDATAROOT")) == NULL) {
816  printf("OCDATAROOT environment variable is not defined.\n");
817  exit(1);
818  }
819 
820  /* for each band, */
821  /* get max possible subframes */
822  if (want_verbose)
823  printf("\nLoading subframe correction tables:\n");
824  for (iband = 0; iband < scan.nbands; iband++) {
825  nrad = 0;
826  isds = BANDLIST[iband].sdsindex;
827  if (isds == RSB_250)
828  nsf = 4;
829  else if (isds == RSB_500)
830  nsf = 2;
831  else
832  nsf = 1;
833 
834  /* if a hi-res band, */
835  if (nsf > 1) {
836 
837  /* get radiance factors */
838  atts = l1b[isds].sds.atts;
839  ical = BANDLIST[iband].bandindex;
840  scale = (double) ((float*) atts[RAD_SCALE].data)[ical];
841  offset = (double) ((float*) atts[RAD_OFFSET].data)[ical];
842 
843  strcpy(sensorName, sensorId2SensorName(sensorID));
844  lowcase(sensorName);
845 
846  /* open subframe correction file */
847  sprintf(file, "%s/%s/%s/cal/subframecor_%s_%d.dat",
848  dataroot,
849  sensorId2SensorDir(sensorID),
850  subsensorId2SubsensorDir(sensorId2SubsensorId(sensorID)),
851  sensorName,
852  BANDLIST[iband].wavelength);
853  if (want_verbose)
854  printf(" %s\n", file);
855  if ((fp = fopen(file, "r")) == NULL) {
856  fprintf(stderr,
857  "-E- %s line %d: unable to open %s for reading\n",
858  __FILE__, __LINE__, file);
859  exit(1);
860  }
861 
862  /* find number of radiance levels and allocate space */
863  while (fgets(line, 80, fp))
864  if (strncmp(line, "Number of radiance levels", 25) == 0) {
865  sscanf(line, "Number of radiance levels: %d", &nrad);
866  TRYMEM(__FILE__, __LINE__,
867  (sfcorr[iband].SI_val =
868  malloc(nrad * sizeof (*sfcorr[iband].SI_val))));
869  TRYMEM(__FILE__, __LINE__,
870  (sfcorr[iband].corr =
871  malloc(nrad * 4 * sizeof (*sfcorr[iband].corr))));
872  break;
873  }
874  while (fgets(line, 80, fp)) /* skip over rest of header */
875  if (strncmp(line, "/end_header", 11) == 0)
876  break;
877 
878  /* read radiances and correction factors */
879  for (irad = 0; irad < nrad; irad++) {
880  if (fgets(line, 80, fp) == NULL) {
881  fprintf(stderr,
882  "-E- %s line %d: File %s contains only %d data records\n",
883  __FILE__, __LINE__, file, irad);
884  exit(1);
885  }
886  sscanf(line, "%lf %lf %lf %lf %lf", &radiance,
887  &rad_corr[0], &rad_corr[1], &rad_corr[2], &rad_corr[3]);
888 
889  /* Scale reference radiances and corrections to SI
890  space so we don't have to scale each pixel later.
891  Also invert factor so we later multiply instead of
892  divide.
893  */
894  SI = radiance / scale + offset;
895  sfcorr[iband].SI_val[irad] = SI;
896  for (isf = 0; isf < 4; isf++) {
897  radiance2 = radiance / rad_corr[isf];
898  SI2 = radiance2 / scale + offset;
899  sfcorr[iband].corr[4 * irad + isf] = SI2 / SI;
900  }
901 
902  }
903  fclose(fp);
904 
905  } // hi-res band
906  sfcorr[iband].nrad = nrad;
907  sfcorr[iband].nsf = nsf;
908 
909  } // iband
910 
911  if (want_verbose)
912  printf("Subframe destriping corrections enabled.\n\n");
913 }
914 
915 void subframe_correction(const int32_t sensorID,
916  const int32_t iband,
917  const modis_sds mds,
918  double *data) {
919  /* Applies subframe destriping to an entire full-width scan at once,
920  so we don't need to correct for any extract offset. */
921  sfcorr_table lut;
922  int32_t i;
923  int32_t idet, iframe, isf, ip;
924  double x; // input value (scaled radiance)
925  double y; // interpolated subframe correction
926 
927  static sfcorr_table *sfcorr = NULL;
928 
929  /* Initialize subframe correction tables */
930  if (sfcorr == NULL) {
931  TRYMEM(__FILE__, __LINE__,
932  (sfcorr = malloc(scan.nbands * sizeof (*sfcorr))));
933  read_sfcorrtables(sfcorr, sensorID);
934  }
935  lut = sfcorr[iband];
936 
937  /* Destripe hi-res bands at native resolution only */
938  if ((lut.nsf > 1) && (lut.nrad > 0) && (lut.nsf == mds.nsf)) {
939  for (idet = 0; idet < mds.ndets; idet++)
940  for (iframe = 0; iframe < mds.nframes; iframe++)
941  for (isf = 0; isf < mds.nsf; isf++) {
942  ip = (idet * mds.nframes + iframe) * mds.nsf + isf;
943 
944  if (data[ip] < MIN_BAD_SI) {
945  x = data[ip];
946 
947  /* find closest lookup table indices */
948  /* TO DO: implement more efficient lookup */
949  for (i = 0; i < lut.nrad - 1; i++) {
950  if (x < lut.SI_val[i])
951  break;
952  }
953 
954  /* use end value if out-of-bounds */
955  if (i == 0 || i == (lut.nrad - 1))
956  y = lut.corr[4 * i + isf];
957 
958  /* otherwise, interpolate between two closest values */
959  else {
960  double x0 = lut.SI_val[i - 1];
961  double x1 = lut.SI_val[i];
962  double y0 = lut.corr[4 * (i - 1) + isf];
963  double y1 = lut.corr[4 * i + isf];
964  y = y0 + (x - x0) * (y1 - y0) / (x1 - x0);
965  }
966 
967  /* Apply the destriping correction */
968  data[ip] *= (float) y;
969  }
970  } // isf, iframe, idet loops
971  } // eligible band
972 }
973 
974 /***********************************************************************/
975 
976 typedef struct {
977  int32_t ndets;
978  int32_t *prev;
979  int32_t *next;
980 } fill_table;
981 
982 void init_detfill(fill_table *fill) {
983  int32_t ndets; // number of detectors
984  int32_t isds, iband, idet, i0, x0, x1;
985  char dead[EOSMETALEN] = "";
986 
987  /* Read "Dead Detector List" global attribute */
988  if (read_att(file_1km.id, "Dead Detector List", dead)) {
989  fprintf(stderr,
990  "-E- %s line %d: Error reading \"Dead Detector List\" from %s.\n",
991  __FILE__, __LINE__, file_1km.file);
992  exit(EXIT_FAILURE);
993  }
994 
995  /* Load values into dead detector table */
996  for (iband = 0; iband < scan.nbands; iband++) {
997  i0 = BANDLIST[iband].flagindex;
998 
999  /* Find number of detectors at native resolution */
1000  isds = BANDLIST[iband].sdsindex;
1001  if (isds == RSB_250)
1002  ndets = 40;
1003  else if (isds == RSB_500)
1004  ndets = 20;
1005  else
1006  ndets = 10;
1007  fill[iband].ndets = ndets;
1008 
1009  /* Allocate space for arrays */
1010  TRYMEM(__FILE__, __LINE__,
1011  (fill[iband].prev = malloc(ndets * sizeof (*fill->prev))));
1012  TRYMEM(__FILE__, __LINE__,
1013  (fill[iband].next = malloc(ndets * sizeof (*fill->next))));
1014 
1015  /* Find previous and next good detectors */
1016  /*
1017  These arrays hold indices for the nearest previous and next
1018  valid detector for each detector. Since detector arrays are
1019  small, opt for code clarity over efficiency.
1020 
1021  Special values:
1022  prev=idet=next: detector is alive
1023  prev = -1 : no previous good detector
1024  next = ndets: no next good detector
1025  */
1026  for (idet = 0; idet < ndets; idet++) {
1027  x0 = idet; // previous
1028  while ((x0 > -1) && dead[i0 + x0]) x0--;
1029  fill[iband].prev[idet] = x0;
1030  x1 = idet; // next
1031  while ((x1 < ndets) && dead[i0 + x1]) x1++;
1032  fill[iband].next[idet] = x1;
1033  }
1034 
1035  } // iband loop
1036 }
1037 
1038 void fill_dead_detectors(const int32_t iband, const modis_sds mds, double *data) {
1039 
1040  int32_t idet, ndets, ipix, npix;
1041  size_t nbytes; // bytes per line (detector)
1042 
1043  int32_t *prev, *next; // previous and next valid detector for each det
1044  int32_t x, x0, x1; // index into data array for current, prev, next det
1045  double y0, y1; // previous and next valid value for current detector
1046  double dx, dy; // change in det, value over valid interval
1047 
1048  static fill_table *fill = NULL;
1049 
1050  /* Initialize previous and next valid detector array table */
1051  if (fill == NULL) {
1052  TRYMEM(__FILE__, __LINE__, (fill = malloc(scan.nbands * sizeof (*fill))));
1053  init_detfill(fill);
1054  }
1055 
1056  /* Fill dead detectors at native resolution only */
1057  ndets = fill[iband].ndets;
1058  if (ndets == mds.ndets) {
1059 
1060  prev = fill[iband].prev;
1061  next = fill[iband].next;
1062  npix = mds.nframes * mds.nsf; // # values per line (detector)
1063  nbytes = npix * sizeof (*data);
1064 
1065  /* Substitute data values as appropriate for each detector: */
1066  for (idet = 0; idet < ndets; idet++) {
1067  if (prev[idet] != next[idet]) { // continue for bad dets only
1068 
1069  /* Between two good detectors: interpolate */
1070  if ((prev[idet] > -1) && (next[idet] < ndets)) {
1071  dx = next[idet] - prev[idet];
1072 
1073  for (ipix = 0; ipix < npix; ipix++) {
1074  x = idet * npix + ipix;
1075  x0 = prev[idet] * npix + ipix;
1076  x1 = next[idet] * npix + ipix;
1077  y0 = data[x0];
1078  y1 = data[x1];
1079 
1080  /* if both values are valid, interpolate */
1081  if ((y0 < MIN_BAD_SI) && (y1 < MIN_BAD_SI)) {
1082  dy = y1 - y0;
1083  data[x] = y0 + (idet - prev[idet]) * (dy / dx);
1084  }
1085  /* otherwise copy adjacent valid value */
1086  else if (y0 < MIN_BAD_SI) data[x] = data[x0];
1087  else if (y1 < MIN_BAD_SI) data[x] = data[x1];
1088 
1089  /* if none, copy maximum adjacent flag */
1090  else
1091  data[x] = MAX(data[x0], data[x1]);
1092  }
1093  }
1094  /* No previous good detector: copy next value or flag */
1095  else if ((prev[idet] == -1) && (next[idet] < ndets))
1096  memcpy(&data[idet * npix], &data[next[idet] * npix], nbytes);
1097 
1098  /* No next good detector: copy previous value or flag */
1099  else if ((prev[idet] > -1) && (next[idet] == ndets))
1100  memcpy(&data[idet * npix], &data[prev[idet] * npix], nbytes);
1101 
1102  /* Default case: no good detectors = no change. */
1103 
1104  } // dead detector
1105  } // idet loop
1106  } // band is at native resolution
1107 }
1108 
1109 /***********************************************************************/
1110 
1111 double interp_bilin(double val[4], double dx, double dy) {
1112  return (1 - dy) * (val[0] * (1 - dx) + val[1] * dx)
1113  + dy * (val[2] * (1 - dx) + val[3] * dx);
1114 }
1115 
1116 double interp_modis_EVdata(double val[4], double dx, double dy) {
1117  int32_t i, ngood = 0;
1118  double sum_val = 0;
1119  double max_val = 0;
1120  double result;
1121 
1122  /* interpolation: account for flagged values */
1123  for (i = 0; i < 4; i++) {
1124  if (val[i] < MIN_BAD_SI) {
1125  sum_val = sum_val + val[i];
1126  ngood++;
1127  }
1128  if (val[i] > max_val)
1129  max_val = val[i];
1130  }
1131 
1132  /* If we have four non-flagged samples, bilinearly interpolate. */
1133  /* If all samples are flagged, set to the max flag value. */
1134  /* If we have 1-3 samples, use average of unflagged values. */
1135  switch (ngood) {
1136  case 4:
1137  result = interp_bilin(val, dx, dy);
1138  break;
1139  case 0:
1140  result = max_val;
1141  break;
1142  default:
1143  /* TO DO: should be WEIGHTED avg of good values */
1144  result = sum_val / ngood;
1145  break;
1146  }
1147  return result;
1148 }
1149 
1150 double interp_modis_Longitude(double val[4], double dx, double dy) {
1151  int32_t i;
1152  double min_val, max_val;
1153 
1154  /* longitude interpolation: account for dateline */
1155  min_val = max_val = val[0];
1156  for (i = 1; i < 4; i++) {
1157  if (val[i] < min_val)
1158  min_val = val[i];
1159  else if (val[i] > max_val)
1160  max_val = val[i];
1161  }
1162  if (max_val > 180.0 + min_val) {
1163  for (i = 0; i < 4; i++) {
1164  if (val[i] < 0.0)
1165  val[i] += 360.0;
1166  }
1167  }
1168  return interp_bilin(val, dx, dy);
1169 }
1170 
1171 double (*interp_modis_var)(double val[4], double dx, double dy);
1172 
1173 int modis_interp(const modis_sds s,
1174  const int32_t sdsband,
1175  const int32_t newres,
1176  double *olddata,
1177  double *newdata) {
1178 
1179  int32_t ny = s.ndets;
1180  int32_t nx = s.nframes * s.nsf;
1181  int32_t nvals = ny * nx;
1182  int32_t oldres = s.resolution;
1183 
1184  int32_t nxnew, nynew;
1185  int32_t ixnew, iynew;
1186 
1187  int32_t x0, y0, x1, y1; // 4 points in old array
1188  double fx, fy;
1189  double dx, dy;
1190  double val[4];
1191  int32_t i;
1192  double factor, xoffset, yoffset;
1193  int32_t maxi = 2; /* 1 duplicates end pixel; 2 extrapolates past end */
1194 
1195  int32_t sdsoffset = nvals * sdsband; // offset within SDS
1196  short alloc_old;
1197 
1198  /* Cast input data to double precision */
1199  alloc_old = (olddata == NULL);
1200  if (alloc_old) {
1201  TRYMEM(__FILE__, __LINE__, (olddata = malloc(nvals * sizeof (*olddata))));
1202 
1203  switch (s.sds.ntype) {
1204  case DFNT_FLOAT32:
1205  LOOPI olddata[i] = (double) ((float*) s.sds.data)[sdsoffset + i];
1206  break;
1207  case DFNT_INT16:
1208  LOOPI olddata[i] = (double) ((int16_t*) s.sds.data)[sdsoffset + i];
1209  break;
1210  case DFNT_UINT16:
1211  LOOPI olddata[i] = (double) ((uint16_t*) s.sds.data)[sdsoffset + i];
1212  break;
1213  default:
1214  fprintf(stderr,
1215  "-E- %s line %d: Cannot interpolate type code: %d\n",
1216  __FILE__, __LINE__, s.sds.ntype);
1217  exit(EXIT_FAILURE);
1218  break;
1219  }
1220  } // alloc_old
1221 
1222  /* Set interp function pointer according to SDS */
1223  if (strncmp(s.sds.name, "EV_", 3) == 0)
1224  interp_modis_var = interp_modis_EVdata;
1225  else if (strcmp(s.sds.name, "Longitude") == 0)
1226  interp_modis_var = interp_modis_Longitude;
1227  else
1228  interp_modis_var = interp_bilin;
1229 
1230  /* Calculate offset between different pixel resolutions */
1231  factor = (double) newres / (double) oldres;
1232  xoffset = 0; // frame centers are aligned
1233  yoffset = 0.5 * (1.0 - factor); // detector edges are aligned
1234  /* yoffset = 0; // for testing only */
1235 
1236  /* Size of new array */
1237  nxnew = nx * oldres / newres;
1238  nynew = ny * oldres / newres;
1239 
1240  /* Calculate index into original data */
1241  for (iynew = 0; iynew < nynew; iynew++) {
1242  fy = factor * iynew - yoffset;
1243  y0 = MIN(MAX((int32_t) floor(fy), 0), ny - maxi);
1244  y1 = MIN(y0 + 1, ny - 1);
1245  dy = fy - (double) y0;
1246 
1247  for (ixnew = 0; ixnew < nxnew; ixnew++) {
1248  fx = factor * ixnew - xoffset;
1249  x0 = MIN(MAX((int32_t) floor(fx), 0), nx - maxi);
1250  x1 = MIN(x0 + 1, nx - 1);
1251  dx = fx - (double) x0;
1252 
1253  /* Do bilinear interpolation */
1254  val[0] = olddata[nx * y0 + x0];
1255  val[1] = olddata[nx * y0 + x1];
1256  val[2] = olddata[nx * y1 + x0];
1257  val[3] = olddata[nx * y1 + x1];
1258  newdata[nxnew * iynew + ixnew] = interp_modis_var(val, dx, dy);
1259  }
1260  }
1261 
1262  if (alloc_old)
1263  free(olddata);
1264 
1265  return SUCCESS;
1266 }
1267 
1268 /***********************************************************************/
1269 
1275 int load_modis_scan(const int32_t iscan,
1276  const int32_t resolution,
1277  int32_t sensorID) {
1278  modis_sds s;
1279  int32_t isds, iband;
1280  int32_t i, i1, i2, nold;
1281  static int32_t nvals = 0; /* pixels per variable at output resolution */
1282  double *olddata = NULL; /* at native resolution */
1283  double *newdata = NULL; /* at output resolution */
1284 
1285  /* Allocate space to hold all data for one scan */
1286  if (nvals == 0) {
1287  nvals = l1b[RSB_250].ndets * l1b[RSB_250].nframes * l1b[RSB_250].nsf;
1288  alloc_scan(nvals);
1289  }
1290 
1291  /* Allocate space to hold interpolation buffer for one variable */
1292  TRYMEM(__FILE__, __LINE__, (newdata = malloc(nvals * sizeof (*newdata))));
1293 
1294  /* ----- Geolocation values ----- */
1295  for (isds = 0; isds < GEO_NUM_SDS; isds++)
1296  read_sds_1scan(&geo[isds].sds,
1297  iscan,
1298  geo[isds].scandim,
1299  geo[isds].ndets);
1300 
1301  if (geo[GEO_LON].resolution == resolution) { /* native resolution */
1302  /* memcpy works only when src and dst have same data type */
1303  memcpy(scan.lon, geo[GEO_LON].sds.data, nvals * sizeof (float));
1304  memcpy(scan.lat, geo[GEO_LAT].sds.data, nvals * sizeof (float));
1305 
1306  /* must copy each element when datatypes differ */
1307  LOOPI{
1308  scan.hgt[i] = (float) ((int16_t *) geo[GEO_HGT].sds.data)[i];
1309  scan.solz[i] = (float) ((int16_t *) geo[GEO_SOLZ].sds.data)[i] * 0.01;
1310  scan.sola[i] = (float) ((int16_t *) geo[GEO_SOLA].sds.data)[i] * 0.01;
1311  scan.senz[i] = (float) ((int16_t *) geo[GEO_SENZ].sds.data)[i] * 0.01;
1312  scan.sena[i] = (float) ((int16_t *) geo[GEO_SENA].sds.data)[i] * 0.01;}
1313 
1314  } else { /* must interpolate */
1315  modis_interp(geo[GEO_LON], 0, resolution, NULL, newdata);
1316  LOOPI scan.lon[i] = (float) newdata[i];
1317  modis_interp(geo[GEO_LAT], 0, resolution, NULL, newdata);
1318  LOOPI scan.lat[i] = (float) newdata[i];
1319  modis_interp(geo[GEO_HGT], 0, resolution, NULL, newdata);
1320  LOOPI scan.hgt[i] = (float) newdata[i];
1321  modis_interp(geo[GEO_SOLZ], 0, resolution, NULL, newdata);
1322  LOOPI scan.solz[i] = (float) newdata[i] * 0.01;
1323  modis_interp(geo[GEO_SOLA], 0, resolution, NULL, newdata);
1324  LOOPI scan.sola[i] = (float) newdata[i] * 0.01;
1325  modis_interp(geo[GEO_SENZ], 0, resolution, NULL, newdata);
1326  LOOPI scan.senz[i] = (float) newdata[i] * 0.01;
1327  modis_interp(geo[GEO_SENA], 0, resolution, NULL, newdata);
1328  LOOPI scan.sena[i] = (float) newdata[i] * 0.01;
1329  }
1330 
1331  /* ----- Radiometric values ----- */
1332  for (isds = 0; isds < L1B_NUM_SDS; isds++)
1333  read_sds_1scan(&l1b[isds].sds,
1334  iscan,
1335  l1b[isds].scandim,
1336  l1b[isds].ndets);
1337 
1338  for (iband = 0; iband < scan.nbands; iband++) {
1339  s = l1b[BANDLIST[iband].sdsindex];
1340  nold = s.ndets * s.nframes * s.nsf;
1341  i1 = nold * BANDLIST[iband].bandindex; // offset within SDS
1342  i2 = nvals * iband; // offset within scan.allbands
1343 
1344  /* Cast SI values to double precision, in native spatial resolution */
1345  TRYMEM(__FILE__, __LINE__, (olddata = malloc(nold * sizeof (*olddata))));
1346  for (i = 0; i < nold; i++)
1347  olddata[i] = (double) ((uint16_t*) s.sds.data)[i1 + i];
1348 
1349  /* Apply subframe correction to destripe BEFORE interpolation */
1350  if ((sensorID == MODISA) && (resolution < 1000))
1351  subframe_correction(sensorID, iband, s, olddata);
1352 
1353  /* Fill any dead detectors with neighboring values */
1354  fill_dead_detectors(iband, s, olddata);
1355 
1356  /* Interpolate to higher resolution, as needed */
1357  if (s.resolution != resolution) {
1358  modis_interp(l1b[BANDLIST[iband].sdsindex],
1359  BANDLIST[iband].bandindex,
1360  resolution, olddata, newdata);
1361 
1362  /* Copy adjusted data into scan structure */
1363  LOOPI scan.allbands[i2++] = (float) newdata[i]; // interpolated
1364  } else
1365  LOOPI scan.allbands[i2++] = (float) olddata[i]; // native resolution
1366 
1367  free(olddata);
1368  }
1369 
1370  /* ----- Store preloaded values for later ----- */
1371  scan.iscan = iscan;
1372  scan.taisec = ((double*) ref[GEO_TAISEC].sds.data)[iscan];
1373  scan.mside = (int32_t) ((int16_t*) ref[GEO_MSIDE].sds.data)[iscan];
1374  for (i = 0; i < 3; i++) {
1375  scan.angles[i] = ((double*) ref[GEO_ANGLES].sds.data)[3 * iscan + i];
1376  scan.mnorm[i] = ((double*) ref[GEO_MNORM].sds.data)[3 * (3 * iscan + i)];
1377  }
1378  free(newdata);
1379  return SUCCESS;
1380 }
1381 
1388 int set_l1rec_scanvals(const int32_t iscan, l1str *l1rec) {
1389  double esdist;
1390  static double lastTAIsec = -1;
1391  int16_t badatt, badmside;
1392  double TAIsec, dsec;
1393  int16_t year, day;
1394  int32_t i;
1395 
1396  /* Mirror Side */
1397  l1rec->mside = scan.mside;
1398  badmside = (l1rec->mside < 0);
1399  if (badmside) /* Fix mirror side if set to -1 */
1400  l1rec->mside = 0;
1401 
1402  /* Check for non-nominal roll, pitch, or yaw */
1403  badatt =
1404  (fabs(scan.angles[0]) > MAX_ATTERR) ||
1405  (fabs(scan.angles[1]) > MAX_ATTERR) ||
1406  (fabs(scan.angles[2]) > MAX_ATTERR);
1407 
1408  /* Flag each pixel of scan for geolocation errors */
1409  for (i = 0; i < l1rec->npix; i++) {
1410  l1rec->pixnum[i] = spixl + i;
1411  l1rec->navfail[i] = badmside;
1412  l1rec->navwarn[i] = badatt;
1413  }
1414 
1415  /* Scan Start Time */
1416  TAIsec = scan.taisec;
1417  if (TAIsec < 0) { /* Workaround for scan time errors */
1418  fprintf(stderr, "-W- %s: Bad time for scan %d; ", __FILE__, iscan);
1419  if (iscan == 0) {
1420  fprintf(stderr, "using granule start time.\n");
1421  TAIsec = unix_to_tai93(file_geo.starttime);
1422  } else {
1423  fprintf(stderr, "extrapolating from last good scan.\n");
1424  TAIsec = lastTAIsec + SCAN_TIME_INTERVAL;
1425  }
1426  }
1427  l1rec->scantime = tai93_to_unix(TAIsec);
1428  lastTAIsec = TAIsec;
1429 
1430  /* Earth-sun distance correction for this scan */
1431  unix2yds(l1rec->scantime, &year, &day, &dsec);
1432  int32_t msec = (int32_t) (dsec * 1000.0);
1433  int32_t yr = (int32_t) year;
1434  int32_t dy = (int32_t) day;
1435  esdist = esdist_(&yr, &dy, &msec);
1436  l1rec->fsol = pow(1.0 / esdist, 2);
1437 
1438  return SUCCESS;
1439 }
1440 
1448 int readl1_modis(filehandle *l1file,
1449  const int32_t line,
1450  l1str *l1rec, int lonlat) {
1451  static int32_t prev_scan = -1;
1452  int32_t iband, iscan, idet, ndets;
1453  size_t ip, ipb;
1454  size_t detoffset, dataptr;
1455  size_t nbytes;
1456  size_t irsb, iteb, ical;
1457  att_struct *atts;
1458  float scale, offset, SI_val, value;
1459 
1460  /* Scan and detector numbers for this line */
1461  ndets = l1b[RSB_250].ndets;
1462  iscan = line / ndets;
1463  idet = line % ndets;
1464 
1465  /* Load all values for next scan as needed */
1466  if (iscan != prev_scan) {
1467  prev_scan = iscan;
1468  load_modis_scan(iscan, l1_input->resolution, l1file->sensorID);
1469  }
1470 
1471  /* The l1rec structure is reinitialized */
1472  /* each line, so always reload basic info. */
1473  l1rec->npix = l1file->npix;
1474  l1rec->detnum = (int32_t) idet;
1475  set_l1rec_scanvals(iscan, l1rec);
1476  unix_to_tai93(file_geo.starttime);
1477 
1478  /* "scan" structure data pointer offset due */
1479  /* to line number and extract start pixel */
1480  detoffset = idet * scan.npix + spixl;
1481 
1482  /* ----- Geolocation values ----- */
1483  nbytes = l1rec->npix * sizeof (float);
1484  memcpy(l1rec->lon, &scan.lon[detoffset], nbytes);
1485  memcpy(l1rec->lat, &scan.lat[detoffset], nbytes);
1486  memcpy(l1rec->solz, &scan.solz[detoffset], nbytes);
1487 
1488  // lonlat only needs time, lon, lat and solar zenith
1489  if (lonlat)
1490  return LIFE_IS_GOOD;
1491 
1492  memcpy(l1rec->height, &scan.hgt[detoffset], nbytes);
1493  memcpy(l1rec->solz, &scan.solz[detoffset], nbytes);
1494  memcpy(l1rec->sola, &scan.sola[detoffset], nbytes);
1495  memcpy(l1rec->senz, &scan.senz[detoffset], nbytes);
1496  memcpy(l1rec->sena, &scan.sena[detoffset], nbytes);
1497 
1498  /* Don't need to check for navigation errors
1499  - it's done by readl1() */
1500 
1501  /* Compute polarization frame rotation angles */
1502  compute_alpha(l1rec->lon, l1rec->lat,
1503  l1rec->senz, l1rec->sena,
1504  scan.mnorm, l1rec->npix, l1rec->alpha);
1505 
1506  /* ----- Radiometric values ----- */
1507  irsb = iteb = 0;
1508  for (iband = 0; iband < scan.nbands; iband++) {
1509 
1510  /* calibration info */
1511  atts = l1b[BANDLIST[iband].sdsindex].sds.atts;
1512  ical = BANDLIST[iband].bandindex;
1513 
1514  /* offset within scan.allbands */
1515  dataptr = iband * scan.nvals + detoffset;
1516 
1517  /*** Thermal Emissive Bands ***/
1518  if (BANDLIST[iband].sdsindex == TEB_1KM) {
1519  scale = ((float*) atts[RAD_SCALE].data)[ical];
1520  offset = ((float*) atts[RAD_OFFSET].data)[ical];
1521 
1522  for (ip = 0; ip < l1rec->npix; ip++) {
1523  value = 0.0; /* init to fill value */
1524 
1525  /* Apply calibration for valid values */
1526  SI_val = scan.allbands[dataptr + ip];
1527  if (SI_val < MIN_BAD_SI) {
1528 
1529  /* Convert radiance from Watts/m^2/micrometer/steradian */
1530  /* to mW/cm2/um/sr, as expected by l2gen */
1531  value = (SI_val - offset) * scale / 10.0;
1532 
1533  /* Terra adjustments */
1534  if (l1rec->l1file->sensorID == MODIST) {
1535  value *= (1.0 - radoff[iteb][idet / (ndets / 10)]);
1536  if (l1rec->mside == 1)
1537  value *= (1.0 + mfact[iteb]);
1538  } /* Aqua adjustments */
1539  else {
1540  /* channel 22 detector balancing */
1541  if (ical == 1)
1542  value /= l1_input->ch22detcor[line % 10];
1543  /* channel 23 detector balancing */
1544  if (ical == 2)
1545  value /= l1_input->ch23detcor[line % 10];
1546  }
1547  }
1548 
1549  /* populate l1rec->Ltir */
1550  ipb = ip * l1rec->l1file->nbandsir + iteb; /* band varies fastest */
1551  l1rec->Ltir[ipb] = value;
1552 
1553  } // ip
1554  iteb++;
1555  }// TEB
1556 
1557  /*** Cirrus Band ***/
1558  else if (BANDLIST[iband].sdsindex == CIR_1KM) {
1559  scale = ((float*) atts[REFL_SCALE].data)[ical];
1560  offset = ((float*) atts[REFL_OFFSET].data)[ical];
1561 
1562  for (ip = 0; ip < l1rec->npix; ip++) {
1563  value = BAD_FLT; /* init to fill value */
1564 
1565  /* Apply calibration for valid values */
1566  SI_val = scan.allbands[dataptr + ip];
1567  if (SI_val < MIN_BAD_SI)
1568 
1569  /* Normalize reflectance by solar zenith angle */
1570  value = (SI_val - offset) * scale
1571  / cos(l1rec->solz[ip] / RADEG);
1572 
1573  /* populate l1rec->rho_cirrus */
1574  l1rec->rho_cirrus[ip] = value;
1575 
1576  } // ip
1577  }// CIR
1578 
1579  /*** Reflective Solar Bands ***/
1580  else {
1581  l1rec->Fo[irsb] = Fobar[irsb] * l1rec->fsol;
1582  scale = ((float*) atts[REFL_SCALE].data)[ical];
1583  offset = ((float*) atts[REFL_OFFSET].data)[ical];
1584 
1585  for (ip = 0; ip < l1rec->npix; ip++) {
1586  value = BAD_FLT; /* init to fill value */
1587 
1588  /* skip night data for visible bands */
1589  if (SOLZNIGHT < l1rec->solz[ip])
1590  continue;
1591 
1592  /* Apply calibration for valid values */
1593  SI_val = scan.allbands[dataptr + ip];
1594  if (SI_val < MIN_BAD_SI)
1595 
1596  /* Convert from reflectance to radiance */
1597  value = (SI_val - offset) * scale * l1rec->Fo[irsb] / PI;
1598 
1599  /* Flag any sentinel values */
1600  else
1601  switch ((int32_t) SI_val) {
1602  case 65527: /* not in earth view: SECTOR_ROTATION_SI */
1603  l1rec->navfail[ip] = 1;
1604  break;
1605  case 65529: /* saturation: TEB_OR_RSB_GT_MAX_SI */
1606  case 65533: /* saturation: SATURATED_DETECTOR_SI */
1607  value = 1000.0;
1608  if (irsb < 13)
1609  l1rec->hilt[ip] = 1;
1610  break;
1611  default: /* all others */
1612  break;
1613  }
1614 
1615  /* populate l1rec->Lt */
1616  ipb = ip * l1rec->l1file->nbands + irsb; /* band varies fastest */
1617  l1rec->Lt[ipb] = value;
1618 
1619  } // ip
1620  irsb++;
1621  } // RSB
1622 
1623  } // iband
1624 
1625  /* Convert IR bands to brightness temperature */
1626  radiance2bt(l1rec, l1_input->resolution);
1627 
1628  return SUCCESS;
1629 }
1630 
1631 /***********************************************************************/
1632 
1637 int closel1_modis() {
1638  int32_t i;
1639 
1640  /* Free memory allocated for global variables */
1641  for (i = 0; i < REF_NUM_SDS; i++)
1642  free_sds(&(ref[i].sds));
1643  for (i = 0; i < GEO_NUM_SDS; i++)
1644  free_sds(&(geo[i].sds));
1645  for (i = 0; i < L1B_NUM_SDS; i++)
1646  free_sds(&(l1b[i].sds));
1647 
1648  /* Free any memory allocated for input files */
1649  if (file_geo.sds) {
1650  for (i = 0; i < file_geo.n_sds; i++)
1651  free_sds(&(file_geo.sds[i]));
1652  free(file_geo.sds);
1653  }
1654  if (file_1km.sds) {
1655  for (i = 0; i < file_1km.n_sds; i++)
1656  free_sds(&(file_1km.sds[i]));
1657  free(file_1km.sds);
1658  }
1659  if (file_hkm.sds) {
1660  for (i = 0; i < file_hkm.n_sds; i++)
1661  free_sds(&(file_hkm.sds[i]));
1662  free(file_hkm.sds);
1663  }
1664  if (file_qkm.sds) {
1665  for (i = 0; i < file_qkm.n_sds; i++)
1666  free_sds(&(file_qkm.sds[i]));
1667  free(file_qkm.sds);
1668  }
1669 
1670  /* End HDF access */
1671  if (file_geo.id) SDend(file_geo.id);
1672  if (file_1km.id) SDend(file_1km.id);
1673  if (file_hkm.id) SDend(file_hkm.id);
1674  if (file_qkm.id) SDend(file_qkm.id);
1675  file_geo.id = file_1km.id = file_hkm.id = file_qkm.id = 0;
1676 
1677  /* TO DO: free scan structure */
1678 
1679  return SUCCESS;
1680 }
1681 
1682 /***********************************************************************/
1683 
1684 int readl1_lonlat_modis(filehandle *l1file, int32_t line, l1str *l1rec) {
1685  static float *lonscan = NULL; /* storage buffer arrays */
1686  static float *latscan = NULL; /* (at output resolution) */
1687 
1688  int32_t ndets = l1b[RSB_250].ndets; /* output resolution dimensions*/
1689  int32_t npixl = l1b[RSB_250].nframes * l1b[RSB_250].nsf;
1690  int32_t nvals = ndets * npixl;
1691  int32_t resolution = l1_input->resolution;
1692 
1693  double *newdata = NULL; /* interpolation buffer */
1694  static int32_t prev_scan = -1;
1695  int32_t i, iscan, idet;
1696  size_t detoffset, nbytes;
1697 
1698  /* Allocate space to hold lon/lat data for one scan at output resolution */
1699  if (lonscan == NULL) {
1700  TRYMEM(__FILE__, __LINE__, (lonscan = malloc(nvals * sizeof (*lonscan))));
1701  TRYMEM(__FILE__, __LINE__, (latscan = malloc(nvals * sizeof (*latscan))));
1702  }
1703 
1704  /* Scan and detector indices for this line */
1705  iscan = line / ndets;
1706  idet = line % ndets;
1707 
1708  /* Load next scan as needed */
1709  if (iscan != prev_scan) {
1710  prev_scan = iscan;
1711 
1712  /* Read values */
1713  read_sds_1scan(&geo[GEO_LON].sds,
1714  iscan,
1715  geo[GEO_LON].scandim,
1716  geo[GEO_LON].ndets);
1717  read_sds_1scan(&geo[GEO_LAT].sds,
1718  iscan,
1719  geo[GEO_LAT].scandim,
1720  geo[GEO_LAT].ndets);
1721 
1722  /* If native resolution, copy values into storage buffers */
1723  if (geo[GEO_LON].resolution == resolution) {
1724  memcpy(lonscan, geo[GEO_LON].sds.data, nvals * sizeof (*lonscan));
1725  memcpy(latscan, geo[GEO_LAT].sds.data, nvals * sizeof (*latscan));
1726  }
1727  /* Otherwise, interpolate first */
1728  else {
1729  TRYMEM(__FILE__, __LINE__,
1730  (newdata = malloc(nvals * sizeof (*newdata))));
1731  modis_interp(geo[GEO_LON], 0, resolution, NULL, newdata);
1732  LOOPI lonscan[i] = (float) newdata[i];
1733  modis_interp(geo[GEO_LAT], 0, resolution, NULL, newdata);
1734  LOOPI latscan[i] = (float) newdata[i];
1735  free(newdata);
1736  }
1737 
1738  } // loaded new scan
1739 
1740  /* Copy specified line from storage buffers into l1rec structure */
1741  detoffset = idet * npixl + spixl;
1742  nbytes = l1rec->npix * sizeof (float);
1743  memcpy(l1rec->lon, &lonscan[detoffset], nbytes);
1744  memcpy(l1rec->lat, &latscan[detoffset], nbytes);
1745 
1746  return SUCCESS;
1747 }
GEO_SDS
GEO_SDS
Definition: l1_hmodis_hdf.c:56
GEO_SCALE
@ GEO_SCALE
Definition: l1_modis.c:78
l1file
int32 l1file(int32 sdfid, int32 *nsamp, int32 *nscans, int16 *dtynum)
Definition: l1stat_chk.c:586
index
an array had not been initialized Several spelling and grammar corrections were which is read from the appropriate MCF the above metadata values were hard coded A problem calculating the average background DN for SWIR bands when the moon is in the space view port was corrected The new algorithm used to calculate the average background DN for all reflective bands when the moon is in the space view port is now the same as the algorithm employed by the thermal bands For non SWIR changes in the averages are typically less than Also for non SWIR the black body DNs remain a backup in case the SV DNs are not available For SWIR the changes in computed averages were larger because the old which used the black body suffered from contamination by the micron leak As a consequence of the if SV DNs are not available for the SWIR the EV pixels will not be the granule time is used to identify the appropriate tables within the set given for one LUT the first two or last two tables respectively will be used for the interpolation If there is only one LUT in the set of it will be treated as a constant LUT The manner in which Earth View data is checked for saturation was changed Previously the raw Earth View DNs and Space View DNs were checked against the lookup table values contained in the table dn_sat The change made is to check the raw Earth and Space View DNs to be sure they are less than the maximum saturation value and to check the Space View subtracted Earth View dns against a set of values contained in the new lookup table dn_sat_ev The metadata configuration and ASSOCIATEDINSTRUMENTSHORTNAME from the MOD02HKM product The same metatdata with extensions and were removed from the MOD021KM and MOD02OBC products ASSOCIATEDSENSORSHORTNAME was set to MODIS in all products These changes are reflected in new File Specification which users may consult for exact the pow functions were eliminated in Emissive_Cal and Emissive bands replaced by more efficient code Other calculations throughout the code were also made more efficient Aside from a few round off there was no difference to the product The CPU time decreased by about for a day case and for a night case A minor bug in calculating the uncertainty index for emissive bands was corrected The frame index(0-based) was previously being used the frame number(1-based) should have been used. There were only a few minor changes to the uncertainty index(maximum of 1 digit). 3. Some inefficient arrays(Sigma_RVS_norm_sq) were eliminated and some code lines in Preprocess_L1A_Data were moved into Process_OBCEng_Emiss. There were no changes to the product. Required RAM was reduced by 20 MB. Now
MAX
#define MAX(A, B)
Definition: swl0_utils.h:25
value
int32 value
Definition: Granule.c:1235
sensorId2SensorDir
const char * sensorId2SensorDir(int sensorId)
Definition: sensorInfo.c:315
MIN
#define MIN(x, y)
Definition: rice.h:169
modis_sds::nscans
int32_t nscans
Definition: l1_hmodis_hdf.c:148
GEO_ANGLES
@ GEO_ANGLES
Definition: l1_modis.c:43
hgt
float * hgt
Definition: l1_modis.c:757
bandinfo::sdsindex
int16_t sdsindex
Definition: l1_hmodis_hdf.c:568
modis_sds::resolution
int16_t resolution
Definition: l1_hmodis_hdf.c:146
SUCCESS
#define SUCCESS
Definition: ObpgReadGrid.h:15
sola
float * sola
Definition: l1_modis.c:759
GEO_LON
@ GEO_LON
Definition: l1_modis.c:57
sfcorr_table
Definition: l1_hmodis_hdf.c:795
j
int j
Definition: decode_rs.h:73
MAX_ATTERR
#define MAX_ATTERR
Definition: l1_modis.c:537
day
int32_t day
Definition: atrem_corl1v2.h:308
allbands
float * allbands
Definition: l1_modis.c:762
sfcorr_table::nrad
int16_t nrad
Definition: l1_hmodis_hdf.c:796
status
int status
Definition: l1_czcs_hdf.c:32
fill_table
Definition: l1_hmodis_hdf.c:976
modis_sds::ndets
int32_t ndets
Definition: l1_hmodis_hdf.c:149
lowcase
char * lowcase(char *instr)
Definition: lowcase.c:10
REFL_SCALE
@ REFL_SCALE
Definition: l1_modis.c:113
TRYMEM
#define TRYMEM(file, line, memstat)
Definition: l1_modis.c:23
SCAN_TIME_INTERVAL
#define SCAN_TIME_INTERVAL
Definition: l1_modis.c:20
modis_file::resolution
int16_t resolution
Definition: l1_hmodis_hdf.c:134
lon
float * lon
Definition: l1_modis.c:755
modis_file::file
char file[FILENAME_MAX]
Definition: l1_hmodis_hdf.c:131
viirs_ler_luts::ref
real ref
Definition: viirs_ler_luts.f95:60
lat
float * lat
Definition: l1_modis.c:756
FAIL
#define FAIL
Definition: ObpgReadGrid.h:18
modis.modis_utils.geofile
geofile
Definition: modis_utils.py:368
angles
double angles[3]
Definition: l1_modis.c:753
NULL
#define NULL
Definition: decode_rs.h:63
RSB_250
@ RSB_250
Definition: l1_modis.c:87
senz
float * senz
Definition: l1_modis.c:760
scan
MOD_PR01 Production producing one five minute granule of output data in each run It can be configured to produce as many as three five minute granules per run Each execution with one construction record and one date file for each dataset In normal these are created by which splits them out of the hour datasets For LANCE they are created by which merges all session MODIS L0 datasets overlapping the requested time and extracts from the merged data those packets which fall within that time period Each scan of data is stored in the L1A granule that covers the start time of that scan
Definition: MOD_PR01_pr.txt:19
resolution
this program makes no use of any feature of the SDP Toolkit that could generate such a then geolocation is calculated at that resolution
Definition: HISTORY.txt:188
modis_sds::scandim
int32_t scandim
Definition: l1_hmodis_hdf.c:152
l1rec
read l1rec
Definition: HOWTO_Add_a_sensor.txt:72
openl1_modis
int openl1_modis(filehandle *l1file)
Definition: l1_modis.c:614
RSB_1KM
@ RSB_1KM
Definition: l1_modis.c:89
L1B_COEFFS
L1B_COEFFS
Definition: l1_hmodis_hdf.c:112
alloc_scan
int alloc_scan(const int32_t nvals)
Definition: l1_modis.c:770
modis_sds::sds
sds_struct sds
Definition: l1_hmodis_hdf.c:145
wavelength
PARAM_TYPE_NONE Default value No parameter is buried in the product name name_prefix is case insensitive string compared to the product name PARAM_TYPE_VIS_WAVE The visible wavelength bands from the sensor are buried in the product name The product name is compared by appending wavelength
Definition: HOWTO_Add_a_product.txt:24
init_sds_byname
int init_sds_byname(int32_t fileid, const char *sdsname, sds_struct *sds)
Definition: hdf4_utils.c:92
msec
int32 * msec
Definition: l1_czcs_hdf.c:31
MODIST
#define MODIST
Definition: sensorDefs.h:18
esdist_
double esdist_(int32_t *year, int32_t *day, int32_t *msec)
sfcorr_table::SI_val
double * SI_val
Definition: l1_hmodis_hdf.c:798
REF_NUM_SDS
@ REF_NUM_SDS
Definition: l1_modis.c:45
LOOPI
#define LOOPI
Definition: l1_modis.c:538
load_modis_scan
int load_modis_scan(const int32_t iscan, const int32_t resolution, int32_t sensorID)
Definition: l1_modis.c:1275
sensorId2SubsensorId
int sensorId2SubsensorId(int sensorId)
Definition: sensorInfo.c:438
modis_file::n_sds
int32_t n_sds
Definition: l1_hmodis_hdf.c:137
interp_bilin
double interp_bilin(double val[4], double dx, double dy)
Definition: l1_modis.c:1111
file
no change in intended resolving MODur00064 Corrected handling of bad ephemeris attitude resolving resolving GSFcd00179 Corrected handling of fill values for[Sensor|Solar][Zenith|Azimuth] resolving MODxl01751 Changed to validate LUT version against a value retrieved from the resolving MODxl02056 Changed to calculate Solar Diffuser angles without adjustment for estimated post launch changes in the MODIS orientation relative to incidentally resolving defects MODxl01766 Also resolves MODxl01947 Changed to ignore fill values in SCI_ABNORM and SCI_STATE rather than treating them as resolving MODxl01780 Changed to use spacecraft ancillary data to recognise when the mirror encoder data is being set by side A or side B and to change calculations accordingly This removes the need for seperate LUTs for Side A and Side B data it makes the new LUTs incompatible with older versions of the and vice versa Also resolves MODxl01685 A more robust GRing algorithm is being which will create a non default GRing anytime there s even a single geolocated pixel in a granule Removed obsolete messages from seed file
Definition: HISTORY.txt:413
names::if
character(len=1000) if
Definition: names.f90:13
fill_table::prev
int32_t * prev
Definition: l1_hmodis_hdf.c:978
L1B_SDS
L1B_SDS
Definition: l1_hmodis_hdf.c:86
mnorm
double mnorm[3]
Definition: l1_modis.c:754
LIFE_IS_GOOD
#define LIFE_IS_GOOD
Definition: passthebuck.h:4
PI
#define PI
Definition: l3_get_org.c:6
CIR_1KM
@ CIR_1KM
Definition: l1_modis.c:90
closel1_modis
int closel1_modis()
Definition: l1_modis.c:1637
modis_sds::nsf
int32_t nsf
Definition: l1_hmodis_hdf.c:151
modis_sds
Definition: l1_hmodis_hdf.c:144
iscan
int32_t iscan
Definition: l1_modis.c:745
compute_alpha
void compute_alpha(float lon[], float lat[], float senz[], float sena[], double mnorm[3], int npix, float alpha[])
Definition: compute_alpha.c:8
unix_to_tai93
double unix_to_tai93(double unixtime)
Definition: yds2tai.c:11
set_l1rec_scanvals
int set_l1rec_scanvals(const int32_t iscan, l1str *l1rec)
Definition: l1_modis.c:1388
sena
float * sena
Definition: l1_modis.c:761
modis_file::shortname
char shortname[FILENAME_MAX]
Definition: l1_hmodis_hdf.c:132
lonlat
subroutine lonlat(alon, alat, xlon, ylat)
Definition: lonlat.f:2
GEO_TAISEC
@ GEO_TAISEC
Definition: l1_modis.c:41
read_sfcorrtables
void read_sfcorrtables(sfcorr_table *sfcorr, int32_t sensorID)
Definition: l1_modis.c:802
tai93_to_unix
double tai93_to_unix(double tai93)
Definition: yds2tai.c:16
L1B_NUM_COEFFS
@ L1B_NUM_COEFFS
Definition: l1_modis.c:117
hico.value_locate.x
list x
Definition: value_locate.py:64
sfcorr_table::nsf
int16_t nsf
Definition: l1_hmodis_hdf.c:797
strdup
char * strdup(const char *)
readl1_lonlat_modis
int readl1_lonlat_modis(filehandle *l1file, int32_t line, l1str *l1rec)
Definition: l1_modis.c:1684
color_dtdb.y
y
Definition: color_dtdb.py:430
modis_sds::nbands
int32_t nbands
Definition: l1_hmodis_hdf.c:147
GEO_NUM_SDS
@ GEO_NUM_SDS
Definition: l1_modis.c:64
fill_table::next
int32_t * next
Definition: l1_hmodis_hdf.c:979
bandinfo::flagindex
int16_t flagindex
Definition: l1_hmodis_hdf.c:570
color_dtdb.scale
int scale
Definition: color_dtdb.py:459
nbands
int32_t nbands
Definition: l1_modis.c:748
modis_file::id
int32_t id
Definition: l1_hmodis_hdf.c:130
l1_input
l1_input_t * l1_input
Definition: l1_options.c:9
open_modis_geofile
int open_modis_geofile(modis_file *geofile)
Definition: l1_modis.c:406
GEO_HGT
@ GEO_HGT
Definition: l1_modis.c:59
want_verbose
int want_verbose
Definition: genutils_globals.c:3
unix2yds
void unix2yds(double usec, short *year, short *day, double *secs)
RAD_SCALE
@ RAD_SCALE
Definition: l1_modis.c:115
filename
char filename[FILENAME_MAX]
Definition: atrem_corl1.h:122
init_detfill
void init_detfill(fill_table *fill)
Definition: l1_modis.c:982
generalauxtype::double
integer, parameter double
Definition: GeneralAuxType.f90:27
fix_mac_rpath.line
line
Definition: fix_mac_rpath.py:56
DFNT_INT16
HDF4 data type of the output SDS Default is DFNT_FLOAT32 Common types used DFNT_INT16
Definition: HOWTO_Add_a_product.txt:67
taisec
double taisec
Definition: l1_modis.c:751
read_sds_1scan
int read_sds_1scan(sds_struct *sds, int32_t iscan, int32_t scandim, int32_t per_scan)
Definition: l1_modis.c:489
get_hdfeos_meta
int get_hdfeos_meta(int32_t sd_id, char *attribute, char *name, char *data)
Definition: hdf4_utils.c:326
color_dtdb.result
result
Definition: color_dtdb.py:197
RADEG
#define RADEG
Definition: czcs_ctl_pt.c:5
GEO_SENA
@ GEO_SENA
Definition: l1_modis.c:63
SI
instead the metadata field ProcessingEnvinronment is filled in from the output of a call to the POSIX compliant function uname from within the L1B code A small bug in L1B_Tables an incorrect comparison of RVS coefficients for TEBs to RVS coefficients for RSBs was being made This was replaced with a comparison between TEB coefficients This error never resulted in an incorrect RVS correction but did lead to recalculating the coefficients for each detector in a thermal band even if the coefficients were the same for all detectors To reduce to overall size of the reflective LUT HDF fill values were eliminated from all LUTs previously dimensioned where and where NUM_TIMES is the number of time dependent table pieces In Preprocess a small error where the trailing dropped scan counter was incremented when the leading dropped scan counter should have been was fixed This counter is internal only and is not yet used for any chiefly to casting of were added to make it LINUX compatible Output of code run on LINUX machines displays differences of at most scaled sector incalculable values of the Emissive calibration factor and incalculable values of SV or BB averages was moved outside the loop over frames in Emissive_Cal c since none of these quantities are frame dependent Initialization of b1 and XMS values in Preprocess c routine Process_OBCENG_Emiss was moved inside the detector loops The code was altered so that if up to five scans are dropped between the leading middle or middle trailing the leading or trailing granule will still be used in emissive calibration to form a cross granule average QA bits and are set for a gap between the leading middle and middle trailing granules respectively This may in rare instances lead to a change in emissive calibration coefficients for scans at the beginning or end of a granule A small bug in the Band correction algorithm was corrected an uncertainty value was being checked against an upper bound whereas the proper quantity to be checked was the corresponding SI
Definition: HISTORY.txt:595
seadasutils.aquarius_utils.start
start
Definition: aquarius_utils.py:27
modis_file::starttime
double starttime
Definition: l1_hmodis_hdf.c:133
data
no change in intended resolving MODur00064 Corrected handling of bad ephemeris attitude data
Definition: HISTORY.txt:356
nvals
int32_t nvals
Definition: l1_modis.c:746
modis_file::nlines
int32_t nlines
Definition: l1_hmodis_hdf.c:135
GEO_NUM_COEFFS
@ GEO_NUM_COEFFS
Definition: l1_modis.c:79
MIN_BAD_SI
#define MIN_BAD_SI
Definition: l1_modis.c:536
sdslist::scandim
int32_t scandim
Definition: l1_hmodis_hdf.c:35
color_dtdb.edges
edges
Definition: color_dtdb.py:331
fill_dead_detectors
void fill_dead_detectors(const int32_t iband, const modis_sds mds, double *data)
Definition: l1_modis.c:1038
basename
#define basename(s)
Definition: l0chunk_modis.c:29
init_geo
int init_geo(const char filename[FILENAME_MAX])
Definition: l1_modis.c:444
RSB_500
@ RSB_500
Definition: l1_modis.c:88
interp_modis_Longitude
double interp_modis_Longitude(double val[4], double dx, double dy)
Definition: l1_modis.c:1150
BAD_FLT
#define BAD_FLT
Definition: jplaeriallib.h:19
modpath_1km
int modpath_1km(const char *oldpath, const char *newchars, char *newpath)
Definition: l1_modis.c:182
init_l1b
int init_l1b(const char filename[FILENAME_MAX], int32_t *max_resolution)
Definition: l1_modis.c:287
GEO_MSIDE
@ GEO_MSIDE
Definition: l1_modis.c:42
sfcorr_table::corr
double * corr
Definition: l1_hmodis_hdf.c:799
sensorId2SensorName
const char * sensorId2SensorName(int sensorId)
Definition: sensorInfo.c:273
radiance2bt
void radiance2bt(l1str *l1rec, int resolution)
Definition: brightness.c:170
L1B_NUM_SDS
@ L1B_NUM_SDS
Definition: l1_modis.c:92
fabs
#define fabs(a)
Definition: misc.h:93
sdslist::name
char * name
Definition: l1_hmodis_hdf.c:36
REFL_OFFSET
@ REFL_OFFSET
Definition: l1_modis.c:114
fopen_warn
void fopen_warn(const char *filename, const char *file, const int32_t line)
Definition: hdf4_utils.c:261
open_modis_l1bfile
int open_modis_l1bfile(modis_file *l1bfile)
Definition: l1_modis.c:230
bandinfo::bandindex
int16_t bandindex
Definition: l1_hmodis_hdf.c:569
free_sds
void free_sds(sds_struct *sds)
Definition: hdf4_utils.c:161
interp_modis_EVdata
double interp_modis_EVdata(double val[4], double dx, double dy)
Definition: l1_modis.c:1116
modis_file::npixls
int32_t npixls
Definition: l1_hmodis_hdf.c:136
readall_sds
int readall_sds(sds_struct *sds)
Definition: hdf4_utils.c:138
s
data_t s[NROOTS]
Definition: decode_rs.h:75
esdist
real *8 function esdist(iyr, iday, msec)
Definition: esdist.f:3
GEO_MNORM
@ GEO_MNORM
Definition: l1_modis.c:44
read_att
#define read_att(obj_id, attname, valptr)
Definition: hdf4utils.h:177
rdsensorinfo
int32_t rdsensorinfo(int32_t, int32_t, const char *, void **)
Definition: rdsensorinfo.c:69
subframe_correction
void subframe_correction(const int32_t sensorID, const int32_t iband, const modis_sds mds, double *data)
Definition: l1_modis.c:915
GEO_SOLZ
@ GEO_SOLZ
Definition: l1_modis.c:60
GEO_LAT
@ GEO_LAT
Definition: l1_modis.c:58
subsensorId2SubsensorDir
const char * subsensorId2SubsensorDir(int subsensorId)
Definition: sensorInfo.c:329
offset
l2prod offset
Definition: HOWTO_Add_a_product.txt:70
DFNT_FLOAT32
HDF4 data type of the output SDS Default is DFNT_FLOAT32 Common types used DFNT_FLOAT32
Definition: HOWTO_Add_a_product.txt:67
GEO_SENZ
@ GEO_SENZ
Definition: l1_modis.c:62
fill_table::ndets
int32_t ndets
Definition: l1_hmodis_hdf.c:977
modis_file
Definition: l1_hmodis_hdf.c:129
modis_file::sds
sds_struct * sds
Definition: l1_hmodis_hdf.c:138
TEB_1KM
@ TEB_1KM
Definition: l1_modis.c:91
mside
int32_t mside
Definition: l1_modis.c:752
REF_SDS
REF_SDS
Definition: l1_hmodis_hdf.c:40
fix_mac_rpath.id
string id
Definition: fix_mac_rpath.py:28
solz
float * solz
Definition: l1_modis.c:758
i
int i
Definition: decode_rs.h:71
strcpy
How many dimensions is the output array Default is Not sure if anything above will work correctly strcpy(l2prod->title, "no title yet")
load_att_byname
int load_att_byname(int32_t obj_id, const char *attname, att_struct *att)
Definition: hdf4_utils.c:31
sensorID
int32_t sensorID[MAXNFILES]
Definition: l2bin.cpp:91
interp_modis_var
double(* interp_modis_var)(double val[4], double dx, double dy)
Definition: l1_modis.c:1171
MODISA
#define MODISA
Definition: sensorDefs.h:19
EOSMETALEN
#define EOSMETALEN
Definition: filetype.c:23
modis_sds::nframes
int32_t nframes
Definition: l1_hmodis_hdf.c:150
modis_interp
int modis_interp(const modis_sds s, const int32_t sdsband, const int32_t newres, double *olddata, double *newdata)
Definition: l1_modis.c:1173
readl1_modis
int readl1_modis(filehandle *l1file, const int32_t line, l1str *l1rec, int lonlat)
Definition: l1_modis.c:1448
GEO_SOLA
@ GEO_SOLA
Definition: l1_modis.c:61
npix
int32_t npix
Definition: l1_modis.c:747
GEO_COEFFS
GEO_COEFFS
Definition: l1_hmodis_hdf.c:77
isodate2unix
double isodate2unix(const char *isodate)
Definition: unix2isodate.c:61
RAD_OFFSET
@ RAD_OFFSET
Definition: l1_modis.c:116


ABOUT

Our Team

Contact Us

DATA

Get Started

Find Data

Use Data

TOOLS

SeaBASS

SeaDAS

Giovanni Data Visualizer

Worldview Mapping Interface

Mission Quality Monitor

HELP

How to Use Earthdata Search

Level 3 & 4 Browser Help

File Search Help

Search & Download Help

Earthdata Forum

Join Mailing List

How to Cite Data

CONNECT WITH US
COMPLIANCE

Web Privacy Policy

Data & Information Policy

Communications Policy

Accessibility

Freedom of Information Act

USA.gov

Responsible NASA Official: Sean Bailey
Curator: OceanColor Webmaster