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_octs_netcdf.cpp
Go to the documentation of this file.
1 
2 #include <stdio.h>
3 #include <stdlib.h>
4 #include <string.h>
5 #include <ctype.h>
6 #include <hdf4utils.h>
7 #include "l1.h"
8 #include "filehdr_struc.h"
9 extern "C" {
10 
11  #include "l1_octs.h"
12 }
13 #include "l1_octs_netcdf.h"
14 #include <netcdf>
15 #include <iostream>
16 
17 #include <hdf.h>
18 #include <mfhdf.h>
19 
20 #define RADIUS 6378.137 /* Earth radius in km */
21 #define FL 1.0/298.257 /* Earth flatening factor */
22 #define NR 560 /* # rows */
23 #define NC 30 /* max # columns */
24 
25 #define MAXOCLIN 6700 /* max # lines */
26 #define MAXOCPIX 2218 /* max # pixels */
27 #define MAXOCARR 10000
28 #define NOCBANDS 8
29 
30 
31 using namespace std;
32 using namespace netCDF;
33 using namespace netCDF::exceptions;
34 
35 static int32_t msec_start;
36 static int32_t msec[MAXOCLIN];
37 static float tilt[MAXOCLIN];
38 static int16_t year, day, nline, npix, sline;
39 static int32_t spix;
40 static int16_t maxday;
41 static float solz1, sola1, senz1, sena1;
42 static float *lon, *lat, *senz, *sena, *solz, *sola;
43 static int32_t ncol;
44 static int32_t scansPerScene;
45 static int32_t linesPerScan;
46 static int16_t *tilt_flag;
47 static int16_t *gain;
48 static float *inst_temp;
49 static int16_t *miss_qual;
50 static int16_t num_tab;
51 static int32_t extr_line_offset = 0;
52 static int16_t *start_scan;
53 static int16_t last_table = 0;
54 static int16_t sample_table[3][8][2][400][2];
55 
56 /* ------------------------------------------------------ */
57 /* CalcViewAngle() - calculates the solar and sensor view */
58 /* angles from one geodetic longitude and */
59 /* latitude, position vector, and solar unit */
60 /* vector. */
61 /* */
62 /* INPUT DESCRIPTION */
63 /* ----- ----------- */
64 /* lon longitude of pixel in degrees */
65 /* lat latitude of pixel in degrees */
66 /* pos[3] vector from earth center to spacecraft in */
67 /* ECR, km */
68 /* usun[3] unit vector from earth center to sun in */
69 /* ECR */
70 /* */
71 /* OUTPUT DESCRIPTION */
72 /* ------ ----------- */
73 /* solz solar zenith angle of pixel in degrees */
74 /* sola solar azimuth angle of pixel in degrees */
75 /* senz sensor zenith angle of pixel in degrees */
76 /* sena sensor azimuth angle of pixel in degrees */
77 /* */
78 
79 /* ------------------------------------------------------ */
80 static int CalcViewAngle_nc(float lon1, float lat1, float pos[3],
81  float usun[3]) {
82  float rmat[3][3]; /* rotation matrix */
83  float rlon, rlat;
84  float up[3], upxy, ea[3], no[3];
85  float phi, R;
86  float gvec[3], scvec[3], senl[3], sunl[3];
87  int i, j;
88 
89  rlon = lon1 * PI / 180.0;
90  rlat = lat1 * PI / 180.0;
91 
92  /* First, we must define the axes (in ECR space) of a
93  pixel-local coordinate system, with the z-axis along
94  the geodetic pixel normal, x-axis pointing east, and
95  y-axis pointing north. */
96  up[0] = cos(rlat) * cos(rlon);
97  up[1] = cos(rlat) * sin(rlon);
98  up[2] = sin(rlat);
99  upxy = sqrt(up[0] * up[0 ] + up[1] * up[1]);
100  ea[0] = -up[1] / upxy;
101  ea[1] = up[0] / upxy;
102  ea[2] = 0.0;
103  /* calculate cross product of up and ea */
104  no[0] = up[1] * ea[2] - ea[1] * up[2];
105  no[1] = up[2] * ea[0] - ea[2] * up[0];
106  no[2] = up[0] * ea[1] - ea[0] * up[1];
107 
108 
109  /* Now we need the pixel-to-spacecraft vector in ECR,
110  Compute geocentric pixel location vector in km and subtract
111  from spacecraft position. */
112  /* geocentric latitude */
113  phi = atan(tan(rlat)*(1 - FL)*(1 - FL));
114  /* dist to Earth surface */
115  R = RADIUS * (1 - FL) / sqrt(1 - (2 - FL) * FL * (cos(phi) * cos(phi)));
116  gvec[0] = R * cos(phi) * cos(rlon);
117  gvec[1] = R * cos(phi) * sin(rlon);
118  gvec[2] = R * sin(phi);
119  for (i = 0; i < 3; i++) {
120  scvec[i] = pos[i] - gvec[i];
121  }
122 
123  /* Now we can transform the pixel-to-spacecraft and Sun
124  vectors into the local frame. */
125  for (i = 0; i < 3; i++) {
126  rmat[0][i] = ea[i];
127  rmat[1][i] = no[i];
128  rmat[2][i] = up[i];
129  }
130  for (i = 0; i < 3; i++) {
131  senl[i] = 0;
132  sunl[i] = 0;
133  for (j = 0; j < 3; j++) {
134  senl[i] = senl[i] + rmat[i][j] * scvec[j];
135  sunl[i] = sunl[i] + rmat[i][j] * usun[j];
136  }
137  }
138 
139  /* Compute the solar zenith and azimuth */
140  solz1 = RADEG * atan(sqrt(sunl[0] * sunl[0] + sunl[1] * sunl[1]) / sunl[2]);
141  if (solz1 < 0.0) solz1 += 180.0;
142 
143  if (solz1 > 0.05) {
144  sola1 = RADEG * (atan2(sunl[0], sunl[1]));
145  } else {
146  sola1 = 0.0;
147  }
148  if (sola1 < 0.0) {
149  sola1 = sola1 + 360.0;
150  }
151 
152  /* Compute the sensor zenith and azimuth */
153  senz1 = RADEG * atan(sqrt(senl[0] * senl[0] + senl[1] * senl[1]) / senl[2]);
154  if (senz1 < 0.0) senz1 += 180.0;
155  if (senz1 > 0.05) {
156  sena1 = RADEG * (atan2(senl[0], senl[1]));
157  } else {
158  sena1 = 0.0;
159  }
160  if (sena1 < 0.0) {
161  sena1 = sena1 + 360.0;
162  }
163  return (0);
164 }
165 
166 
167 
168 /* ------------------------------------------------------ */
169 /* navigation() - generates navigation info interpolated */
170 /* for each pixel, using info from a NASDA */
171 /* OCTS HDF file */
172 /* */
173 
174 /* ------------------------------------------------------ */
175 extern "C" void navigation_netcdf(NcFile &dataFile) {
176  int16_t det;
177 
178  int16_t *pxl;
179  float *xctl;
180  float *inlon, *inlat;
181  float *insolz, *insola;
182  float *insenz, *insena;
183  float *usun, *pos;
184  int32_t *row;
185  int32_t *indx;
186  float *yctl;
187  float *in1, *in2;
188 
189  float iusun[3], ipos[3];
190  float usun_sum;
191  int32_t i, j, k, l;
192  int32_t eline, nlon, nlat;
193  int32_t ilat[MAXOCLIN], ilon[MAXOCLIN];
194  float out1[MAXOCLIN], out2[MAXOCLIN];
195  int32_t jout;
196  float *spl_aux;
197  float *x_ctl_ll, *y_ctl_ll, *z_ctl_ll;
198  float *x_ctl_sol, *y_ctl_sol, *z_ctl_sol;
199  float *x_ctl_sen, *y_ctl_sen, *z_ctl_sen;
200  float * in_ptr[3][3], *out_ptr[3][2];
201 
202  inlon = (float *) calloc(ncol*scansPerScene, sizeof (float));
203  inlat = (float *) calloc(ncol*scansPerScene, sizeof (float));
204  insolz = (float *) calloc(ncol*scansPerScene, sizeof (float));
205  insola = (float *) calloc(ncol*scansPerScene, sizeof (float));
206  insenz = (float *) calloc(ncol*scansPerScene, sizeof (float));
207  insena = (float *) calloc(ncol*scansPerScene, sizeof (float));
208  usun = (float *) calloc(3 * scansPerScene, sizeof (float));
209  pos = (float *) calloc(3 * scansPerScene, sizeof (float));
210  row = (int32_t *) calloc(scansPerScene, sizeof (int));
211  pxl = (int16_t *) calloc(ncol, sizeof (int16_t));
212  xctl = (float *) calloc(ncol, sizeof (float));
213  indx = (int32_t *) calloc(scansPerScene, sizeof (int));
214  yctl = (float *) calloc(scansPerScene, sizeof (float));
215  in1 = (float *) calloc(scansPerScene, sizeof (float));
216  in2 = (float *) calloc(scansPerScene, sizeof (float));
217 
218  /* read the data sets */
219  NcVar tempVar = dataFile.getVar("det");
220  tempVar.getVar(&det);
221  tempVar = dataFile.getVar("pxl"); /* control pt. pixels */
222  tempVar.getVar(pxl);
223  tempVar = dataFile.getVar("lon"); /* geodetic longitude */
224  tempVar.getVar(inlon);
225  tempVar = dataFile.getVar("lat"); /* geodetic latitude */
226  tempVar.getVar(inlat);
227  tempVar = dataFile.getVar("sun_ref"); /* geocentric ECR solar refercence vector */
228  tempVar.getVar(usun);
229  tempVar = dataFile.getVar("orb_vec"); /* geocentric ECR S/C position */
230  tempVar.getVar(pos);
231 
232  det = det - 1;
233  for (i = 0; i < ncol; i++) {
234  pxl[i] = pxl[i] - 1;
235  }
236 
237  /* convert sun vector to unit vector */
238  for (i = 0; i < scansPerScene; i++) {
239  usun_sum = 0;
240  for (j = 0; j < 3; j++) {
241  usun_sum = usun_sum + usun[3 * i + j] * usun[3 * i + j];
242  }
243  for (j = 0; j < 3; j++) {
244  usun[3 * i + j] = usun[3 * i + j] / sqrt(usun_sum);
245  }
246  }
247 
248  /* define control point grid */
249  if (want_verbose)
250  printf("scansPerScene,ncol = %d %d \n", scansPerScene, ncol);
251  for (i = 0; i < ncol; i++) {
252  xctl[i] = (float) pxl[i];
253  }
254  for (i = 0; i < scansPerScene; i++) {
255  row[i] = i * linesPerScan + det;
256  }
257 
258  /* define output grid */
259  eline = sline + nline - 1;
260  nlon = npix; /* # of output pixs */
261  if (nline < linesPerScan * scansPerScene) {
262  nlat = nline; /* # of output lines */
263  } else nlat = linesPerScan*scansPerScene;
264 
265  for (i = 0; i < nlon; i++) {
266  ilon[i] = i + spix; /* pix #'s of output grid */
267  }
268  for (i = 0; i < nlat; i++) {
269  ilat[i] = i + sline; /* line #'s of output grid */
270  }
271 
272  /* create index of output lines to control point lines */
273  k = 0;
274  for (i = 0; i < nlat; i++) {
275  for (j = 0; j < scansPerScene; j++) {
276  if (ilat[i] == row[j]) {
277  indx[k] = i;
278  k++;
279  }
280  }
281  }
282 
283  for (i = 0; i < nlat / linesPerScan; i++) {
284  yctl[i] = (float) indx[i];
285  }
286 
287 
288  /* compute solar and sensor zenith and azimuth at input
289  control points from info provided */
290  if (want_verbose)
291  printf("Computing view angles at input control points\n");
292  for (i = 0; i < scansPerScene; i++) {
293  for (j = 0; j < ncol; j++) {
294  for (k = 0; k < 3; k++) {
295  ipos[k] = pos[3 * i + k];
296  iusun[k] = usun[3 * i + k];
297  }
298  CalcViewAngle_nc(inlon[i * ncol + j], inlat[i * ncol + j], ipos, iusun);
299  insolz[i * ncol + j] = solz1;
300  insola[i * ncol + j] = sola1;
301  insenz[i * ncol + j] = senz1;
302  insena[i * ncol + j] = sena1;
303  }
304  }
305 
306 
307  /* Compute unit vectors from lon/lat of control points */
308  x_ctl_ll = (float *) calloc(ncol*scansPerScene, sizeof (float));
309  y_ctl_ll = (float *) calloc(ncol*scansPerScene, sizeof (float));
310  z_ctl_ll = (float *) calloc(ncol*scansPerScene, sizeof (float));
311 
312  x_ctl_sol = (float *) calloc(ncol*scansPerScene, sizeof (float));
313  y_ctl_sol = (float *) calloc(ncol*scansPerScene, sizeof (float));
314  z_ctl_sol = (float *) calloc(ncol*scansPerScene, sizeof (float));
315 
316  x_ctl_sen = (float *) calloc(ncol*scansPerScene, sizeof (float));
317  y_ctl_sen = (float *) calloc(ncol*scansPerScene, sizeof (float));
318  z_ctl_sen = (float *) calloc(ncol*scansPerScene, sizeof (float));
319 
320 
321  for (i = 0; i < scansPerScene; i++) {
322  for (j = 0; j < ncol; j++) {
323  inlon[i * ncol + j] = inlon[i * ncol + j] / RADEG;
324  inlat[i * ncol + j] = inlat[i * ncol + j] / RADEG;
325 
326  x_ctl_ll[i * ncol + j] = cos(inlat[i * ncol + j]) * cos(inlon[i * ncol + j]);
327  y_ctl_ll[i * ncol + j] = cos(inlat[i * ncol + j]) * sin(inlon[i * ncol + j]);
328  z_ctl_ll[i * ncol + j] = sin(inlat[i * ncol + j]);
329 
330 
331  insola[i * ncol + j] = insola[i * ncol + j] / RADEG;
332  insolz[i * ncol + j] = insolz[i * ncol + j] / RADEG;
333 
334  x_ctl_sol[i * ncol + j] = cos(insolz[i * ncol + j]) * cos(insola[i * ncol + j]);
335  y_ctl_sol[i * ncol + j] = cos(insolz[i * ncol + j]) * sin(insola[i * ncol + j]);
336  z_ctl_sol[i * ncol + j] = sin(insolz[i * ncol + j]);
337 
338 
339  insena[i * ncol + j] = insena[i * ncol + j] / RADEG;
340  insenz[i * ncol + j] = insenz[i * ncol + j] / RADEG;
341 
342  x_ctl_sen[i * ncol + j] = cos(insenz[i * ncol + j]) * cos(insena[i * ncol + j]);
343  y_ctl_sen[i * ncol + j] = cos(insenz[i * ncol + j]) * sin(insena[i * ncol + j]);
344  z_ctl_sen[i * ncol + j] = sin(insenz[i * ncol + j]);
345  }
346  }
347 
348  in_ptr[0][0] = x_ctl_ll;
349  in_ptr[0][1] = y_ctl_ll;
350  in_ptr[0][2] = z_ctl_ll;
351  in_ptr[1][0] = x_ctl_sol;
352  in_ptr[1][1] = y_ctl_sol;
353  in_ptr[1][2] = z_ctl_sol;
354  in_ptr[2][0] = x_ctl_sen;
355  in_ptr[2][1] = y_ctl_sen;
356  in_ptr[2][2] = z_ctl_sen;
357 
358  out_ptr[0][0] = lon;
359  out_ptr[0][1] = lat;
360  out_ptr[1][0] = sola;
361  out_ptr[1][1] = solz;
362  out_ptr[2][0] = sena;
363  out_ptr[2][1] = senz;
364 
365 
366  /* we now have all the info at each control point, so we
367  can interpolate to all pixels, all lines */
368  /* interpolate angles across each control point line */
369  if (want_verbose)
370  printf("Interpolating rows for longitude/azimuth\n");
371  spl_aux = (float *) calloc(ncol, sizeof (float));
372 
373  for (i = 0; i < scansPerScene; i++) {
374  jout = row[i] - sline;
375  if ((row[i] >= sline) && (row[i] <= eline)) {
376  for (l = 0; l < 3; l++) {
377  spline(xctl, in_ptr[l][0] + i*ncol, ncol, 1e30, 1e30, spl_aux);
378  for (j = 0; j < nlon; j++)
379  splint(xctl, in_ptr[l][0] + i * ncol, spl_aux, ncol,
380  (float) ilon[j], out_ptr[l][0] + jout * npix + j);
381 
382  spline(xctl, in_ptr[l][1] + i*ncol, ncol, 1e30, 1e30, spl_aux);
383  for (j = 0; j < nlon; j++)
384  splint(xctl, in_ptr[l][1] + i * ncol, spl_aux, ncol,
385  (float) ilon[j], out_ptr[l][1] + jout * npix + j);
386  }
387  }
388  }
389  free(spl_aux);
390 
391 
392  /* fill missing lines by interpolating columns */
393  if (want_verbose)
394  printf("Interpolating columns for longitude/azimuth\n");
395  spl_aux = (float *) calloc(nlat / linesPerScan, sizeof (float));
396 
397  for (i = 0; i < nlon; i++) {
398  for (l = 0; l < 3; l++) {
399  for (k = 0; k < nlat / linesPerScan; k++) {
400  in1[k] = *(out_ptr[l][0] + indx[k] * npix + i);
401  in2[k] = *(out_ptr[l][1] + indx[k] * npix + i);
402  }
403  spline(yctl, in1, nlat / linesPerScan, 1e30, 1e30, spl_aux);
404  for (j = 0; j < nlat; j++)
405  splint(yctl, in1, spl_aux, nlat / linesPerScan, (float) j, (float *) &out1[j]);
406 
407  spline(yctl, in2, nlat / linesPerScan, 1e30, 1e30, spl_aux);
408  for (j = 0; j < nlat; j++)
409  splint(yctl, in2, spl_aux, nlat / linesPerScan, (float) j, (float *) &out2[j]);
410 
411  for (j = 0; j < nlat; j++) {
412  *(out_ptr[l][0] + j * npix + i) = atan2(out2[j], out1[j]) * RADEG;
413  if (l >= 1 && *(out_ptr[l][0] + j * npix + i) < 0) {
414  *(out_ptr[l][0] + j * npix + i) += 360;
415  }
416  }
417  }
418  }
419  free(spl_aux);
420 
421 
422  if (want_verbose)
423  printf("Interpolating rows for latitude/zenith\n");
424  spl_aux = (float *) calloc(ncol, sizeof (float));
425 
426  for (i = 0; i < scansPerScene; i++) {
427  jout = row[i] - sline;
428  if ((row[i] >= sline) && (row[i] <= eline)) {
429  for (l = 0; l < 3; l++) {
430  spline(xctl, in_ptr[l][2] + i*ncol, ncol, 1e30, 1e30, spl_aux);
431  for (j = 0; j < nlon; j++)
432  splint(xctl, in_ptr[l][2] + i * ncol, spl_aux, ncol,
433  (float) ilon[j], out_ptr[l][1] + jout * npix + j);
434  }
435  }
436  }
437  free(spl_aux);
438 
439 
440  /* fill missing lines by interpolating columns */
441  if (want_verbose)
442  printf("Interpolating columns for latitude/zenith\n");
443  spl_aux = (float *) calloc(nlat / linesPerScan, sizeof (float));
444 
445  for (i = 0; i < nlon; i++) {
446  for (l = 0; l < 3; l++) {
447  for (k = 0; k < nlat / linesPerScan; k++) {
448  in1[k] = *(out_ptr[l][1] + indx[k] * npix + i);
449  }
450  spline(yctl, in1, nlat / linesPerScan, 1e30, 1e30, spl_aux);
451  for (j = 0; j < nlat; j++)
452  splint(yctl, in1, spl_aux, nlat / linesPerScan, (float) j, (float *) &out1[j]);
453 
454  for (j = 0; j < nlat; j++) {
455  *(out_ptr[l][1] + j * npix + i) = asin(out1[j]) * RADEG;
456  }
457  }
458  }
459  free(spl_aux);
460 
461 
462  free(x_ctl_ll);
463  free(y_ctl_ll);
464  free(z_ctl_ll);
465  free(x_ctl_sol);
466  free(y_ctl_sol);
467  free(z_ctl_sol);
468  free(x_ctl_sen);
469  free(y_ctl_sen);
470  free(z_ctl_sen);
471 
472  free(inlon);
473  free(inlat);
474  free(insolz);
475  free(insola);
476  free(insenz);
477  free(insena);
478  free(usun);
479  free(pos);
480  free(row);
481  free(pxl);
482  free(xctl);
483  free(indx);
484  free(yctl);
485  free(in1);
486  free(in2);
487 }
488 
489 /* ------------------------------------------------------ */
490 /* openl1_octs() - opens a OCTS HDF level-1 file */
491 /* for reading, if not already opened. */
492 /* Reads the global attributes. */
493 /* Calculates the subscene info, if */
494 /* necessary, ane gets the navigation */
495 /* info. */
496 /* */
497 /* CALLS: getHDFattr() */
498 /* navigation() */
499 /* rdSDS() */
500 /* */
501 
502 /* ------------------------------------------------------ */
503 extern "C" int32_t openl1_octs_netcdf(filehandle *l1file) {
504  int32_t i, j;
505  int32_t pixPerScan;
506  int32_t linesPerScene;
507  int32_t msec_temp1[MAXOCARR], msec_temp2[MAXOCARR];
508  char buf[32];
509 
510  // Reading some Global Attributes and Dimentions
511  try {
512  NcFile dataFile(l1file->name, NcFile::read);
513  char tempTimeCoverage[27];
514 
515  // Start Time
516  dataFile.getAtt("time_coverage_start").getValues((void*) &tempTimeCoverage);
517 
518  // temp to hold the value and then reassign it to the 16 bit int
519  int32_t tempDay, tempYear;
520  isodate2ydmsec(tempTimeCoverage, &tempYear, &tempDay, &msec_start);
521  day = tempDay;
522  year = tempYear;
523 
524  /* get attributes pertaining to the full scene */
525  pixPerScan = dataFile.getDim("nsamp").getSize();
526  scansPerScene = dataFile.getDim("rec").getSize();
527  ncol = dataFile.getDim("pxls").getSize();
528  linesPerScan = 2;
529  dataFile.getAtt("orbit_node_longitude").getValues((void*) &l1file->orbit_node_lon);
530  dataFile.getAtt("orbit_number").getValues((void*) &l1file->orbit_number);
531  dataFile.getAtt("node_crossing_time").getValues((char*) &buf);
532  reform_octs_time(buf);
533  l1file->node_crossing_time = zulu2unix(buf);
534 
535  linesPerScene = scansPerScene * linesPerScan;
536 
537  /* Adjust spix if this is an extract file */
538  if (!dataFile.getAtt("extract_pixel_start").isNull()) {
539  dataFile.getAtt("extract_pixel_start").getValues((void*) &spix);
540  dataFile.getAtt("extract_line_start").getValues((void*) &extr_line_offset);
541  printf("File is Level-1A extract starting on line %d.\n", extr_line_offset + 1);
542  } else {
543  extr_line_offset = 0;
544  spix = 0;
545  }
546  sline = 0;
547  nline = linesPerScene;
548  npix = pixPerScan;
549 
550 
551  /* Check that number of scene lines not greater than allocation */
552  if (nline > MAXOCLIN) {
553  printf("-E- %s: Number of scene lines: %d, greater than array allocation: %d.\n",
554  "openl1_octs", nline, MAXOCLIN);
555  dataFile.close();
556  exit(FATAL_ERROR);
557  }
558 
559 
560  /* Check that number of pixels per scan not greater than allocation */
561  if (npix > MAXOCPIX) {
562  printf("-E- %s: Number of scan pixels: %d, greater than array allocation: %d.\n",
563  "openl1_octs", npix, MAXOCPIX);
564  dataFile.close();
565  exit(FATAL_ERROR);
566  }
567 
568 
569  lon = (float *) calloc(npix*nline, sizeof (float));
570  lat = (float *) calloc(npix*nline, sizeof (float));
571  senz = (float *) calloc(npix*nline, sizeof (float));
572  sena = (float *) calloc(npix*nline, sizeof (float));
573  solz = (float *) calloc(npix*nline, sizeof (float));
574  sola = (float *) calloc(npix*nline, sizeof (float));
575 
576  l1file->npix = npix;
577  l1file->nscan = nline;
578  strcpy(l1file->spatialResolution, "3.5 km");
579 
580  /* compute and store the view angles for the scene
581  or subscene */
582  navigation_netcdf(dataFile);
583 
584  /* get the time for each line of full scene */
585  // Get a variable
586  NcVar tempVar = dataFile.getVar("msec");
587  tempVar.getVar(msec_temp1);
588 
589  for (i = 0; i < scansPerScene; i++) {
590  for (j = 0; j < linesPerScan; j++) {
591  msec_temp2[i * linesPerScan + j] = msec_temp1[i];
592  }
593  }
594  msec_start = msec_temp2[0];
595  j = 0;
596  for (i = sline; i < (sline + nline); i++) {
597  msec[j] = msec_temp2[i];
598  j++;
599  }
600  // Get the tilt variable
601  tempVar = dataFile.getVar("tilt");
602  tempVar.getVar(tilt);
603 
604  maxday = 365 + abs(LeapCheck(year));
605 
606  /* Get tilt flag, instrument temperature, and gain index */
607  tilt_flag = (int16_t *) calloc(scansPerScene, sizeof (int16_t));
608  tempVar = dataFile.getVar("tilt_flag");
609  tempVar.getVar(tilt_flag);
610  inst_temp = (float *) calloc(4 * scansPerScene, sizeof (float));
611  tempVar = dataFile.getVar("inst_temp");
612  tempVar.getVar(inst_temp);
613  gain = (int16_t *) calloc(8 * scansPerScene, sizeof (int16_t));
614  tempVar = dataFile.getVar("gain");
615  tempVar.getVar(gain);
616  miss_qual = (int16_t *) calloc(scansPerScene, sizeof (int16_t));
617  tempVar = dataFile.getVar("miss_qual");
618  tempVar.getVar(miss_qual);
619 
620  /* Read GAC subsampling table */
621  tempVar = dataFile.getVar("num_tables");
622  tempVar.getVar(&num_tab);
623  start_scan = (int16_t *) calloc(num_tab, sizeof (int16_t));
624  tempVar = dataFile.getVar("start_line");
625  tempVar.getVar(start_scan);
626 
627  for (i = 0; i < num_tab; i++) {
628  if ((extr_line_offset + sline) / linesPerScan < start_scan[i]) {
629  last_table = i - 1;
630  break;
631  }
632  }
633 
634 
635  tempVar = dataFile.getVar("samp_table");
636 
637  vector<size_t> samp_start(tempVar.getDimCount());
638  vector<size_t> samp_edges(tempVar.getDimCount());
639 
640  samp_start[0] = last_table;
641  samp_edges[0] = 1;
642  samp_start[1] = 0;
643  samp_edges[1] = 3;
644  samp_start[2] = 0;
645  samp_edges[2] = 8;
646  samp_start[3] = 0;
647  samp_edges[3] = 2;
648  samp_start[4] = 0;
649  samp_edges[4] = 400;
650  samp_start[5] = 0;
651  samp_edges[5] = 2;
652 
653  tempVar.getVar(samp_start, samp_edges, &sample_table);
654 
655 
656  return (0);
657  }
658  catch (NcException& e) {
659  cout << "-E- Error in reading input NetCDF: " << e.what() << endl;
660  exit(1);
661  }
662 }
663 
664 /* ------------------------------------------------------ */
665 /* readl1_octs() - reads a OCTS HDF level-1 record. */
666 /* */
667 
668 /* ------------------------------------------------------ */
669 extern "C" int32_t readl1_octs_netcdf(filehandle *l1file, int32_t recnum, l1str *l1rec) {
670  int32_t i, status;
671  uint16_t dataarr[NOCBANDS][MAXOCPIX];
672  float scanarr[MAXOCPIX][NOCBANDS];
673  int32_t tilt_deg;
674  int32_t scan;
675  int32_t ip, ib, iw;
676  int32_t nwave = l1rec->l1file->nbands;
677  int32_t *bindx = l1rec->l1file->bindx;
678  char *cal_path = l1_input->calfile;
679  static int32_t FirstCall = 1;
680  int32_t current_scan = recnum / linesPerScan;
681 
682  /* load scan times */
683  /* if the day changed between the start of the scene
684  and the start of the subscene, adjust the day/year
685  accordingly */
686  static int32_t dayIncrimented = 0;
687  if (msec[recnum] < msec_start) {
688  if(!dayIncrimented) {
689  dayIncrimented = 1;
690  day = day + 1;
691  if (day > maxday) {
692  year = year + 1;
693  day = day - maxday;
694  }
695  }
696  }
697  if (msec[recnum] >= 86400000L) {
698  msec[recnum] = msec[recnum] - 86400000L;
699  if(!dayIncrimented) {
700  dayIncrimented = 1;
701  day = day + 1;
702  if (day > maxday) {
703  year = year + 1;
704  day = day - maxday;
705  }
706  }
707  }
708  l1rec->scantime = yds2unix(year, day, (double) (msec[recnum] / 1.e3));
709 
710  /* load standard navigation */
711  memcpy(l1rec->lon, &lon [recnum * npix], sizeof (float)*npix);
712  memcpy(l1rec->lat, &lat [recnum * npix], sizeof (float)*npix);
713  memcpy(l1rec->solz, &solz[recnum * npix], sizeof (float)*npix);
714  memcpy(l1rec->sola, &sola[recnum * npix], sizeof (float)*npix);
715  memcpy(l1rec->senz, &senz[recnum * npix], sizeof (float)*npix);
716  memcpy(l1rec->sena, &sena[recnum * npix], sizeof (float)*npix);
717 
718  /* Read the l1a data */
719  try {
720  NcFile dataFile(l1file->name, NcFile::read);
721 
722  NcVar l1aDataVar = dataFile.getVar("l1a_data");
723 
724  vector<size_t> start(l1aDataVar.getDimCount());
725  vector<size_t> edges(l1aDataVar.getDimCount());
726  start[1] = recnum;
727  start[2] = 0;
728  edges[0] = 1;
729  edges[1] = 1;
730  edges[2] = npix;
731 
732  l1aDataVar.getVar(start, edges, &dataarr);
733 
734 
735  for (i = 0; i < NOCBANDS; i++) {
736  start[0] = i;
737  l1aDataVar.getVar(start, edges, &dataarr[i][0]);
738  }
739 
740  scan = (recnum + extr_line_offset) / linesPerScan;
741 
742  l1rec->tilt = tilt[current_scan];
743 
744  tilt_deg = 0;
745  if (tilt_flag[current_scan] == 1) tilt_deg = +20;
746  if (tilt_flag[current_scan] == 2) tilt_deg = 0;
747  if (tilt_flag[current_scan] == 4) tilt_deg = -20;
748 
749  /* Read the next subsample table if necessary */
750  if (FirstCall == 1 && recnum > 0) {
751  for (i = 0; i < num_tab; i++) {
752  if (scan / linesPerScan < start_scan[i]) {
753  last_table = i;
754  FirstCall = 0;
755  break;
756  }
757  }
758  }
759  NcVar samp_tableVar = dataFile.getVar("samp_table");
760 
761  vector<size_t> samp_start(samp_tableVar.getDimCount());
762  vector<size_t> samp_edges(samp_tableVar.getDimCount());
763 
764  samp_start[0] = last_table;
765  samp_edges[0] = 1;
766  samp_start[1] = 0;
767  samp_edges[1] = 3;
768  samp_start[2] = 0;
769  samp_edges[2] = 8;
770  samp_start[3] = 0;
771  samp_edges[3] = 2;
772  samp_start[4] = 0;
773  samp_edges[4] = 400;
774  samp_start[5] = 0;
775  samp_edges[5] = 2;
776 
777  if (scan >= start_scan[last_table + 1]) {
778  last_table++;
779  samp_start[0] = last_table;
780 
781  samp_tableVar.getVar(samp_start, samp_edges, &sample_table);
782  }
783 
784  }
785  catch (NcException& e) {
786  cout << "-E- Error in reading in readl1a_octs_netcdf: " << e.what() << endl;
787  exit(1);
788  }
789 
790  /* Calibrate the l1a data */
791  status = get_octs_cal(cal_path, year, day, msec, recnum, npix, spix, tilt_deg,
792  gain, inst_temp, sample_table, scansPerScene,
793  dataarr, scanarr);
794 
795  if (status < 0) {
796  fprintf(stderr,
797  "-E- %s line %d: Error applying calibration table \"%s\".\n",
798  __FILE__, __LINE__, cal_path);
799  exit(status);
800  }
801 
802  /* copy the scan array over all bands, and npix pixels */
803  for (ip = 0; ip < npix; ip++) {
804 
805  // if solz NaN set navfail
806  if (isnan(l1rec->solz[ip]))
807  l1rec->navfail[ip] = 1;
808 
809  for (iw = 0; iw < nwave; iw++) {
810  ib = bindx[iw];
811  l1rec->Lt [ip * nwave + ib] = scanarr[ip][iw];
812 
813  /* Also set flags */
814  if (gain[iw * scansPerScene + current_scan] > 0)
815  l1rec->navwarn[ip] = 1;
816  if (miss_qual[current_scan] > 0)
817  l1rec->navwarn[ip] = 1;
818  if (dataarr[iw][ip] > 1022)
819  l1rec->hilt[ip] = 1;
820  if (scanarr[ip][iw] <= 0.0)
821  l1rec->navfail[ip] = 1;
822  }
823  }
824 
825  l1rec->npix = l1file->npix;
826 
827  return (0);
828 }
829 
830 
831 /* ------------------------------------------------------ */
832 /* closel1_octs() - closes the level 1 OCTS HDF file */
833 /* */
834 
835 /* ------------------------------------------------------ */
836 extern "C" int32_t closel1_octs_netcdf(filehandle *l1file) {
837  free(lon);
838  free(lat);
839  free(senz);
840  free(sena);
841  free(solz);
842  free(sola);
843 
844  if (l1file->format == FT_OCTSL1A) {
845  free(start_scan);
846  free(tilt_flag);
847  free(inst_temp);
848  }
849 
850  /* End access to the HDF file */
851  SDend(l1file->sd_id);
852 
853  return (0);
854 }
l1file
int32 l1file(int32 sdfid, int32 *nsamp, int32 *nscans, int16 *dtynum)
Definition: l1stat_chk.c:586
bindx
const int bindx[3]
Definition: DbLutNetcdf.cpp:28
j
int j
Definition: decode_rs.h:73
day
int32_t day
Definition: atrem_corl1v2.h:308
status
int status
Definition: l1_czcs_hdf.c:32
yds2unix
double yds2unix(int16_t year, int16_t day, double secs)
Definition: yds2unix.c:7
L
#define L(lambda, T)
Definition: PreprocessP.h:185
gain
int16 * gain
Definition: l1_czcs_hdf.c:33
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
l1rec
read l1rec
Definition: HOWTO_Add_a_sensor.txt:72
spline
subroutine spline(s, x, y, n, in, t, il, iu, vl, vu, e, u)
Definition: phs.f:1348
pos
float32 * pos
Definition: l1_czcs_hdf.c:35
filehdr_struc.h
msec
int32 * msec
Definition: l1_czcs_hdf.c:31
FT_OCTSL1A
@ FT_OCTSL1A
Definition: filetype.h:36
navigation_netcdf
void navigation_netcdf(NcFile &dataFile)
Definition: l1_octs_netcdf.cpp:175
readl1_octs_netcdf
int32_t readl1_octs_netcdf(filehandle *l1file, int32_t recnum, l1str *l1rec)
Definition: l1_octs_netcdf.cpp:669
PI
#define PI
Definition: l3_get_org.c:6
NOCBANDS
#define NOCBANDS
Definition: l1_octs_netcdf.cpp:28
recnum
read recnum
Definition: HOWTO_Add_a_sensor.txt:72
zulu2unix
double zulu2unix(char *zulu)
Definition: zulu2unix.c:3
MAXOCPIX
#define MAXOCPIX
Definition: l1_octs_netcdf.cpp:26
FL
#define FL
Definition: l1_octs_netcdf.cpp:21
RADIUS
#define RADIUS
Definition: l1_octs_netcdf.cpp:20
LeapCheck
int LeapCheck(int yr)
Definition: l1_octs.c:1126
l1_input
l1_input_t * l1_input
Definition: l1_options.c:9
FATAL_ERROR
#define FATAL_ERROR
Definition: swl0_parms.h:5
want_verbose
int want_verbose
Definition: genutils_globals.c:3
MAXOCLIN
#define MAXOCLIN
Definition: l1_octs_netcdf.cpp:25
isodate2ydmsec
void isodate2ydmsec(char *date, int32_t *year, int32_t *day, int32_t *msec)
Definition: date2ydmsec.c:20
RADEG
#define RADEG
Definition: czcs_ctl_pt.c:5
seadasutils.aquarius_utils.start
start
Definition: aquarius_utils.py:27
color_dtdb.edges
edges
Definition: color_dtdb.py:331
spix
int32 spix
Definition: l1_czcs_hdf.c:21
openl1_octs_netcdf
int32_t openl1_octs_netcdf(filehandle *l1file)
Definition: l1_octs_netcdf.cpp:503
splint
subroutine splint(xa, ya, y2a, n, x, y)
Definition: atrem_app_refl_f90_cubeio.f:3270
R
#define R
Definition: make_L3_v1.1.c:96
get_octs_cal
int32_t get_octs_cal(char *file, int16_t year, int16_t day, int32_t msec[MAXOCLIN], int32_t recnum, int16_t npix, int32_t spix, int32_t tilt, int16_t gainset[MAXOCLIN], float inst_temp[MAXOCLIN], int16_t sample_table[3][8][2][400][2], int32_t scansPerScene, uint16_t l1acnts[NOCBANDS][MAXOCPIX], float l1brads[MAXOCPIX][NOCBANDS])
Definition: l1_octs.c:95
MAXOCARR
#define MAXOCARR
Definition: l1_octs_netcdf.cpp:27
tilt
int16_t * tilt
Definition: l2bin.cpp:80
abs
#define abs(a)
Definition: misc.h:90
i
int i
Definition: decode_rs.h:71
l1_octs_netcdf.h
strcpy
How many dimensions is the output array Default is Not sure if anything above will work correctly strcpy(l2prod->title, "no title yet")
l1_octs.h
closel1_octs_netcdf
int32_t closel1_octs_netcdf(filehandle *l1file)
Definition: l1_octs_netcdf.cpp:836
k
int k
Definition: decode_rs.h:73
npix
int npix
Definition: get_cmp.c:28
reform_octs_time
void reform_octs_time(char *time)
Definition: l1_octs.c:1167


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