ocssw
V2022
|
make_L3_v1.1.c
Go to the documentation of this file.
6 cc -Ae +O2 +Ofastaccess +Onolimit make_L3.c -o make_L3 ./LIB/lib.a -lf -lm -I /usr/local/HDF4.0/include -L/usr/local/HDF4.0/lib -lmfhdf -ldf -ljpeg -lz
8 cc -O2 -fullwarn -woff 1209,1210 -64 make_L3.c -o make_L3 ./LIB/lib.a -I./LIB -I$HDFINC -L$HDFLIB -lmfhdf -ldf -lz -ljpeg -lm -lmalloc
9 or cc -O2 -fullwarn -woff 1209,1210 -64 make_L3.c -o make_L3 ./LIB/lib.a -I./LIB -I$HDFINC -L/home/jack/libhdf64bit.4.1r3/ -lmfhdf -ldf -lz -lm -ljpeg
79 #define WAVELENGTH { "412", "443", "490", "510", "555", "670", "765", "865"} /* Band wavelength */
80 #define IRRAD {171.07, 189.91, 194.33, 188.24, 185.93, 152.14, 123.55, 100.00} /* TOA solar irradiance */
81 #define TAURAY {0.3139, 0.2341, 0.1561, 0.1324, 0.0942, 0.0439, 0.0257, 0.0155} /* TOA molecular reflectance */
82 #define OZONE {.00103, .00400, .02536, .04200, .09338, .04685, .00837, .00485} /* Coef. for ozone absorption */
83 #define WVAPOR {0.0000, 0.0000, 0.0000, 0.0000, 0.0008, 0.0043, 0.0007, 0.0057} /* Coef. for water vapor abs. */
84 #define OXYGEN {0.0000, 0.0000, 0.0000, 0.0000, 0.0000, 0.0000, 0.0000, 0.0000} /* Coef. for oxygen absorption */
103 #define CROSS_PRODUCT(x1, y1, x2, y2, x3, y3) ( ((x2) - (x1)) * ((y3) - (y1)) - ((y2) - (y1)) * ((x3) - (x1)) )
151 int get_lonlat(double xproj, double yproj, double *lon, double *lat, int projtype, double lon_center);
152 int get_xyproj(double lon, double lat, double *xproj, double *yproj, int projtype, double lon_center);
157 void process_pixel(int ipix, int iline_out, int ipix_out, char **wl, INPUT input, OUTPUT *proj, PROCESSINFO process);
159 int read_file_block(char *filename, char **wl, int offset, int startline, int buflines, OUTPUT *proj);
161 void resample_geo_scan(int Npixels, float *inlon, float *inlat, double *outlon, double *outlat);
164 int update_file(char *filename, char **wl, int offset, int startline, int buflines, OUTPUT *proj, char verbose, char gzip);
325 if (get_attributes(input.sd_id, prefix, &input.Nscans, &input.Npixels, &input.startpix, &input.subsampling, &input.year, &input.doy, &input.orbno, input.title) == -1) {
392 if (update_file(proj.filename, wl, 0, proj.startline, proj.buflines, &proj, TRUE, process.gzip) == -1)
448 printf("Calibration -> Geolocation -> Computation of surface reflectance (atmospheric correction) -> Projection ...\n");
450 printf("Calibration -> Geolocation -> Computation of top-of-the-atmosphere reflectance (no atmospheric correction) -> Projection ...\n");
459 proj.filelist = (char *) realloc(proj.filelist, (strlen(proj.filelist) + strlen(inputfile) + 8) * sizeof (char));
461 proj.filelist = (char *) malloc((strlen(inputfile) + 6) * sizeof (char)); /* +1 for initial "\n", +3 for "0: ", */
552 if ((input.level == L1A && (input.corruption = compute_l1b(input.sd_id, iscan, input.rho, CALIBTABLE, input.newfile)) == -1) ||
685 if (update_file(proj.filename, wl, 0, proj.startline, proj.buflines, &proj, TRUE, process.gzip) == -1)
765 scan = (double) ((input->startpix - 1 + input->subsampling * ipix - (NpixelsLAC - 1) / 2.) * IFOV);
780 a = coef[0] * cosa[ipix] * cosa[ipix] + coef[1] * cosa[ipix] * sina[ipix] + coef[2] * sina[ipix] * sina[ipix];
874 void resample_geo_scan(int Npixels, float *inlon, float *inlat, double *outlon, double *outlat) {
902 int read_file_block(char *filename, char **wl, int offset, int startline, int buflines, OUTPUT *proj) {
975 int update_file(char *filename, char **wl, int offset, int startline, int buflines, OUTPUT *proj, char verbose, char gzip) {
991 write_sds(sd_id, SDSname, proj->refl[ib] + offset, startline, buflines, proj->h, proj->w, DFNT_INT16, &fill_int16, sf_refl, verbose, gzip);
997 write_sds(sd_id, SDSname, proj->reflTOA[ib] + offset, startline, buflines, proj->h, proj->w, DFNT_INT16, &fill_int16, sf_refl, verbose, gzip);
1001 write_sds(sd_id, SDSname, proj->ifno + offset, startline, buflines, proj->h, proj->w, DFNT_INT16, &fill_int16, 0, verbose, gzip);
1004 write_sds(sd_id, SDSname, proj->senz + offset, startline, buflines, proj->h, proj->w, DFNT_INT16, &fill_int16, sf_angle, verbose, gzip);
1008 write_sds(sd_id, SDSname, proj->solz + offset, startline, buflines, proj->h, proj->w, DFNT_INT16, &fill_int16, sf_angle, verbose, gzip);
1010 write_sds(sd_id, SDSname, proj->phi + offset, startline, buflines, proj->h, proj->w, DFNT_INT16, &fill_int16, sf_angle, verbose, gzip);
1014 write_sds(sd_id, SDSname, proj->smoke + offset, startline, buflines, proj->h, proj->w, DFNT_INT16, &fill_int16, sf_smoke, verbose, gzip);
1018 write_sds(sd_id, SDSname, proj->ndvi + offset, startline, buflines, proj->h, proj->w, DFNT_INT16, &fill_int16, sf_vi, verbose, gzip);
1022 write_sds(sd_id, SDSname, proj->evi + offset, startline, buflines, proj->h, proj->w, DFNT_INT16, &fill_int16, sf_vi, verbose, gzip);
1026 write_sds(sd_id, SDSname, proj->gemi + offset, startline, buflines, proj->h, proj->w, DFNT_INT16, &fill_int16, sf_vi, verbose, gzip);
1029 if (SDsetattr(sd_id, "File list", DFNT_CHAR8, (int) strlen(proj->filelist) + 1, proj->filelist) == -1)
1044 memcpy(proj->refl[ib] + (offset2 + irow) * Ncols, proj->refl[ib] + (offset1 + irow) * Ncols, Ncols * sizeof (int16));
1046 memcpy(proj->reflTOA[ib] + (offset2 + irow) * Ncols, proj->reflTOA[ib] + (offset1 + irow) * Ncols, Ncols * sizeof (int16));
1048 memcpy(proj->ifno + (offset2 + irow) * Ncols, proj->ifno + (offset1 + irow) * Ncols, Ncols * sizeof (int16));
1050 memcpy(proj->senz + (offset2 + irow) * Ncols, proj->senz + (offset1 + irow) * Ncols, Ncols * sizeof (int16));
1052 memcpy(proj->solz + (offset2 + irow) * Ncols, proj->solz + (offset1 + irow) * Ncols, Ncols * sizeof (int16));
1053 memcpy(proj->phi + (offset2 + irow) * Ncols, proj->phi + (offset1 + irow) * Ncols, Ncols * sizeof (int16));
1056 memcpy(proj->smoke + (offset2 + irow) * Ncols, proj->smoke + (offset1 + irow) * Ncols, Ncols * sizeof (int16));
1058 memcpy(proj->ndvi + (offset2 + irow) * Ncols, proj->ndvi + (offset1 + irow) * Ncols, Ncols * sizeof (int16));
1060 memcpy(proj->evi + (offset2 + irow) * Ncols, proj->evi + (offset1 + irow) * Ncols, Ncols * sizeof (int16));
1062 memcpy(proj->gemi + (offset2 + irow) * Ncols, proj->gemi + (offset1 + irow) * Ncols, Ncols * sizeof (int16));
1066 void process_pixel(int ipix, int iline_out, int ipix_out, char **wl, INPUT input, OUTPUT *proj, PROCESSINFO process) {
1091 if (update_file(proj->filename, wl, 0, old_startline, offset, proj, FALSE, process.gzip) == -1)
1094 if (read_file_block(proj->filename, wl, (proj->buflines - offset) * proj->w, proj->startline + proj->buflines - offset,
1100 if (update_file(proj->filename, wl, (proj->buflines - offset) * proj->w, old_startline + proj->buflines - offset,
1108 if (update_file(proj->filename, wl, 0, old_startline, proj->buflines, proj, FALSE, process.gzip) == -1)
1171 int get_xyproj(double lon, double lat, double *xproj, double *yproj, int projtype, double lon_center) {
1178 case PLATECARREE: *xproj = lon * RAD2DEG; /* PLATECARREE was not initialized with RAD2DEG factor */
1181 case SINUSOIDAL: *xproj = lon * cos(lat) * RAD2DEG; /* SINUSOIDAL was not initialized with RAD2DEG factor */
1195 int get_lonlat(double xproj, double yproj, double *lon, double *lat, int projtype, double lon_center) {
1269 proj->evi[idx_out] = EVI(input.rho[BAND2][ipix], input.rho[BAND6][ipix], input.rho[BAND8][ipix]);
1287 return (int16) (rho[BAND1][ipix] * (rho[BAND2][ipix] - rho[BAND4][ipix]) * (rho[BAND5][ipix] - rho[BAND6][ipix]) * 10000);
1375 if ((*data = (void *) malloc(size * (sizeoftype = DFKNTsize(num_type)))) == NULL) { /* Allocate a new array */
1379 printf(" (%.1f Mbytes allocated: %d items of %d bytes)\n", size * sizeoftype / 1e6, size, sizeoftype);
1395 new_cmd_line = (char *) realloc(new_cmd_line, (count + strlen(cmd_line) + 1 + 1) * sizeof (char));
1491 proj->cmd_line = (char *) realloc(proj->cmd_line, (strlen(proj->cmd_line) + strlen(arg) + 1 + 1) * sizeof (char));
1730 (int) ((1e6 * process->maxmem - NLPSURF * NPPSURF * DFKNTsize(DFNT_INT16)) / proj->w / nbytes));
1737 printf("Output projection has to be bufferized: %d lines out of %d\n", proj->buflines, proj->h);
1739 /* bufstep is the shift (in # of lines) to apply to the buffer when pixel is out of the current buffer */
1778 case MINBLUE: sprintf(process->compositename, "minimum band %d (%s nm)", BLUEBAND + 1, wl[BLUEBAND]);
1799 fprintf(stderr, "Cannot make composite without band %d (%s nm)\n", BLUEBAND + 1, wl[BLUEBAND]);
1819 fprintf(stderr, "Cannot compute EVI without bands %s, %s and %s\n", wl[BAND2], wl[BAND6], wl[BAND8]);
1854 alloc_new_array((void *) &proj->refl[ib], proj->buflines * proj->w, DFNT_INT16, process.maxmem);
1864 alloc_new_array((void *) &proj->reflTOA[ib], proj->buflines * proj->w, DFNT_INT16, process.maxmem);
1934 if (SDsetattr(proj->sd_id, "Projection", DFNT_CHAR8, (int) strlen(proj->projname) + 1, proj->projname) == -1)
1941 if (SDsetattr(proj->sd_id, "Maximum sensor zenith angle (degrees)", DFNT_FLOAT32, 1, &process.maxsenz) == -1)
1942 fprintf(stderr, "Cannot write attribute \"%s\" in file %s\n", "Maximum sensor zenith angle", proj->filename);
1946 fprintf(stderr, "Cannot write attribute \"%s\" in file %s\n", "Compositing criterion", proj->filename);
1950 fprintf(stderr, "Cannot write attribute \"%s\" in file %s\n", "Subset Scan Range", proj->filename);
1952 fprintf(stderr, "Cannot write attribute \"%s\" in file %s\n", "Subset Pixel Range", proj->filename);
2038 if (fabs(input->tilt[iscan] - input->tilt[iscan + 1 - 2 * (iscan == input->Nscans - 1)]) > 0.1)
2254 fprintf(stderr, "SDS %s has the wrong dimensions (%dx%d)\n", L2sds[k].sds_name, input->Npixels, input->Nscans);
2274 ((float32 *) L2sds[k].data)[j] = L2sds[k].slope * ((float32 *) L2sds[k].data)[j] + L2sds[k].intercept;
2275 else if (k > 2 * NBANDS && k < L2FIELDS - 6) /* TOA refl, longitude and latitude: float32 -> float32 */
2277 ((float32 *) L2sds[k].data)[j] = L2sds[k].slope * ((float32 *) L2sds[k].data)[j] + L2sds[k].intercept;
2283 input->ndvi[j] = floor((L2sds[L2FIELDS - 6].slope * input->ndvi[j] + L2sds[L2FIELDS - 6].intercept) / sf_vi + 0.5);
2286 input->evi[j] = floor((L2sds[L2FIELDS - 5].slope * input->evi[j] + L2sds[L2FIELDS - 5].intercept) / sf_vi + 0.5);
2289 input->solz[j] = (L2sds[L2FIELDS - 4].slope * input->solz[j] + L2sds[L2FIELDS - 4].intercept) / sf_angle + 0.5;
void check_process_parameters(PROCESSINFO *process, char **wl)
Definition: make_L3_v1.1.c:1762
void red(int iz1, int iz2, int jz1, int jz2, int jm1, int jm2, int jmf, int ic1, int jc1, int jcf, int kc, float ***c, float **s)
Definition: make_L3_v1.1.c:124
void make_psurf(signed short *psurf, char *file, int Nlines, int Npixels, float psurf0)
Definition: make_psurf.c:6
void parse_command_line(int argc, char **argv, OUTPUT *proj, PROCESSINFO *process)
Definition: make_L3_v1.1.c:1480
int read_sds_block(int32 sds_id, int iline, int buflines, void *dataP)
Definition: read_write.c:19
These are used to scale the SD before writing it to the HDF4 file The default is and which means the product is not scaled at all Since the product is usually stored as a float inside of this is a way to write the float out as a integer l2prod min
Definition: HOWTO_Add_a_product.txt:76
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
void resample_geo_scan(int Npixels, float *inlon, float *inlat, double *outlon, double *outlat)
Definition: make_L3_v1.1.c:874
int molwforint(double r, double center_long, double false_east, double false_north)
Definition: proj_molwfor.c:37
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 and name_suffix ie aph_412_giop where prod_ix will be set to PARAM_TYPE_IR_WAVE same search method as PARAM_TYPE_VIS_WAVE except only wavelength above are looped through but prod_ix is still based ie aph_2_giop for the second band
Definition: HOWTO_Add_a_product.txt:42
int robinvint(double r, double center_long, double false_east, double false_north)
Definition: proj_robinv.c:46
HDF4 data type of the output SDS Default is DFNT_FLOAT32 Common types used DFNT_INT32
Definition: HOWTO_Add_a_product.txt:67
int write_sds(int32 sd_id, char *SDSname, void *data, int startline, int buflines, int Nl, int Np, int32 numtype, void *fillvalue, float64 scale_factor, char verbose, char gzip)
Definition: read_write.c:75
int update_file(char *filename, char **wl, int offset, int startline, int buflines, OUTPUT *proj, char verbose, char gzip)
Definition: make_L3_v1.1.c:975
int get_attributes(int32 sd_id, char *prefix, int *Nlines, int *Npixels, int *startpix, int *subsampling, int *year, int *doy, int *orbno, char *title)
Definition: get_attributes.c:4
int alloc_new_array(void **data, int size, int32 num_type, int maxmem)
Definition: make_L3_v1.1.c:1372
void set_buffer_parameters(OUTPUT *proj, PROCESSINFO *process)
Definition: make_L3_v1.1.c:1695
void init_output_file(OUTPUT *proj, INPUT input, PROCESSINFO process)
Definition: make_L3_v1.1.c:1917
int setlinebuf(FILE *stream)
Definition: make_L3_v1.1.c:107
int robforint(double r, double center_long, double false_east, double false_north)
Definition: proj_robfor.c:43
void update_pixel(INPUT input, int16 ndvi, int ipix, int idx_out, int ifno, OUTPUT *proj)
Definition: make_L3_v1.1.c:1223
char * strdup(const char *)
void process_pixel(int ipix, int iline_out, int ipix_out, char **wl, INPUT input, OUTPUT *proj, PROCESSINFO process)
Definition: make_L3_v1.1.c:1066
void memcpy_block(int offset2, int offset1, int Nrows, int Ncols, OUTPUT *proj)
Definition: make_L3_v1.1.c:1036
int get_xyproj(double lon, double lat, double *xproj, double *yproj, int projtype, double lon_center)
Definition: make_L3_v1.1.c:1171
HDF4 data type of the output SDS Default is DFNT_FLOAT32 Common types used DFNT_INT16
Definition: HOWTO_Add_a_product.txt:67
no change in intended resolving MODur00064 Corrected handling of bad ephemeris attitude data
Definition: HISTORY.txt:356
int get_geolocation_scan(int iscan, INPUT *input, PROCESSINFO process)
Definition: make_L3_v1.1.c:748
int molwinvint(double r, double center_long, double false_east, double false_north)
Definition: proj_molwinv.c:35
Extra metadata that will be written to the HDF4 file l2prod rank
Definition: HOWTO_Add_a_product.txt:80
int get_lonlat(double xproj, double yproj, double *lon, double *lat, int projtype, double lon_center)
Definition: make_L3_v1.1.c:1195
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 and name_suffix ie aph_412_giop where prod_ix will be set to PARAM_TYPE_IR_WAVE same search method as PARAM_TYPE_VIS_WAVE except only wavelength above are looped through but prod_ix is still based ie aph_2_giop for the second and prod_ix set to PARAM_TYPE_INT name_prefix is compared with the beginning of the product name If name_suffix is not empty the it must match the end of the product name The characters right after the prefix are read as an integer and prod_ix is set to that number strncpy(l2prod->name_prefix, "myprod", UNITLEN)
HDF4 data type of the output SDS Default is DFNT_FLOAT32 Common types used DFNT_FLOAT32
Definition: HOWTO_Add_a_product.txt:67
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 as required for compatibility with version of the SDP toolkit Corrected test output file names to end in per delivery and then split off a new MYD_PR03 pcf file for Aqua Added AssociatedPlatformInstrumentSensor to the inventory metadata in MOD01 mcf and MOD03 mcf Created new versions named MYD01 mcf and MYD03 where AssociatedPlatformShortName is rather than Terra The program itself has been changed to read the Satellite Instrument validate it against the input L1A and LUT and to use it determine the correct files to retrieve the ephemeris and attitude data from Changed to produce a LocalGranuleID starting with MYD03 if run on Aqua data Added the Scan Type file attribute to the Geolocation copied from the L1A and attitude_angels to radians rather than degrees The accumulation of Cumulated gflags was moved from GEO_validate_earth_location c to GEO_locate_one_scan c
Definition: HISTORY.txt:464
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] angles
Definition: HISTORY.txt:366
int compute_l1b(int32 sd_id, int iline, float *l1b_data[], char *calibtable, char newfile)
Definition: compute_l1b.c:9
Definition: hybrid.c:27
How many dimensions is the output array Default is Not sure if anything above will work correctly strcpy(l2prod->title, "no title yet")
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
int haminvint(double r, double center_long, double false_east, double false_north)
Definition: proj_haminv.c:40
int32 read_sds(l1info_struct l1info, char *arr_name, int32 *exp_ntyp, void *array)
Definition: read_sds.c:5
int project_scan(int iscan, char **wl, INPUT *input, OUTPUT *proj, PROCESSINFO process)
Definition: make_L3_v1.1.c:1962
Definition: make_L3_v1.1.c:134
int read_file_block(char *filename, char **wl, int offset, int startline, int buflines, OUTPUT *proj)
Definition: make_L3_v1.1.c:902
int hamforint(double r, double center_long, double false_east, double false_north)
Definition: proj_hamfor.c:41
void init_output_projection(OUTPUT *proj, PROCESSINFO process)
Definition: make_L3_v1.1.c:1848