OB.DAAC Logo
NASA Logo
Ocean Color Science Software

ocssw V2022
l1asubset_modis.c
Go to the documentation of this file.
1 /*
2  Modification history:
3  Programmer Organization Date Description of change
4  -------------- ------------ -------- ---------------------
5  Joel Gales Futuretech 03/28/03 Original Development
6  Joel Gales Futuretech 10/18/07 Rescaled pixels with
7  value less than 4095
8  are set to 4095.
9  Joel Gales Futuretech 10/26/07 Back out prev change
10  Joel Gales Futuretech 10/26/07 Put back prev change
11  put check if < 0 before
12  setting to 4095
13  */
14 
15 
16 #include <stdio.h>
17 #include <stdlib.h>
18 #include <math.h>
19 #include <time.h>
20 #include <string.h>
21 #include <unistd.h>
22 #include "hdf.h"
23 #include "mfhdf.h"
24 
25 #include <timeutils.h>
26 
27 #define basename(s) (strrchr((s), '/') == NULL ? (s) : strrchr((s), '/') + 1)
28 
29 
30 #define VERSION "1.45"
31 
32 void usage(char *progname) {
33  printf("This is version %s of %s (compiled on %s %s)\n",
34  VERSION, progname, __DATE__, __TIME__);
35 
36  printf("\nUsage: %s ifile ofile\n", progname);
37  printf(" ifile = input filename\n");
38  printf(" ofile = output filename\n");
39  exit(0);
40 }
41 
42 int main(int argc, char *argv[]) {
43 
44  int32 i;
45  int32 j;
46  int32 k;
47  int32 l;
48  int32 n;
49  int32 sds_id;
50  int32 i32;
51  int32 dims[8];
52  int32 rank;
53  int32 nelem;
54  int32 dtype;
55  int32 nattrs;
56  int32 status;
57 
58  int32 vdid;
59  int32 vdid_w;
60  int32 nflds;
61  int32 nrec;
62  int32 nread;
63  int32 nwrite;
64  int32 rescale;
65  int32 last_bad = 0;
66  int32 daynotneeded[] = {11, 12, 13, -1}; /* 0-based */
67  int32 ngtnotneeded[] = {1, 4, 10, 13, 14, 15, 16, -1}; /* 0-based */
68 
69  int32 nlines;
70  int32 npix;
71 
72  int16 idata_HR[16];
73  int16 sz;
74  int16 *data_1km;
75 
76  float32 *median_250_500;
77  float32 sumx, sumy, sumx2, sumxy, a, b, scale_parm[2];
78  float32 f32;
79  float32 rescale_parm[8];
80 
81  static char buffer[2 * 2048];
82  static char buf2[2048];
83 
84  char sat[8];
85 
86  char *data;
87  char *cptr;
88 
89  int32 HDFfid_r;
90  int32 sd_id_r;
91  int32 HDFfid_w;
92  int32 sd_id_w;
93  int32 sds_id_w;
94  int32 ndatasets;
95  int32 nglobal_attr;
96  int32 count;
97  int32 dim_id_r;
98  int32 dim_id_w;
99  int32 start[3] = {0, 0, 0};
100  int32 start_w[3] = {0, 0, 0};
101  int32 maxsize;
102  int32 *lonebuf;
103  int32 coremeta = 0;
104 
105  char *tmp_str;
106 
107  int compare_int16(const void *, const void *);
108  int compare_float32(const void *, const void *);
109 
110  FILE *fp;
111 
112  struct tm *tmnow;
113  time_t tnow;
114 
115  char prodtime[80];
116 
117 
118  if (argc == 1) {
119  usage(argv[0]);
120  exit(0);
121  }
122 
123  time(&tnow);
124  tmnow = gmtime(&tnow);
125  strftime(prodtime, 80, "%Y-%m-%dT%XZ", tmnow);
126 
127  HDFfid_r = Hopen(argv[1], DFACC_READ, 0);
128  status = Vstart(HDFfid_r);
129  sd_id_r = SDstart(argv[1], DFACC_RDONLY);
130 
131 
132  /* Create output HDF file */
133  HDFfid_w = Hopen(argv[2], DFACC_CREATE, 0);
134  status = Vstart(HDFfid_w);
135  sd_id_w = SDstart(argv[2], DFACC_RDWR);
136 
137  printf("argc: %d\n", argc);
138 
139  /*
140  rescale = 0 (Perform rescaling)
141  rescale = 1 (Perform rescaling with granule-based parameters)
142  rescale = -1 (No rescaling)
143  */
144  rescale = 0;
145  if (argc >= 4) {
146  rescale = atoi(argv[3]);
147 
148  }
149 
150  /* Determine number of datasets in input files */
151  SDfileinfo(sd_id_r, &ndatasets, &nglobal_attr);
152 
153  /* For each dataset ... */
154  for (j = 0; j < ndatasets; j++) {
155  sds_id = SDselect(sd_id_r, j);
156  status = SDgetinfo(sds_id, buffer, &rank, dims, &dtype, &nattrs);
157 
158  if (strcmp(buffer, "SD_250m") == 0) continue;
159  if (strcmp(buffer, "SD_500m") == 0) continue;
160 
161  if (strcmp(buffer, "SRCA_250m") == 0) continue;
162  if (strcmp(buffer, "SRCA_500m") == 0) continue;
163 
164  if (strcmp(buffer, "BB_250m") == 0) continue;
165  if (strcmp(buffer, "BB_500m") == 0) continue;
166 
167  if (strcmp(buffer, "SV_250m") == 0) continue;
168  if (strcmp(buffer, "SV_500m") == 0) continue;
169 
170  if (strcmp(buffer, "EV_250m") == 0) continue;
171  if (strcmp(buffer, "EV_500m") == 0) continue;
172 
173 
174  /*
175  printf("Name :%s rank: %d type: %d\n", buffer, rank, dtype);
176  printf("dims: %d %d %d\n", dims[0], dims[1], dims[2]);
177  */
178 
179  sds_id_w = SDcreate(sd_id_w, buffer, dtype, rank, dims);
180  if (sds_id_w == -1) {
181  printf("Field: %s cannot be created\n", buffer);
182  exit(-1);
183  }
184 
185  for (i = 0; i < nattrs; i++) {
186  status = SDreadattr(sds_id, i, (VOIDP) buf2);
187  status = SDattrinfo(sds_id, i, buffer, &dtype, &count);
188  status = SDsetattr(sds_id_w, buffer, dtype, count, (VOIDP) buf2);
189  }
190 
191  for (i = 0; i < rank; i++) {
192  dim_id_r = SDgetdimid(sds_id, i);
193  dim_id_w = SDgetdimid(sds_id_w, i);
194  SDdiminfo(dim_id_r, buffer, &count, &dtype, &nattrs);
195  SDsetdimname(dim_id_w, buffer);
196  }
197 
198  SDendaccess(sds_id_w);
199  SDendaccess(sds_id);
200  } /* End Create SD */
201 
202 
203  /* Copy Datasets */
204  /* ------------- */
205  n = 0;
206  for (j = 0; j < ndatasets; j++) {
207  sds_id = SDselect(sd_id_r, j);
208  status = SDgetinfo(sds_id, buffer, &rank, dims, &dtype, &nattrs);
209 
210  if (strcmp(buffer, "SD_250m") == 0) continue;
211  if (strcmp(buffer, "SD_500m") == 0) continue;
212 
213  if (strcmp(buffer, "SRCA_250m") == 0) continue;
214  if (strcmp(buffer, "SRCA_500m") == 0) continue;
215 
216  if (strcmp(buffer, "BB_250m") == 0) continue;
217  if (strcmp(buffer, "BB_500m") == 0) continue;
218 
219  if (strcmp(buffer, "SV_250m") == 0) continue;
220  if (strcmp(buffer, "SV_500m") == 0) continue;
221 
222  if (strcmp(buffer, "EV_250m") == 0) continue;
223  if (strcmp(buffer, "EV_500m") == 0) continue;
224 
225  sds_id_w = SDselect(sd_id_w, n++);
226 
227  /*
228  printf("% 3d Name :%s rank: %d type: %d\n", j, buffer, rank, dtype);
229  for (k=0; k<rank; k++) printf("dim %d: %d\n", k, dims[k]);
230  */
231 
232  nelem = dims[0];
233  for (k = 1; k < rank; k++) nelem *= dims[k];
234 
235  data = (char *) calloc(nelem, DFKNTsize(dtype));
236 
237  status = SDreaddata(sds_id, start, NULL, dims, (VOIDP) data);
238 
239 
240 #if 0
241  /* Check whether last scan is bad */
242  /* ------------------------------ */
243  if (strcmp(buffer, "Frame count array") == 0 && last_bad == 0) {
244  memcpy(&i16, &data[(6 * 202 + 1)*2], 2);
245  if (i16 == 0) last_bad = 1;
246  printf("last_bad: %d\n", last_bad);
247  }
248 #endif
249 
250 
251  /* Flag unneeded day bands */
252  /* ----------------------- */
253  if (strcmp(buffer, "EV_1km_day") == 0) {
254  printf("Masking Unneeded Day Bands\n");
255  for (k = 0; k < dims[0]; k++) {
256  l = 0;
257  while (daynotneeded[l] != -1) {
258  if (k == 0) printf("%d\n", daynotneeded[l]);
259  memset(&data[k * 2 * dims[1] * dims[2] + daynotneeded[l]*2 * dims[2]],
260  255, 2 * dims[2]);
261  l++;
262  }
263  }
264  }
265 
266  /* Flag unneeded night bands */
267  /* ------------------------- */
268  if (strcmp(buffer, "EV_1km_night") == 0) {
269  printf("Masking Unneeded Night Bands\n");
270  for (k = 0; k < dims[0]; k++) {
271  l = 0;
272  while (ngtnotneeded[l] != -1) {
273  if (k == 0) printf("%d\n", ngtnotneeded[l]);
274  memset(&data[k * 2 * dims[1] * dims[2] + ngtnotneeded[l]*2 * dims[2]],
275  255, 2 * dims[2]);
276  l++;
277  }
278  }
279  }
280 
281 
282  status = SDwritedata(sds_id_w, start_w, NULL, dims, (VOIDP) data);
283  if (status == -1) {
284  printf("write status: %d\n\n", status);
285  exit(-1);
286  }
287 
288  free(data);
289 
290  SDendaccess(sds_id);
291  SDendaccess(sds_id_w);
292  }
293 
294 
295  if (SDselect(sd_id_r, SDnametoindex(sd_id_r, "EV_250m")) == -1 ||
296  SDselect(sd_id_r, SDnametoindex(sd_id_r, "EV_500m")) == -1) {
297  printf("No 250/500 EV field found.\n");
298  rescale = -1;
299  }
300 
301  if (rescale == -1) {
302  printf("No rescaling performed\n");
303  goto skip_rescale;
304  }
305 
306 
307  /* Read rescaling parameters */
308  /* ------------------------- */
309  if (rescale == 0) {
310  if ((tmp_str = getenv("OCDATAROOT")) == NULL) {
311  printf("OCDATAROOT environment variable is not defined.\n");
312  return (1);
313  }
314 
315  strcpy(buf2, tmp_str);
316  strcat(buf2, "/modisa/l1a_scaling.dat");
317 
318  if ((fp = fopen(buf2, "r")) == NULL) {
319  printf("Error: Unable to open rescaling parameter file: %s.\n", buf2);
320  return (1);
321  }
322 
323  for (i = 0; i < 4; i++) {
324  if (fscanf(fp, "%f %f\n",
325  &rescale_parm[2 * i], &rescale_parm[2 * i + 1]) == EOF) {
326  printf("Error: Insufficent number of scaling parameters: %s.\n",
327  buf2);
328  return (1);
329  }
330  }
331  }
332 
333 
334  /* Rescale "4095" 1km B pixels (B10) using 500m B pixels (B3) */
335  /* ---------------------------------------------------------- */
336  printf("Rescale saturated 1km B pixels (B10) using 500m B pixels (B3)\n");
337  sz = 2;
338  sds_id_w = SDselect(sd_id_w, SDnametoindex(sd_id_w, "EV_1km_day"));
339  status = SDgetinfo(sds_id_w, buffer, &rank, dims, &dtype, &nattrs);
340  start[1] = 2;
341  dims[1] = 1;
342 
343  nelem = dims[0];
344  for (k = 1; k < rank; k++) nelem *= dims[k];
345 
346  nlines = dims[0];
347  npix = dims[2];
348 
349  data_1km = (int16 *) calloc(nelem, sz);
350  median_250_500 = (float32 *) calloc(nelem, sizeof (float32));
351 
352  status = SDreaddata(sds_id_w, start, NULL, dims, (VOIDP) data_1km);
353 
354  sds_id = SDselect(sd_id_r, SDnametoindex(sd_id_r, "EV_500m"));
355  status = SDgetinfo(sds_id, buffer, &rank, dims, &dtype, &nattrs);
356 
357  dims[0] = 2;
358  start[1] = 0;
359  dims[1] = 1;
360 
361  nelem = dims[0];
362  for (k = 1; k < rank; k++) nelem *= dims[k];
363 
364  data = (char *) calloc(nelem, sz);
365 
366  sumx = 0;
367  sumy = 0;
368  sumx2 = 0;
369  sumxy = 0;
370  n = 0;
371  for (i = 0; i < nlines; i++) {
372  start[0] = dims[0] * i;
373  status = SDreaddata(sds_id, start, NULL, dims, (VOIDP) data);
374 
375  if (i % 400 == 0) printf("i: %d\n", i);
376 
377  for (j = 0; j < npix; j++) {
378 
379  memcpy(&idata_HR[0], &data[sz * 2 * (j + 0 * npix)], 2 * sz);
380  memcpy(&idata_HR[2], &data[sz * 2 * (j + 1 * npix)], 2 * sz);
381 
382  qsort((void *) idata_HR, dims[0] * dims[0], sizeof (int16),
383  compare_int16);
384 
385  for (k = 0; k < 4; k++) {
386  if (idata_HR[k] > 0) {
387 
388  if ((k % 2) == 0)
389  median_250_500[i * npix + j] =
390  0.5 * (idata_HR[(2 + k) / 2] + idata_HR[(4 + k) / 2]);
391  else
392  median_250_500[i * npix + j] = (float32) idata_HR[(3 + k) / 2];
393 
394  if (data_1km[i * npix + j] > 0 && data_1km[i * npix + j] < 4095) {
395  n++;
396  sumx += median_250_500[i * npix + j];
397  sumy += data_1km[i * npix + j];
398  sumx2 += median_250_500[i * npix + j] * median_250_500[i * npix + j];
399  sumxy += median_250_500[i * npix + j] * data_1km[i * npix + j];
400  }
401  break;
402  }
403  }
404  }
405  }
406  a = (sumxy - sumx * sumy / n) / (sumx2 - sumx * sumx / n);
407  b = (sumy / n) - a * (sumx / n);
408  scale_parm[0] = b;
409  scale_parm[1] = a;
410  status = SDsetattr(sd_id_w, "B_10_scale_parm_computed", DFNT_FLOAT32, 2,
411  (VOIDP) scale_parm);
412 
413  if (rescale == 0) {
414  b = rescale_parm[0];
415  a = rescale_parm[1];
416  status = SDsetattr(sd_id_w, "B_10_scale_parm", DFNT_FLOAT32, 2,
417  (VOIDP) & rescale_parm[0]);
418  } else {
419  status = SDsetattr(sd_id_w, "B_10_scale_parm", DFNT_FLOAT32, 2,
420  (VOIDP) scale_parm);
421  }
422  printf("a: %f b: %f\n", a, b);
423 
424  for (i = 0; i < nlines; i++)
425  for (j = 0; j < npix; j++) {
426  if (data_1km[i * npix + j] == 4095) {
427  data_1km[i * npix + j] = a * median_250_500[i * npix + j] + b;
428  if (data_1km[i * npix + j] >= 0 && data_1km[i * npix + j] < 4095)
429  data_1km[i * npix + j] = 4095;
430  }
431  if (data_1km[i * npix + j] < 0) data_1km[i * npix + j] = 32767;
432  }
433 
434  status = SDgetinfo(sds_id_w, buffer, &rank, dims, &dtype, &nattrs);
435  start[0] = 0;
436  start[1] = 2;
437  dims[1] = 1;
438  status = SDwritedata(sds_id_w, start, NULL, dims, (VOIDP) data_1km);
439 
440  free(data);
441  free(data_1km);
442  free(median_250_500);
443 
444  SDendaccess(sds_id);
445  SDendaccess(sds_id_w);
446 
447 
448 
449  /* Rescale "4095" 1km G pixels (B12) using 500m G pixels (B4) */
450  /* ---------------------------------------------------------- */
451  printf("Rescale saturated 1km G pixels (B12) using 500m G pixels (B4)\n");
452  sz = 2;
453  sds_id_w = SDselect(sd_id_w, SDnametoindex(sd_id_w, "EV_1km_day"));
454  status = SDgetinfo(sds_id_w, buffer, &rank, dims, &dtype, &nattrs);
455  start[1] = 4;
456  dims[1] = 1;
457 
458  nelem = dims[0];
459  for (k = 1; k < rank; k++) nelem *= dims[k];
460 
461  nlines = dims[0];
462  npix = dims[2];
463 
464  data_1km = (int16 *) calloc(nelem, sz);
465  median_250_500 = (float32 *) calloc(nelem, sizeof (float32));
466 
467  status = SDreaddata(sds_id_w, start, NULL, dims, (VOIDP) data_1km);
468 
469  sds_id = SDselect(sd_id_r, SDnametoindex(sd_id_r, "EV_500m"));
470  status = SDgetinfo(sds_id, buffer, &rank, dims, &dtype, &nattrs);
471 
472  dims[0] = 2;
473  start[1] = 1;
474  dims[1] = 1;
475 
476  nelem = dims[0];
477  for (k = 1; k < rank; k++) nelem *= dims[k];
478 
479  data = (char *) calloc(nelem, sz);
480 
481  sumx = 0;
482  sumy = 0;
483  sumx2 = 0;
484  sumxy = 0;
485  n = 0;
486  for (i = 0; i < nlines; i++) {
487  start[0] = dims[0] * i;
488  status = SDreaddata(sds_id, start, NULL, dims, (VOIDP) data);
489 
490  if (i % 400 == 0) printf("i: %d\n", i);
491 
492  for (j = 0; j < npix; j++) {
493 
494  memcpy(&idata_HR[0], &data[sz * 2 * (j + 0 * npix)], 2 * sz);
495  memcpy(&idata_HR[2], &data[sz * 2 * (j + 1 * npix)], 2 * sz);
496 
497  qsort((void *) idata_HR, dims[0] * dims[0], sizeof (int16),
498  compare_int16);
499 
500  for (k = 0; k < 4; k++) {
501  if (idata_HR[k] > 0) {
502 
503  if ((k % 2) == 0)
504  median_250_500[i * npix + j] =
505  0.5 * (idata_HR[(2 + k) / 2] + idata_HR[(4 + k) / 2]);
506  else
507  median_250_500[i * npix + j] = (float32) idata_HR[(3 + k) / 2];
508 
509  if (data_1km[i * npix + j] > 0 && data_1km[i * npix + j] < 4095) {
510  n++;
511  sumx += median_250_500[i * npix + j];
512  sumy += data_1km[i * npix + j];
513  sumx2 += median_250_500[i * npix + j] * median_250_500[i * npix + j];
514  sumxy += median_250_500[i * npix + j] * data_1km[i * npix + j];
515  }
516  break;
517  }
518  }
519  }
520  }
521  a = (sumxy - sumx * sumy / n) / (sumx2 - sumx * sumx / n);
522  b = (sumy / n) - a * (sumx / n);
523  scale_parm[0] = b;
524  scale_parm[1] = a;
525  status = SDsetattr(sd_id_w, "B_12_scale_parm_computed", DFNT_FLOAT32, 2,
526  (VOIDP) scale_parm);
527 
528  if (rescale == 0) {
529  b = rescale_parm[2];
530  a = rescale_parm[3];
531  status = SDsetattr(sd_id_w, "B_12_scale_parm", DFNT_FLOAT32, 2,
532  (VOIDP) & rescale_parm[2]);
533  } else {
534  status = SDsetattr(sd_id_w, "B_12_scale_parm", DFNT_FLOAT32, 2,
535  (VOIDP) scale_parm);
536  }
537  printf("a: %f b: %f\n", a, b);
538 
539  for (i = 0; i < nlines; i++)
540  for (j = 0; j < npix; j++) {
541  if (data_1km[i * npix + j] == 4095) {
542  data_1km[i * npix + j] = a * median_250_500[i * npix + j] + b;
543  if (data_1km[i * npix + j] >= 0 && data_1km[i * npix + j] < 4095)
544  data_1km[i * npix + j] = 4095;
545  }
546  if (data_1km[i * npix + j] < 0) data_1km[i * npix + j] = 32767;
547  }
548 
549  status = SDgetinfo(sds_id_w, buffer, &rank, dims, &dtype, &nattrs);
550  start[0] = 0;
551  start[1] = 4;
552  dims[1] = 1;
553  status = SDwritedata(sds_id_w, start, NULL, dims, (VOIDP) data_1km);
554 
555  free(data);
556  free(data_1km);
557  free(median_250_500);
558 
559  SDendaccess(sds_id);
560  SDendaccess(sds_id_w);
561 
562 
563  /* Rescale "4095" 1km R pixels (B13) using 250m R pixels (B1) */
564  /* ---------------------------------------------------------- */
565  printf("Rescale saturated 1km R pixels (B13) using 250m R pixels (B1)\n");
566  sz = 2;
567  sds_id_w = SDselect(sd_id_w, SDnametoindex(sd_id_w, "EV_1km_day"));
568  status = SDgetinfo(sds_id_w, buffer, &rank, dims, &dtype, &nattrs);
569  start[1] = 5;
570  dims[1] = 1;
571 
572  nelem = dims[0];
573  for (k = 1; k < rank; k++) nelem *= dims[k];
574 
575  nlines = dims[0];
576  npix = dims[2];
577 
578  data_1km = (int16 *) calloc(nelem, sz);
579  median_250_500 = (float32 *) calloc(nelem, sizeof (float32));
580 
581  status = SDreaddata(sds_id_w, start, NULL, dims, (VOIDP) data_1km);
582 
583  sds_id = SDselect(sd_id_r, SDnametoindex(sd_id_r, "EV_250m"));
584  status = SDgetinfo(sds_id, buffer, &rank, dims, &dtype, &nattrs);
585 
586  dims[0] = 4;
587  start[1] = 0;
588  dims[1] = 1;
589 
590  nelem = dims[0];
591  for (k = 1; k < rank; k++) nelem *= dims[k];
592 
593  data = (char *) calloc(nelem, sz);
594 
595  sumx = 0;
596  sumy = 0;
597  sumx2 = 0;
598  sumxy = 0;
599  n = 0;
600  for (i = 0; i < nlines; i++) {
601  start[0] = dims[0] * i;
602  status = SDreaddata(sds_id, start, NULL, dims, (VOIDP) data);
603 
604  if (i % 400 == 0) printf("i: %d\n", i);
605 
606  for (j = 0; j < npix; j++) {
607 
608  memcpy(&idata_HR[0], &data[sz * 4 * (j + 0 * npix)], 4 * sz);
609  memcpy(&idata_HR[4], &data[sz * 4 * (j + 1 * npix)], 4 * sz);
610  memcpy(&idata_HR[8], &data[sz * 4 * (j + 2 * npix)], 4 * sz);
611  memcpy(&idata_HR[12], &data[sz * 4 * (j + 3 * npix)], 4 * sz);
612 
613  qsort((void *) idata_HR, dims[0] * dims[0], sizeof (int16),
614  compare_int16);
615 
616  for (k = 0; k < dims[0] * dims[0]; k++) {
617  if (idata_HR[k] > 0) {
618 
619  if ((k % 2) == 0)
620  median_250_500[i * npix + j] =
621  0.5 * (idata_HR[(14 + k) / 2] + idata_HR[(16 + k) / 2]);
622  else
623  median_250_500[i * npix + j] = (float32) idata_HR[(15 + k) / 2];
624 
625  if (data_1km[i * npix + j] > 0 && data_1km[i * npix + j] < 4095) {
626  n++;
627  sumx += median_250_500[i * npix + j];
628  sumy += data_1km[i * npix + j];
629  sumx2 += median_250_500[i * npix + j] * median_250_500[i * npix + j];
630  sumxy += median_250_500[i * npix + j] * data_1km[i * npix + j];
631  }
632  break;
633  }
634  }
635  }
636  }
637  a = (sumxy - sumx * sumy / n) / (sumx2 - sumx * sumx / n);
638  b = (sumy / n) - a * (sumx / n);
639  scale_parm[0] = b;
640  scale_parm[1] = a;
641  status = SDsetattr(sd_id_w, "B_13_scale_parm_computed", DFNT_FLOAT32, 2,
642  (VOIDP) scale_parm);
643 
644  if (rescale == 0) {
645  b = rescale_parm[4];
646  a = rescale_parm[5];
647  status = SDsetattr(sd_id_w, "B_13_scale_parm", DFNT_FLOAT32, 2,
648  (VOIDP) & rescale_parm[4]);
649  } else {
650  status = SDsetattr(sd_id_w, "B_13_scale_parm", DFNT_FLOAT32, 2,
651  (VOIDP) scale_parm);
652  }
653  printf("a: %f b: %f\n", a, b);
654 
655  for (i = 0; i < nlines; i++)
656  for (j = 0; j < npix; j++) {
657  if (data_1km[i * npix + j] == 4095) {
658  data_1km[i * npix + j] = a * median_250_500[i * npix + j] + b;
659  if (data_1km[i * npix + j] >= 0 && data_1km[i * npix + j] < 4095)
660  data_1km[i * npix + j] = 4095;
661  }
662  /* if (data_1km[i*npix+j] == -32768) data_1km[i*npix+j] = 32767;*/
663  if (data_1km[i * npix + j] < 0) data_1km[i * npix + j] = 32767;
664  }
665 
666  status = SDgetinfo(sds_id_w, buffer, &rank, dims, &dtype, &nattrs);
667  start[0] = 0;
668  start[1] = 5;
669  dims[1] = 1;
670  status = SDwritedata(sds_id_w, start, NULL, dims, (VOIDP) data_1km);
671 
672  free(data);
673  free(data_1km);
674  free(median_250_500);
675 
676  SDendaccess(sds_id);
677  SDendaccess(sds_id_w);
678 
679 
680 
681  /* Rescale "4095" 1km IR pixels (B16) using 250m IR pixels (B2) */
682  /* ------------------------------------------------------------ */
683  printf("Rescale saturated 1km IR pixels (B16) using 250m IR pixels (B2)\n");
684  sz = 2;
685  sds_id_w = SDselect(sd_id_w, SDnametoindex(sd_id_w, "EV_1km_day"));
686  status = SDgetinfo(sds_id_w, buffer, &rank, dims, &dtype, &nattrs);
687  start[1] = 10;
688  dims[1] = 1;
689 
690  nelem = dims[0];
691  for (k = 1; k < rank; k++) nelem *= dims[k];
692 
693  nlines = dims[0];
694  npix = dims[2];
695 
696  data_1km = (int16 *) calloc(nelem, sz);
697  median_250_500 = (float32 *) calloc(nelem, sizeof (float32));
698 
699  status = SDreaddata(sds_id_w, start, NULL, dims, (VOIDP) data_1km);
700 
701  sds_id = SDselect(sd_id_r, SDnametoindex(sd_id_r, "EV_250m"));
702  status = SDgetinfo(sds_id, buffer, &rank, dims, &dtype, &nattrs);
703 
704  dims[0] = 4;
705  start[1] = 1;
706  dims[1] = 1;
707 
708  nelem = dims[0];
709  for (k = 1; k < rank; k++) nelem *= dims[k];
710 
711  data = (char *) calloc(nelem, sz);
712 
713  sumx = 0;
714  sumy = 0;
715  sumx2 = 0;
716  sumxy = 0;
717  n = 0;
718  for (i = 0; i < nlines; i++) {
719  start[0] = dims[0] * i;
720  status = SDreaddata(sds_id, start, NULL, dims, (VOIDP) data);
721 
722  if (i % 400 == 0) printf("i: %d\n", i);
723 
724  for (j = 0; j < npix; j++) {
725 
726  memcpy(&idata_HR[0], &data[sz * 4 * (j + 0 * npix)], 4 * sz);
727  memcpy(&idata_HR[4], &data[sz * 4 * (j + 1 * npix)], 4 * sz);
728  memcpy(&idata_HR[8], &data[sz * 4 * (j + 2 * npix)], 4 * sz);
729  memcpy(&idata_HR[12], &data[sz * 4 * (j + 3 * npix)], 4 * sz);
730 
731  qsort((void *) idata_HR, dims[0] * dims[0], sizeof (int16),
732  compare_int16);
733 
734  for (k = 0; k < dims[0] * dims[0]; k++) {
735  if (idata_HR[k] > 0) {
736 
737  if ((k % 2) == 0)
738  median_250_500[i * npix + j] =
739  0.5 * (idata_HR[(14 + k) / 2] + idata_HR[(16 + k) / 2]);
740  else
741  median_250_500[i * npix + j] = (float32) idata_HR[(15 + k) / 2];
742 
743  if (data_1km[i * npix + j] > 0 && data_1km[i * npix + j] < 4095) {
744  n++;
745  sumx += median_250_500[i * npix + j];
746  sumy += data_1km[i * npix + j];
747  sumx2 += median_250_500[i * npix + j] * median_250_500[i * npix + j];
748  sumxy += median_250_500[i * npix + j] * data_1km[i * npix + j];
749  }
750  break;
751  }
752  }
753  }
754  }
755  a = (sumxy - sumx * sumy / n) / (sumx2 - sumx * sumx / n);
756  b = (sumy / n) - a * (sumx / n);
757  scale_parm[0] = b;
758  scale_parm[1] = a;
759  status = SDsetattr(sd_id_w, "B_16_scale_parm_computed", DFNT_FLOAT32, 2,
760  (VOIDP) scale_parm);
761 
762  if (rescale == 0) {
763  b = rescale_parm[6];
764  a = rescale_parm[7];
765  status = SDsetattr(sd_id_w, "B_16_scale_parm", DFNT_FLOAT32, 2,
766  (VOIDP) & rescale_parm[6]);
767  } else {
768  status = SDsetattr(sd_id_w, "B_16_scale_parm", DFNT_FLOAT32, 2,
769  (VOIDP) scale_parm);
770  }
771  printf("a: %f b: %f\n", a, b);
772 
773  for (i = 0; i < nlines; i++)
774  for (j = 0; j < npix; j++) {
775  if (data_1km[i * npix + j] == 4095) {
776  data_1km[i * npix + j] = a * median_250_500[i * npix + j] + b;
777  if (data_1km[i * npix + j] >= 0 && data_1km[i * npix + j] < 4095)
778  data_1km[i * npix + j] = 4095;
779  }
780  if (data_1km[i * npix + j] < 0) data_1km[i * npix + j] = 32767;
781  }
782 
783  status = SDgetinfo(sds_id_w, buffer, &rank, dims, &dtype, &nattrs);
784  start[0] = 0;
785  start[1] = 10;
786  dims[1] = 1;
787  status = SDwritedata(sds_id_w, start, NULL, dims, (VOIDP) data_1km);
788 
789  free(data);
790  free(data_1km);
791  free(median_250_500);
792 
793  SDendaccess(sds_id);
794  SDendaccess(sds_id_w);
795 
796 
797 skip_rescale:
798 
799  /* Write lone Vdatas */
800  /* ----------------- */
801  maxsize = VSlone(HDFfid_r, NULL, 0);
802  lonebuf = (int32 *) calloc(maxsize, sizeof (int32));
803  maxsize = VSlone(HDFfid_r, lonebuf, maxsize);
804 
805 
806  /* For each lone vdata ... */
807  /* ----------------------- */
808  for (i = 0; i < maxsize; i++) {
809  vdid = VSattach(HDFfid_r, lonebuf[i], "r");
810  vdid_w = VSattach(HDFfid_w, -1, "w");
811  VSgetname(vdid, buffer);
812  VSsetname(vdid_w, buffer);
813  /* printf("Vlone name: %s\n", buffer);*/
814 
815  VSsetinterlace(vdid_w, VSgetinterlace(vdid));
816 
817 
818  /* Replace "," in fieldname string by 0 */
819  /* ------------------------------------ */
820  nflds = VSgetfields(vdid, buf2);
821  k = strlen(buf2);
822  for (j = 0; j < k; j++) if (buf2[j] == ',') buf2[j] = 0;
823 
824 
825  /* Parse fieldname string and define fields */
826  /* ---------------------------------------- */
827  cptr = &buf2[0];
828  for (j = 0; j < nflds; j++) {
829  VSfdefine(vdid_w, cptr, VFfieldtype(vdid, j), VFfieldorder(vdid, j));
830 
831  while (1) {
832  cptr++;
833  if (*cptr == 0) {
834  cptr++;
835  break;
836  }
837  }
838  }
839 
840  for (j = 0; j < k - 1; j++) if (buf2[j] == 0) buf2[j] = ',';
841  status = VSsetfields(vdid, buf2);
842  status = VSsetfields(vdid_w, buf2);
843 
844  VSQuerycount(vdid, &nrec);
845  /* printf("nrec: %d vsize: %d\n", nrec, VSsizeof(vdid, buf2));*/
846  data = (char *) calloc(nrec * VSsizeof(vdid, buf2), 1);
847 
848  nread = VSread(vdid, (uint8*) data, nrec, VSgetinterlace(vdid));
849  nwrite = VSwrite(vdid_w, (uint8*) data, nrec, VSgetinterlace(vdid_w));
850  /* printf("nread: %d\n", nread);*/
851 
852  if (nread != nwrite) {
853  printf("%d %d %d\n", i, nread, nwrite);
854  }
855 
856  free(data);
857  VSdetach(vdid_w);
858  VSdetach(vdid);
859  }
860  free(lonebuf);
861 
862 
863  /* Write Global Attributes */
864  /* ----------------------- */
865  for (i = 0; i < nglobal_attr; i++) {
866  status = SDattrinfo(sd_id_r, i, buffer, &dtype, &count);
867 
868  if (strcmp(buffer, "CoreMetadata.0") == 0) {
869  coremeta = 1;
870  }
871 
872  if (strcmp(buffer, "Number of scans") == 0) {
873  strcpy(buffer, "Number of Scans");
874  }
875 
876  /* printf("%d %s\n", count, buffer);*/
877  data = (char *) calloc(count, DFKNTsize(dtype));
878 
879  status = SDreadattr(sd_id_r, i, (VOIDP) data);
880 
881  if (strcmp(buffer, "Number of Scans") == 0 && last_bad == 1) {
882  i32 = 202;
883  memcpy(data, &i32, sizeof (i32));
884  }
885 
886 
887  if (strcmp(buffer, "Processing Time") == 0) strcpy(data, ydhmsf(now(), 'L'));
888  status = SDsetattr(sd_id_w, buffer, dtype, count, (VOIDP) data);
889  /* printf("%d %d\n", i, status);*/
890 
891  if (strcmp(buffer, "Satellite") == 0) {
892  strcpy(sat, data);
893  }
894 
895  free(data);
896  }
897 
898 
899  if (coremeta == 0) {
900 
901 #if 0
902  strcpy(buf2, getenv("SIMBIOS_ROOT"));
903  strcat(buf2, "/data/modis/static/");
904  strcat(buf2, sat);
905  strcat(buf2, "/coremeta.txt");
906  printf("%s\n", buf2);
907  fp = fopen(buf2, "r");
908 
909  if (fp == 0x0) {
910  printf("\"coremeta.txt\" cannot be found.\n");
911  } else {
912  buffer[0] = 10;
913  buffer[1] = 0;
914 
915  while (fgets(buf2, 80, fp) != NULL) {
916  strcat(buffer, buf2);
917  }
918  if (feof(fp) != 0)
919  status = SDsetattr(sd_id_w, "CoreMetadata.0", DFNT_CHAR,
920  strlen(buffer) + 1, (VOIDP) buffer);
921  }
922 #endif
923 
924  /* Write CoreMetadata (jmg) */
925  strcpy(buffer, "\n");
926  strcat(buffer, "GROUP = INVENTORYMETADATA\n");
927  strcat(buffer, " GROUPTYPE = MASTERGROUP\n");
928 
929  strcat(buffer, " GROUP = ECSDATAGRANULE\n");
930 
931  strcat(buffer, " OBJECT = LOCALGRANULEID\n");
932  strcat(buffer, " NUM_VAL = 1\n");
933 
934 
935  memset(buf2, 0, sizeof (buf2));
936  strcpy(buf2, getenv("L0File"));
937  strcpy(buf2, basename(buf2));
938 
939  strcat(buffer, " VALUE = ");
940  strcat(buffer, "\"");
941  strcat(buffer, buf2);
942  strcat(buffer, "\"\n");
943  strcat(buffer, " END_OBJECT = LOCALGRANULEID\n");
944 
945  strcat(buffer, " OBJECT = PRODUCTIONDATETIME\n");
946  strcat(buffer, " NUM_VAL = 1\n");
947  strcat(buffer, " VALUE = ");
948  strcat(buffer, "\"");
949  strcat(buffer, prodtime);
950  strcat(buffer, "\"\n");
951  strcat(buffer, " END_OBJECT = PRODUCTIONDATETIME\n");
952 
953  strcat(buffer, " OBJECT = DAYNIGHTFLAG\n");
954  strcat(buffer, " NUM_VAL = 1\n");
955  strcat(buffer, " VALUE = ");
956  memset(buf2, 0, sizeof (buf2));
957  if ((i = SDfindattr(sd_id_r, "DAYNIGHTFLAG")) != -1)
958  status = SDreadattr(sd_id_r, i, buf2);
959  strcat(buffer, "\"");
960  strcat(buffer, buf2);
961  strcat(buffer, "\"\n");
962  strcat(buffer, " END_OBJECT = DAYNIGHTFLAG\n");
963 
964  strcat(buffer, " OBJECT = REPROCESSINGACTUAL\n");
965  strcat(buffer, " NUM_VAL = 1\n");
966  strcat(buffer, " VALUE = ");
967  strcat(buffer, "\"processed once\"");
968  strcat(buffer, "\n");
969  strcat(buffer, " END_OBJECT = REPROCESSINGACTUAL\n");
970 
971  strcat(buffer, " END_GROUP = ECSDATAGRANULE\n");
972 
973  strcat(buffer, " GROUP = ORBITCALCULATEDSPATIALDOMAIN\n");
974 
975  strcat(buffer, " OBJECT = ORBITCALCULATEDSPATIALDOMAINCONTAINER\n");
976  strcat(buffer, " CLASS = \"1\"\n");
977  strcat(buffer, " OBJECT = EQUATORCROSSINGDATE\n");
978  strcat(buffer, " CLASS = \"1\"\n");
979  strcat(buffer, " NUM_VAL = 1\n");
980  memset(buf2, 0, sizeof (buf2));
981  if ((i = SDfindattr(sd_id_r, "EQUATORCROSSINGDATE")) != -1)
982  status = SDreadattr(sd_id_r, i, buf2);
983  strcat(buffer, " VALUE = ");
984  strcat(buffer, "\"");
985  strcat(buffer, buf2);
986  strcat(buffer, "\"\n");
987  strcat(buffer, " END_OBJECT = EQUATORCROSSINGDATE\n");
988 
989  strcat(buffer, " OBJECT = EQUATORCROSSINGTIME\n");
990  strcat(buffer, " CLASS = \"1\"\n");
991  strcat(buffer, " NUM_VAL = 1\n");
992  memset(buf2, 0, sizeof (buf2));
993  if ((i = SDfindattr(sd_id_r, "EQUATORCROSSINGTIME")) != -1)
994  status = SDreadattr(sd_id_r, i, buf2);
995  strcat(buffer, " VALUE = ");
996  strcat(buffer, "\"");
997  strcat(buffer, buf2);
998  strcat(buffer, "\"\n");
999  strcat(buffer, " END_OBJECT = EQUATORCROSSINGTIME\n");
1000 
1001  strcat(buffer, " OBJECT = ORBITNUMBER\n");
1002  strcat(buffer, " CLASS = \"1\"\n");
1003  strcat(buffer, " NUM_VAL = 1\n");
1004  strcat(buffer, " VALUE = ");
1005  memset(buf2, 0, sizeof (buf2));
1006  if ((i = SDfindattr(sd_id_r, "ORBITNUMBER")) != -1) {
1007  status = SDreadattr(sd_id_r, i, &i32);
1008  sprintf(buf2, "%6d\n", i32);
1009  } else sprintf(buf2, "%6d\n", -1);
1010  strcat(buffer, buf2);
1011  strcat(buffer, " END_OBJECT = ORBITNUMBER\n");
1012 
1013  strcat(buffer, " OBJECT = EQUATORCROSSINGLONGITUDE\n");
1014  strcat(buffer, " CLASS = \"1\"\n");
1015  strcat(buffer, " NUM_VAL = 1\n");
1016  strcat(buffer, " VALUE = ");
1017  memset(buf2, 0, sizeof (buf2));
1018  if ((i = SDfindattr(sd_id_r, "EQUATORCROSSINGLONGITUDE")) != -1) {
1019  status = SDreadattr(sd_id_r, i, &f32);
1020  sprintf(buf2, "%8.5f\n", f32);
1021  } else sprintf(buf2, "%8.5f\n", 999.0);
1022  strcat(buffer, buf2);
1023  strcat(buffer, " END_OBJECT = EQUATORCROSSINGLONGITUDE\n");
1024 
1025  strcat(buffer, " END_OBJECT = ORBITCALCULATEDSPATIALDOMAINCONTAINER\n");
1026  strcat(buffer, " END_GROUP = ORBITCALCULATEDSPATIALDOMAIN\n");
1027 
1028  strcat(buffer, " GROUP = PGEVERSIONCLASS\n");
1029  strcat(buffer, " OBJECT = PGEVERSION\n");
1030  strcat(buffer, " NUM_VAL = 1\n");
1031  strcat(buffer, " VALUE = ");
1032  strcat(buffer, "\"4.0.3\"\n");
1033  strcat(buffer, " END_OBJECT = PGEVERSION\n");
1034  strcat(buffer, " END_GROUP = PGEVERSIONCLASS\n");
1035 
1036  strcat(buffer, " GROUP = SPATIALDOMAINCONTAINER\n");
1037  strcat(buffer, " GROUP = HORIZONTALSPATIALDOMAINCONTAINER\n");
1038  strcat(buffer, " GROUP = GPOLYGON\n");
1039  strcat(buffer, " OBJECT = GPOLYGONCONTAINER\n");
1040  strcat(buffer, " CLASS = \"1\"\n");
1041 
1042  strcat(buffer, " GROUP = GRINGPOINT\n");
1043  strcat(buffer, " CLASS = \"1\"\n");
1044  strcat(buffer, " OBJECT = GRINGPOINTLONGITUDE\n");
1045  strcat(buffer, " NUM_VAL= 4\n");
1046  strcat(buffer, " CLASS = \"1\"\n");
1047  strcat(buffer, " VALUE = ");
1048  memset(buf2, 0, sizeof (buf2));
1049  sprintf(buf2, "(%8.5f, %8.5f, %8.5f, %8.5f)\n", -1., -1., -1., -1.);
1050  strcat(buffer, buf2);
1051  strcat(buffer, " END_OBJECT= GRINGPOINTLONGITUDE\n");
1052 
1053  strcat(buffer, " OBJECT = GRINGPOINTLATITUDE\n");
1054  strcat(buffer, " NUM_VAL= 4\n");
1055  strcat(buffer, " CLASS = \"1\"\n");
1056  strcat(buffer, " VALUE = ");
1057  memset(buf2, 0, sizeof (buf2));
1058  sprintf(buf2, "(%8.5f, %8.5f, %8.5f, %8.5f)\n", -1., -1., -1., -1.);
1059  strcat(buffer, buf2);
1060  strcat(buffer, " END_OBJECT= GRINGPOINTLATITUDE\n");
1061 
1062  strcat(buffer, " OBJECT = GRINGPOINTSEQUENCENO\n");
1063  strcat(buffer, " NUM_VAL= 4\n");
1064  strcat(buffer, " CLASS = \"1\"\n");
1065  strcat(buffer, " VALUE = ");
1066  memset(buf2, 0, sizeof (buf2));
1067  sprintf(buf2, "(%4d, %4d, %4d, %4d)\n", 1, 2, 3, 4);
1068  strcat(buffer, buf2);
1069  strcat(buffer, " END_OBJECT= GRINGPOINTSEQUENCENO\n");
1070  strcat(buffer, " END_GROUP = GRINGPOINT\n");
1071 
1072  strcat(buffer, " GROUP = GRING\n");
1073  strcat(buffer, " CLASS = \"1\"\n");
1074  strcat(buffer, " OBJECT = EXCLUSIONGRINGFLAG\n");
1075  strcat(buffer, " NUM_VAL= 1\n");
1076  strcat(buffer, " CLASS = \"1\"\n");
1077  strcat(buffer, " VALUE = ");
1078  memset(buf2, 0, sizeof (buf2));
1079  if ((i = SDfindattr(sd_id_r, "EXCLUSIONGRINGFLAG")) != -1)
1080  status = SDreadattr(sd_id_r, i, buf2);
1081  strcat(buffer, "\"");
1082  strcat(buffer, buf2);
1083  strcat(buffer, "\"\n");
1084  strcat(buffer, " END_OBJECT= EXCLUSIONGRINGFLAG\n");
1085  strcat(buffer, " END_GROUP = GRING\n");
1086  strcat(buffer, " END_OBJECT = GPOLYGONCONTAINER\n");
1087  strcat(buffer, " END_GROUP = GPOLYGON\n");
1088  strcat(buffer, " END_GROUP = HORIZONTALSPATIALDOMAINCONTAINER\n");
1089  strcat(buffer, " END_GROUP = SPATIALDOMAINCONTAINER\n");
1090 
1091  strcat(buffer, " GROUP = RANGEDATETIME\n");
1092 
1093  strcat(buffer, " OBJECT = RANGEENDINGDATE\n");
1094  strcat(buffer, " NUM_VAL = 1\n");
1095  strcat(buffer, " VALUE = ");
1096  memset(buf2, 0, sizeof (buf2));
1097  if ((i = SDfindattr(sd_id_r, "RANGEENDINGDATE")) != -1)
1098  status = SDreadattr(sd_id_r, i, buf2);
1099  strcat(buffer, "\"");
1100  strcat(buffer, buf2);
1101  strcat(buffer, "\"\n");
1102  strcat(buffer, " END_OBJECT = RANGEENDINGDATE\n");
1103 
1104  strcat(buffer, " OBJECT = RANGEENDINGTIME\n");
1105  strcat(buffer, " NUM_VAL = 1\n");
1106  strcat(buffer, " VALUE = ");
1107  memset(buf2, 0, sizeof (buf2));
1108  if ((i = SDfindattr(sd_id_r, "RANGEENDINGTIME")) != -1)
1109  status = SDreadattr(sd_id_r, i, buf2);
1110  strcat(buffer, "\"");
1111  strcat(buffer, buf2);
1112  strcat(buffer, "\"\n");
1113  strcat(buffer, " END_OBJECT = RANGEENDINGTIME\n");
1114 
1115  strcat(buffer, " OBJECT = RANGEBEGINNINGDATE\n");
1116  strcat(buffer, " NUM_VAL = 1\n");
1117  strcat(buffer, " VALUE = ");
1118  memset(buf2, 0, sizeof (buf2));
1119  if ((i = SDfindattr(sd_id_r, "RANGEBEGINNINGDATE")) != -1)
1120  status = SDreadattr(sd_id_r, i, buf2);
1121  strcat(buffer, "\"");
1122  strcat(buffer, buf2);
1123  strcat(buffer, "\"\n");
1124  strcat(buffer, " END_OBJECT = RANGEBEGINNINGDATE\n");
1125 
1126  strcat(buffer, " OBJECT = RANGEBEGINNINGTIME\n");
1127  strcat(buffer, " NUM_VAL = 1\n");
1128  strcat(buffer, " VALUE = ");
1129  memset(buf2, 0, sizeof (buf2));
1130  if ((i = SDfindattr(sd_id_r, "RANGEBEGINNINGTIME")) != -1)
1131  status = SDreadattr(sd_id_r, i, buf2);
1132  strcat(buffer, "\"");
1133  strcat(buffer, buf2);
1134  strcat(buffer, "\"\n");
1135  strcat(buffer, " END_OBJECT = RANGEBEGINNINGTIME\n");
1136 
1137  strcat(buffer, " END_GROUP = RANGEDATETIME\n");
1138 
1139  strcat(buffer, " GROUP = ADDITIONALATTRIBUTES\n");
1140 
1141  strcat(buffer, " OBJECT = ADDITIONALATTRIBUTESCONTAINER\n");
1142  strcat(buffer, " CLASS = \"1\"\n");
1143 
1144  strcat(buffer, " OBJECT = ADDITIONALATTRIBUTENAME\n");
1145  strcat(buffer, " CLASS = \"1\"\n");
1146  strcat(buffer, " NUM_VAL = 1\n");
1147  strcat(buffer, " VALUE = GRANULENUMBER\n");
1148  strcat(buffer, " END_OBJECT = ADDITIONALATTRIBUTENAME\n");
1149 
1150  strcat(buffer, " END_OBJECT = ADDITIONALATTRIBUTESCONTAINER\n");
1151  strcat(buffer, " END_GROUP = ADDITIONALATTRIBUTES\n");
1152 
1153 
1154  strcat(buffer, " GROUP = CollectionDescriptionClass\n");
1155  strcat(buffer, " OBJECT = SHORTNAME\n");
1156  strcat(buffer, " NUM_VAL = 1\n");
1157  strcat(buffer, " VALUE = ");
1158  strcat(buffer, "\"MYD01SS\"\n");
1159  strcat(buffer, " END_OBJECT = SHORTNAME\n");
1160  strcat(buffer, " END_GROUP = CollectionDescriptionClass\n");
1161 
1162  strcat(buffer, "END_GROUP = INVENTORYMETADATA\n");
1163  strcat(buffer, "END\n");
1164  strcat(buffer, "");
1165 
1166  status = SDsetattr(sd_id_w, "CoreMetadata.0", DFNT_CHAR, strlen(buffer) + 1,
1167  (VOIDP) buffer);
1168 
1169  }
1170 
1171 
1172  SDend(sd_id_r);
1173  Hclose(HDFfid_r);
1174  Vend(HDFfid_r);
1175 
1176  SDend(sd_id_w);
1177  Hclose(HDFfid_w);
1178  Vend(HDFfid_w);
1179 
1180  return 0;
1181 }
1182 
1183 int compare_int16(const void * p1, const void * p2) {
1184  int16 a = *((int16 *) p1);
1185  int16 b = *((int16 *) p2);
1186 
1187  if (a > b)
1188  return 1;
1189  else if (a < b)
1190  return -1;
1191  else
1192  return 0;
1193 }
1194 
1195 int compare_float32(const void * p1, const void * p2) {
1196  float32 a = *((float32 *) p1);
1197  float32 b = *((float32 *) p2);
1198 
1199  if (a > b)
1200  return 1;
1201  else if (a < b)
1202  return -1;
1203  else
1204  return 0;
1205 }
1206 
1207 
1208 
1209 
1210 /*
1211 qsort ((void *) data_1km, nelem, sizeof(int16), compare_int16);
1212 for (k=0; k<nelem; k++) if (data_1km[k] > 0) break;
1213 for (l=0; l<nelem; l++) if (data_1km[nelem-l-1] < 32767) break;
1214 printf("%f\n", (float32) data_1km[(k+(nelem-l-1))/2]);
1215 
1216 qsort ((void *) median_250_500, nelem, sizeof(float32), compare_float32);
1217 for (k=0; k<nelem; k++) if (median_250_500[k] > 0) break;
1218 for (l=0; l<nelem; l++) if (median_250_500[nelem-l-1] < 32767) break;
1219 printf("%f\n", median_250_500[(k+(nelem-l-1))/2]);
1220 
1221 f32 = 1. / (res_b[0]*data_1km[(k+(nelem-l-1))/2]/median_250_500[(k+(nelem-l-1))/2]+
1222  1000*res_b[1]/median_250_500[(k+(nelem-l-1))/2] + const_b);
1223 status = SDsetattr(sd_id_w, "B_TC_RS", DFNT_FLOAT32, 1, (VOIDP) &f32);
1224 printf("%f\n",f32);
1225  */
1226 
1227 
1228 
1229 /*
1230 if (masking == 1) {
1231  if (argc == 6) {
1232  scan_incr = atoi(argv[4]);
1233  pixl_incr = atoi(argv[5]);
1234  } else {
1235  printf("Scan and Pixel Increments must be defined for masking.\n");
1236  exit(-1);
1237  }
1238 }
1239  */
1240 
1241 
1242 #if 0
1243 /* Flag skipped pixels and scans (if masking on) */
1244 /* --------------------------------------------- */
1245 if (strncmp(buffer, "EV", 2) == 0 && masking == 1) {
1246  printf("Masking %s\n", buffer);
1247  for (k = 0; k < dims[0]; k++) {
1248  if (k % scan_incr != 0) {
1249  memset(&data[k * 2 * dims[1] * dims[2]], 255, 2 * dims[1] * dims[2]);
1250  } else {
1251 
1252  for (l = 0; l < dims[1]; l++)
1253  for (i = 0; i < dims[2]; i++)
1254  if ((i % pixl_incr) != 0)
1255  memset(&data[k * 2 * dims[1] * dims[2] + l * 2 * dims[2] + i * 2], 255, 2);
1256  }
1257  }
1258 }
1259 #endif
ydhmsf
char * ydhmsf(double dtime, char zone)
Definition: ydhmsf.c:12
cubeio::int16
integer, parameter int16
Definition: cubeio.f90:3
compare_float32
int compare_float32(const void *p1, const void *p2)
Definition: l1asubset_modis.c:1195
j
int j
Definition: decode_rs.h:73
status
int status
Definition: l1_czcs_hdf.c:32
compare_int16
int compare_int16(const void *p1, const void *p2)
Definition: l1asubset_modis.c:1183
NULL
#define NULL
Definition: decode_rs.h:63
tm
float tm[MODELMAX]
Definition: atrem_corl1v2.h:248
usage
void usage(char *progname)
Definition: l1asubset_modis.c:32
VERSION
#define VERSION
Definition: l1asubset_modis.c:30
basename
#define basename(s)
Definition: l1asubset_modis.c:27
seadasutils.aquarius_utils.start
start
Definition: aquarius_utils.py:27
data
no change in intended resolving MODur00064 Corrected handling of bad ephemeris attitude data
Definition: HISTORY.txt:356
b
data_t b[NROOTS+1]
Definition: decode_rs.h:77
dtype
dtype
Definition: DDataset.hpp:31
time
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
now
u5 which has been done in the LOCALGRANULEID metadata should have an extension NRT It is requested to identify the NRT production Changes from v6 which may affect scientific the sector rotation may actually occur during one of the scans earlier than the one where it is first reported As a the b1 values are about the LOCALGRANULEID metadata should have an extension NRT It is requested to identify the NRT to fill pixels affected by dead subframes with a special value Output the metadata of noisy and dead subframe Dead Subframe EV and Detector Quality Flag2 Removed the function call of Fill_Dead_Detector_SI to stop interpolating SI values for dead but also for all downstream products for science test only Changes from v5 which will affect scientific to conform to MODIS requirements Removed the Mixed option from the ScanType in the code because the L1A Scan Type is never Mixed Changed for ANSI C compliance and comments to better document the fact that when the HDF_EOS metadata is stricly the and products are off by and in the track respectively Corrected some misspelling of RCS swir_oob_sending_detector to the Reflective LUTs to enable the SWIR OOB correction detector so that if any of the sending detectors becomes noisy or non near by good detectors from the same sending band can be specified as the substitute in the new look up table Code change for adding an additional dimension of mirror side to the Band_21_b1 LUT to separate the coefficient of the two mirror sides for just like other thermal emissive so that the L1B code can calibrate Band scan to scan with mirror side dependency which leads better calibration result Changes which do not affect scientific when the EV data are not provided in this Crosstalk Correction will not be performed to the Band calibration data Changes which do not affect scientific and BB_500m in L1A Logic was added to turn off the or to spatial aggregation processes and the EV_250m_Aggr1km_RefSB and EV_500m_Aggr1km_RefSB fields were set to fill values when SDSs EV_250m and EV_500m are absent in L1A file Logic was added to skip the processing and turn off the output of the L1B QKM and HKM EV data when EV_250m and EV_500m are absent from L1A In this the new process avoids accessing and reading the and L1A EV skips and writing to the L1B and EV omits reading and subsampling SDSs from geolocation file and writing them to the L1B and omits writing metadata to L1B and EV and skips closing the L1A and L1B EV and SDSs Logic was added to turn off the L1B OBC output when the high resolution OBC SDSs are absent from L1A This is accomplished by skipping the openning the writing of metadata and the closing of the L1B OBC hdf which is Bit in the scan by scan bit QA has been changed Until now
Definition: HISTORY.txt:361
rank
Extra metadata that will be written to the HDF4 file l2prod rank
Definition: HOWTO_Add_a_product.txt:80
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
main
int main(int argc, char *argv[])
Definition: l1asubset_modis.c:42
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")
a
PGE01 indicating that PGE02 PGE01 V6 for and PGE01 V2 for MOD03 were used to produce the granule By convention adopted in all MODIS Terra PGE02 code versions are The fourth digit of the PGE02 version denotes the LUT version used to produce the granule The source of the metadata environment variable ProcessingCenter was changed from a QA LUT value to the Process Configuration A sign used in error in the second order term was changed to a
Definition: HISTORY.txt:424
viirs_ler_luts::p1
real p1
Definition: viirs_ler_luts.f95:37
k
int k
Definition: decode_rs.h:73
npix
int npix
Definition: get_cmp.c:27
f32
float32 f32
Definition: l2bin.cpp:104
count
int count
Definition: decode_rs.h:79


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

Level 1 & 2 Browser Help

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 Policy

Freedom of Information Act

USA.gov

Responsible NASA Official: Sean Bailey
Curator: OceanColor Webmaster