|
ocssw
1.0
|
00001 /* ============================================================================ */ 00002 /* module l1_hmodis_hdf.c - functions to read MODIS HIRES L1B for MSL12 */ 00003 /* */ 00004 /* Written By: B. Franz, NASA/SIMBIOS, January 2003. */ 00005 /* Modified By: J. Gales, NASA/OBPG, August 2005. */ 00006 /* Conversion for 16-band MODIS: B. Franz, J. Gales, 2006. */ 00007 /* */ 00008 /* ============================================================================ */ 00009 00010 #include "l1_hmodis_hdf.h" 00011 #include "l1_modis_hdf.h" 00012 #include "l12_proto.h" 00013 #include "hdf_utils.h" 00014 #include "alloc_2d.h" 00015 00016 #define EOSMETALEN 32768 00017 00018 #define MAXBANDS (NBANDS+NBANDSIR+1) 00019 #define NDET 40 00020 00021 static int have_extract = 0; 00022 00023 static int32 sd_id[3] = {-1}; 00024 static int32 sd_id_g; 00025 static int32 spix = 0; 00026 static int32 sscan = 0; /* start scan in original granule */ 00027 static int32 rscan = 0; /* start scan in extract granule */ 00028 static int32 resolution = 1000; 00029 static int32 ndet = 10; 00030 00031 static int32 np1000 = 1354; 00032 static int32 sp1000 = 0; 00033 static int32 np500 = 2708; 00034 static int32 sp500 = 0; 00035 static int32 np250 = 5416; 00036 static int32 sp250 = 0; 00037 00038 static int32 nd1000 = 10; 00039 static int32 nd500 = 20; 00040 static int32 nd250 = 40; 00041 00042 static float32 **lon1000 = NULL; 00043 static float32 **lat1000 = NULL; 00044 static float32 **hgt1000 = NULL; 00045 static float32 **solz1000 = NULL; 00046 static float32 **sola1000 = NULL; 00047 static float32 **senz1000 = NULL; 00048 static float32 **sena1000 = NULL; 00049 static float32 **fdatabuf[MAXBANDS]; 00050 00051 static int destripe = 0; 00052 00053 /* RSMAS corrections to Terra IR radiances (20,22,23,27,28,29,31,32) */ 00054 static float radoff[NBANDSIR][NDET] = { 00055 {-0.000972, 0.014200, 0.000022, 0.000238,-0.000003, 0.003760, 0.002337, 0.000084, 0.008848, 0.005050}, 00056 {-0.002300, 0.000600, 0.000236,-0.000280,-0.000003, 0.002798, 0.003496, 0.018035, 0.002942, 0.001787}, 00057 {-0.003833, 0.003657, 0.002118, 0.000798,-0.000003, 0.000208, 0.000399, 0.000553, 0.000258, 0.021128}, 00058 {0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0}, 00059 {0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0}, 00060 {0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0}, 00061 {-0.000423,-0.000242,-0.000330,-0.000065,-0.000001,-0.000006, 0.000064, 0.000398, 0.000362, 0.000322}, 00062 {-0.000433,-0.000246,-0.000222,-0.000148,-0.000001, 0.000068, 0.000221, 0.000221, 0.000437, 0.000303} 00063 }; 00064 static float mfact[NBANDSIR] = {0.0,0.0,0.0005,0.0,0.0,0.0,0.0028,0.00278}; 00065 00066 /* RSMAS AQUA modis corrections for some bands and detectors */ 00067 static float ch20d2cor = 0.996666; /* channel 20 detector 2 */ 00068 static float ch20d3cor = 1.0044444; /* channel 20 detector 3 */ 00069 static float ch20d10cor = 1.0033333; /* channel 20 detector 10 */ 00070 static float ch22d1cor = 1.0030075; /* channel 22 detector 1 */ 00071 static float ch22d8cor = 1.0007518; /* channel 22 detector 8 */ 00072 static float ch22d9cor = 1.0015037; /* channel 22 detector 9 */ 00073 static float ch22d10cor = 1.0022556; /* channel 22 detector 10 */ 00074 static float ch23d1cor = 1.0192857; /* channel 23 detector 1 */ 00075 00076 00077 /* ----------------------------------------------------------------------------------- */ 00078 /* modis_data_interp() - interpolates from 1km counts to 250 or 500m pixels. */ 00079 /* */ 00080 /* ipixl = full resolution pixel (frame) number, from zero */ 00081 /* iline = full resolution line (detector) number, from zero */ 00082 /* */ 00083 /* B. Franz, SAIC, January 2006 . */ 00084 /* ----------------------------------------------------------------------------------- */ 00085 void modis_data_interp(int band, int32_t ipixl2, int32_t iline2, int res1, int res2, float32 *data) 00086 { 00087 int32_t ipixl; 00088 int32_t iline; 00089 float fpixl; 00090 float fline; 00091 float dpixl; 00092 float dline; 00093 int nline; 00094 int npixl; 00095 float val[4]; 00096 float sumv; 00097 int cntv; 00098 uint16 maxv; 00099 int i; 00100 00101 /* note: this function presumes that ipixl2=0 is at the center of a 1km pixel, */ 00102 /* as it will always be for files generated via MODl1extract L1A */ 00103 00104 switch (res1+res2) { 00105 case 1500: /* 1000 to 500 */ 00106 nline = nd1000; 00107 npixl = np1000; 00108 fpixl = 0.50*ipixl2; 00109 fline = 0.50*iline2 - 0.25; 00110 break; 00111 case 1250: /* 1000 to 250 */ 00112 nline = nd1000; 00113 npixl = np1000; 00114 fpixl = 0.25*ipixl2; 00115 fline = 0.25*iline2 - 0.375; 00116 break; 00117 case 750: /* 500 to 250 */ 00118 nline = nd500; 00119 npixl = np500; 00120 fpixl = 0.50*ipixl2; 00121 fline = 0.50*iline2 - 0.25; 00122 break; 00123 default: 00124 printf("%s Line %d: Unknown resolution change, %d to %d\n",__FILE__,__LINE__,res1,res2); 00125 exit(1); 00126 } 00127 00128 /* index and weighting into input data block */ 00129 iline = MIN(MAX((int)floor((float)fline),0),nline-2); 00130 dline = (float)fline - (float)iline; 00131 ipixl = MIN(MAX((int)floor((float)fpixl),0),npixl-2); 00132 dpixl = (float)fpixl - (float)ipixl; 00133 00134 /* interpolation: account for flagged values */ 00135 00136 val[0] = fdatabuf[band][iline+0][ipixl+0]; 00137 val[1] = fdatabuf[band][iline+0][ipixl+1]; 00138 val[2] = fdatabuf[band][iline+1][ipixl+0]; 00139 val[3] = fdatabuf[band][iline+1][ipixl+1]; 00140 00141 sumv = 0.0; 00142 cntv = 0; 00143 maxv = 0; 00144 00145 for( i = 0; i < 4; i++ ) { 00146 if (val[i] <= 65500) { 00147 sumv = sumv + val[i]; 00148 cntv++; 00149 } 00150 if( val[i] > maxv ) 00151 maxv = val[i]; 00152 } 00153 00154 /* If we have four non-flagged samples, we bilinearly interpolate. If we */ 00155 /* have 0-3 samples, just set to the max flag value (or use average) */ 00156 00157 switch (cntv) { 00158 case 4: 00159 *data = ( (1.0-dline)*(val[0]*(1.0-dpixl) + val[1]*dpixl) + 00160 dline*(val[2]*(1.0-dpixl) + val[3]*dpixl) ); 00161 break; 00162 case 0: 00163 *data = maxv; 00164 break; 00165 default: 00166 //*data = maxv; 00167 *data = sumv/cntv; // re-enabled this; don't know why it was disabled 00168 break; 00169 } 00170 } 00171 00172 00173 /* ----------------------------------------------------------------------------------- */ 00174 /* modis_geo_interp() - interpolates from 1km geolocation to 250 or 500m pixels. */ 00175 /* */ 00176 /* ipixl = full resolution pixel (frame) number, from zero */ 00177 /* iline = full resolution line (detector) number, from zero */ 00178 /* */ 00179 /* B. Franz, SAIC, January 2006 . */ 00180 /* ----------------------------------------------------------------------------------- */ 00181 void modis_geo_interp(int32_t ipixl2, int32_t iline2, int resolution, 00182 float *lon, float *lat, float *hgt, 00183 float *solz, float *sola, 00184 float *senz, float *sena, int lonlat) 00185 { 00186 float fpixl; 00187 float fline; 00188 double dpixl; 00189 double dline; 00190 int32_t ipixl; 00191 int32_t iline; 00192 double min_val; 00193 double max_val; 00194 double val[4]; 00195 int ival; 00196 00197 /* note: this function presumes that ipixl2=0 is at the center of a 1km pixel, */ 00198 /* as it will always be for files generated via MODl1extract L1A */ 00199 00200 switch (resolution) { 00201 case 250: 00202 fpixl = 0.25*ipixl2; 00203 fline = 0.25*iline2 - 0.375; 00204 break; 00205 case 500: 00206 fpixl = 0.50*ipixl2; 00207 fline = 0.50*iline2 - 0.25; 00208 break; 00209 default: 00210 printf("%s Line %d: Unknown resolution = %d\n",__FILE__,__LINE__,resolution); 00211 exit(1); 00212 } 00213 00214 /* index and weighting into 1km geolocation block */ 00215 iline = MIN(MAX((int)floor((double)fline),0),nd1000-2); 00216 dline = (double)fline - (double)iline; 00217 ipixl = MIN(MAX((int)floor((double)fpixl),0),np1000-2); 00218 dpixl = (double)fpixl - (double)ipixl; 00219 00220 /* longitude interpolation: account for dateline */ 00221 00222 val[0] = (double)lon1000[iline+0][ipixl+0]; 00223 val[1] = (double)lon1000[iline+0][ipixl+1]; 00224 val[2] = (double)lon1000[iline+1][ipixl+0]; 00225 val[3] = (double)lon1000[iline+1][ipixl+1]; 00226 00227 min_val = max_val = val[0]; 00228 00229 for( ival = 1; ival < 4; ival++ ) { 00230 if( val[ival] < min_val ) 00231 min_val = val[ival]; 00232 else if( val[ival] > max_val ) 00233 max_val = val[ival]; 00234 } 00235 00236 if( max_val > 180.0 + min_val ) { 00237 for( ival = 0; ival < 4; ival++ ) { 00238 if( val[ival] < 0.0 ) 00239 val[ival] += 360.0; 00240 } 00241 } 00242 00243 *lon = ((1.0-dline)*(val[0]*(1.0-dpixl) + val[1]*dpixl) + 00244 dline*(val[2]*(1.0-dpixl) + val[3]*dpixl)); 00245 00246 if( *lon >= 180.0 ) 00247 *lon = *lon-360.0; 00248 00249 00250 /* latitude interpolation */ 00251 00252 val[0] = (double)lat1000[iline+0][ipixl+0]; 00253 val[1] = (double)lat1000[iline+0][ipixl+1]; 00254 val[2] = (double)lat1000[iline+1][ipixl+0]; 00255 val[3] = (double)lat1000[iline+1][ipixl+1]; 00256 00257 *lat = ((1.0-dline)*(val[0]*(1.0-dpixl) + val[1]*dpixl) + 00258 dline*(val[2]*(1.0-dpixl) + val[3]*dpixl)); 00259 00260 if(lonlat) 00261 return; 00262 00263 /* height interpolation */ 00264 00265 val[0] = (double)hgt1000[iline+0][ipixl+0]; 00266 val[1] = (double)hgt1000[iline+0][ipixl+1]; 00267 val[2] = (double)hgt1000[iline+1][ipixl+0]; 00268 val[3] = (double)hgt1000[iline+1][ipixl+1]; 00269 00270 *hgt = ((1.0-dline)*(val[0]*(1.0-dpixl) + val[1]*dpixl) + 00271 dline*(val[2]*(1.0-dpixl) + val[3]*dpixl)); 00272 00273 00274 /* solar zenith interpolation */ 00275 00276 val[0] = (double)solz1000[iline+0][ipixl+0]; 00277 val[1] = (double)solz1000[iline+0][ipixl+1]; 00278 val[2] = (double)solz1000[iline+1][ipixl+0]; 00279 val[3] = (double)solz1000[iline+1][ipixl+1]; 00280 00281 *solz = ((1.0-dline)*(val[0]*(1.0-dpixl) + val[1]*dpixl) + 00282 dline*(val[2]*(1.0-dpixl) + val[3]*dpixl)); 00283 00284 00285 /* solar azimuth interpolation */ 00286 00287 val[0] = (double)sola1000[iline+0][ipixl+0]; 00288 val[1] = (double)sola1000[iline+0][ipixl+1]; 00289 val[2] = (double)sola1000[iline+1][ipixl+0]; 00290 val[3] = (double)sola1000[iline+1][ipixl+1]; 00291 00292 *sola = ((1.0-dline)*(val[0]*(1.0-dpixl) + val[1]*dpixl) + 00293 dline*(val[2]*(1.0-dpixl) + val[3]*dpixl)); 00294 00295 00296 /* sensor zenith interpolation */ 00297 00298 val[0] = (double)senz1000[iline+0][ipixl+0]; 00299 val[1] = (double)senz1000[iline+0][ipixl+1]; 00300 val[2] = (double)senz1000[iline+1][ipixl+0]; 00301 val[3] = (double)senz1000[iline+1][ipixl+1]; 00302 00303 *senz = ((1.0-dline)*(val[0]*(1.0-dpixl) + val[1]*dpixl) + 00304 dline*(val[2]*(1.0-dpixl) + val[3]*dpixl)); 00305 00306 00307 /* sensor azimuth interpolation, handle singularity */ 00308 00309 val[0] = (double)sena1000[iline+0][ipixl+0]; 00310 val[1] = (double)sena1000[iline+0][ipixl+1]; 00311 val[2] = (double)sena1000[iline+1][ipixl+0]; 00312 val[3] = (double)sena1000[iline+1][ipixl+1]; 00313 00314 if (*senz < 0.05) 00315 *sena = 0.0; 00316 else 00317 *sena = ((1.0-dline)*(val[0]*(1.0-dpixl) + val[1]*dpixl) + 00318 dline*(val[2]*(1.0-dpixl) + val[3]*dpixl)); 00319 } 00320 00321 00322 /* ----------------------------------------------------------------------------------- */ 00323 /* subframe_calibration() - loads subframe calibration factors and applies. */ 00324 /* */ 00325 /* B. Franz, SAIC, February 2003. */ 00326 /* ----------------------------------------------------------------------------------- */ 00327 void subframe_calibration(int32_t sensorID, int32_t ib, float mrad, float brad, 00328 int32_t spix, int32_t npix, float32 data[]) 00329 { 00330 #define NSFWAVE 7 00331 #define MAXSF 4 00332 #define MAXRAD 100 00333 00334 static int firstCall = 1; 00335 static float sftable[NSFWAVE][1+MAXSF][MAXRAD]; 00336 static int twave [NSFWAVE] = {469,555,645,859,1240,1640,2130}; 00337 static int nsfsub [NSFWAVE] = {2,2,4,4,2,2,2}; 00338 static int nsfrad [NSFWAVE]; 00339 static int32_t iband [NSFWAVE]; 00340 00341 float rad, sfcor; 00342 int32_t irad, ip, iw, ipb, isf; 00343 00344 if (firstCall) { 00345 char *tmp_str; 00346 char line [81]; 00347 char path [FILENAME_MAX] = ""; 00348 char file [FILENAME_MAX] = ""; 00349 FILE *fp = NULL; 00350 00351 if ((tmp_str = getenv("OCDATAROOT")) == NULL) { 00352 printf("OCDATAROOT environment variable is not defined.\n"); 00353 exit(1); 00354 } 00355 00356 strcpy(path,tmp_str); strcat(path,"/"); strcat(path,sensorDir[sensorID]); strcat(path,"/cal/"); 00357 00358 for (iw=0; iw<NSFWAVE; iw++) { 00359 iband[iw] = bindex_get(twave[iw]); 00360 sprintf(file,"%s%s%s%s%3d%s",path,"subframecor_",sensorDir[sensorID],"_",twave[iw],".dat"); 00361 if(want_verbose) 00362 printf("\nLoading subframe correction table from %s",file); 00363 if ( (fp = fopen(file,"r")) == NULL ) { 00364 fprintf(stderr, 00365 "-E- %s line %d: unable to open %s for reading\n",__FILE__,__LINE__,file); 00366 exit(1); 00367 } 00368 00369 while ( fgets( line, 80, fp ) ) { 00370 if ( strncmp(line,"/end_header",11) == 0 ) 00371 break; 00372 } 00373 00374 irad = 0; 00375 while ( fgets( line, 80, fp ) ) { 00376 if (nsfsub[iw] == 2) 00377 sscanf(line,"%f %f %f",&sftable[iw][0][irad], 00378 &sftable[iw][1][irad],&sftable[iw][2][irad]); 00379 else 00380 sscanf(line,"%f %f %f %f %f",&sftable[iw][0][irad], 00381 &sftable[iw][1][irad],&sftable[iw][2][irad],&sftable[iw][3][irad],&sftable[iw][4][irad]); 00382 irad++; 00383 if (irad == MAXRAD) { 00384 fprintf(stderr, 00385 "-E- %s line %d: subframe correction table exceeded allocation.\n", 00386 __FILE__,__LINE__); 00387 exit(1); 00388 } 00389 } 00390 nsfrad[iw] = irad; 00391 fclose(fp); 00392 } 00393 00394 if(want_verbose) 00395 printf("\nSubframe destriping corrections enabled.\n\n"); 00396 00397 firstCall = 0; 00398 } 00399 00400 00401 for (iw=0; iw<NSFWAVE; iw++) 00402 if (ib == iband[iw]) 00403 break; 00404 if (iw == NSFWAVE) { 00405 fprintf(stderr,"-E- %s line %d: attempt to correct non-subframed band %d.\n",__FILE__,__LINE__,ib); 00406 return; 00407 } 00408 00409 for (ip=0; ip<npix; ip++) { 00410 if (data[ip] < 65000) { 00411 rad = (data[ip] - brad)*mrad; 00412 isf = (ip + spix) % nsfsub[iw]; 00413 for (irad=0; irad<nsfrad[iw]-1; irad++) { 00414 if (rad < sftable[iw][0][irad]) break; 00415 } 00416 if (irad == 0 || irad == (nsfrad[iw]-1)) 00417 sfcor = sftable[iw][isf+1][irad]; 00418 else { 00419 sfcor = sftable[iw][isf+1][irad-1] 00420 + (rad - sftable[iw][0][irad-1]) 00421 * (sftable[iw][isf+1][irad]-sftable[iw][isf+1][irad-1]) / 00422 (sftable[iw][0][irad]-sftable[iw][0][irad-1]); 00423 } 00424 data[ip] /= sfcor; 00425 } 00426 } 00427 } 00428 00429 00430 00431 /* ----------------------------------------------------------------------------------- */ 00432 /* openl1_hmodis_hdf() - opens a MODIS L1B file for reading. */ 00433 /* */ 00434 /* B. Franz, SAIC, February 2003. */ 00435 /* ----------------------------------------------------------------------------------- */ 00436 int openl1_hmodis_hdf(filehandle *file) 00437 { 00438 int32 npix; 00439 int32 nscan; 00440 00441 int32 sds_id; 00442 int32 rank; 00443 int32 dims[3]; 00444 int32 type; 00445 int32 numattr; 00446 00447 int32 itemp; 00448 uint32 gflags[8]; 00449 char cdata[32]; 00450 char namebuf [FILENAME_MAX]; 00451 char name250 [FILENAME_MAX]; 00452 char name500 [FILENAME_MAX]; 00453 char name1000[FILENAME_MAX]; 00454 char dirbuf [FILENAME_MAX]; 00455 char temp [FILENAME_MAX]; 00456 char *p = NULL; 00457 int convention = 0; 00458 00459 memset(namebuf, '\0',FILENAME_MAX); 00460 memset(name250, '\0',FILENAME_MAX); 00461 memset(name500, '\0',FILENAME_MAX); 00462 memset(name1000,'\0',FILENAME_MAX); 00463 00464 strcpy(temp, file->name); 00465 strcpy(dirbuf, dirname(temp)); 00466 strcat(dirbuf, "/"); 00467 strcpy(temp, file->name); 00468 strcpy(namebuf, basename(temp)); 00469 00470 if ((p = strstr(namebuf,"QKM")) != NULL) { 00471 if (file->input->resolution == -1) 00472 file->input->resolution = 250; 00473 } else if ((p = strstr(namebuf,"HKM")) != NULL) { 00474 if (file->input->resolution == -1) 00475 file->input->resolution = 500; 00476 } else if ((p = strstr(namebuf,"1KM")) != NULL) { 00477 if (file->input->resolution == -1) 00478 file->input->resolution = 1000; 00479 convention = 1; 00480 } else if ((p = strstr(namebuf,"LAC")) != NULL) { 00481 if (file->input->resolution == -1) 00482 file->input->resolution = 1000; 00483 convention = 2; 00484 } else { 00485 p = namebuf+strlen(namebuf); 00486 if (file->input->resolution == -1) 00487 file->input->resolution = 1000; 00488 convention = 3; 00489 } 00490 00491 strcpy (name250, dirbuf); 00492 strncat(name250, namebuf,p-namebuf); 00493 strcpy (name500, dirbuf); 00494 strncat(name500, namebuf,p-namebuf); 00495 strcpy (name1000,dirbuf); 00496 strncat(name1000,namebuf,p-namebuf); 00497 00498 /* Capture processing resolution */ 00499 resolution = file->input->resolution; 00500 if(want_verbose) 00501 printf("Processing at %d meter resolution.\n",resolution); 00502 00503 /* If not following std naming conventions, assume a 1km file */ 00504 if (p == NULL) { 00505 sd_id[2] = SDstart(file->name, DFACC_RDONLY); 00506 if(sd_id[2] == FAIL){ 00507 fprintf(stderr,"-E- %s line %d: SDstart(%s, %d) failed.\n", 00508 __FILE__,__LINE__,file->name,DFACC_RDONLY); 00509 return(HDF_FUNCTION_ERROR); 00510 } 00511 if(want_verbose) 00512 printf(" 1000-meter file: %s\n",file->name); 00513 00514 } else { 00515 00516 if (resolution < 500) { 00517 if (convention != 3) { 00518 strcat(name250,"QKM"); 00519 strcat(name250,p+3); 00520 } else { 00521 strcat(name250,"_QKM"); 00522 } 00523 if(want_verbose) 00524 printf(" 250-meter file: %s\n",name250); 00525 sd_id[0] = SDstart(name250, DFACC_RDONLY); 00526 if(sd_id[0] == FAIL){ 00527 fprintf(stderr,"-E- %s line %d: SDstart(%s, %d) failed.\n", 00528 __FILE__,__LINE__,name250,DFACC_RDONLY); 00529 return(HDF_FUNCTION_ERROR); 00530 } 00531 } 00532 00533 if (resolution < 1000) { 00534 if (convention != 3) { 00535 strcat(name500,"HKM"); 00536 strcat(name500,p+3); 00537 } else { 00538 strcat(name500,"_HKM"); 00539 } 00540 if(want_verbose) 00541 printf(" 500-meter file: %s\n",name500); 00542 sd_id[1] = SDstart(name500, DFACC_RDONLY); 00543 if(sd_id[1] == FAIL){ 00544 fprintf(stderr,"-E- %s line %d: SDstart(%s, %d) failed.\n", 00545 __FILE__,__LINE__,name500,DFACC_RDONLY); 00546 return(HDF_FUNCTION_ERROR); 00547 } 00548 } 00549 00550 switch (convention) { 00551 case 1: 00552 strcat(name1000,"1KM"); 00553 strcat(name1000,p+3); 00554 break; 00555 case 2: 00556 strcat(name1000,"LAC"); 00557 strcat(name1000,p+3); 00558 break; 00559 case 3: 00560 break; 00561 default: 00562 strcat(name1000,"1KM"); 00563 strcat(name1000,p+3); 00564 sd_id[2] = SDstart(name1000, DFACC_RDONLY); 00565 if(sd_id[2] == FAIL) { 00566 memset(name1000,'\0',FILENAME_MAX); 00567 memcpy(name1000,namebuf,p-namebuf); 00568 strcat(name1000,"LAC"); 00569 strcat(name1000,p+3); 00570 } 00571 break; 00572 } 00573 if(want_verbose) 00574 printf(" 1000-meter file: %s\n",name1000); 00575 sd_id[2] = SDstart(name1000, DFACC_RDONLY); 00576 if(sd_id[2] == FAIL){ 00577 fprintf(stderr,"-E- %s line %d: SDstart(%s, %d) failed.\n", 00578 __FILE__,__LINE__,name1000,DFACC_RDONLY); 00579 return(HDF_FUNCTION_ERROR); 00580 } 00581 00582 /* Make sure this is not an ocean-band sub-setted L1 file */ 00583 if (getHDFattr(sd_id[2],"Rescaled Ocean R","EV_250_Aggr1km_RefSB",(VOIDP)&itemp) == 0) { 00584 printf("\nThis L1B file appears to have been derived from an ocean-band subsetted L1A file,\n"); 00585 printf("so it does not contain any data from the MODIS Land/Cloud bands, at any resolution.\n"); 00586 printf("Set resolution=-1 for standard ocean processing.\n"); 00587 exit(1); 00588 } 00589 } 00590 if(want_verbose) 00591 printf("\n"); 00592 00593 00594 /* Get pixel and scan dimensions */ 00595 switch (resolution) { 00596 case 250: 00597 sds_id = SDselect(sd_id[0],SDnametoindex(sd_id[0],("EV_250_RefSB"))); 00598 if (SDgetinfo(sds_id,NULL,&rank,dims,&type,&numattr) == -1) { 00599 fprintf(stderr,"-E- %s line %d: error getting dimension info.\n", 00600 __FILE__,__LINE__); 00601 return(HDF_FUNCTION_ERROR); 00602 } 00603 npix = dims[2]; 00604 nscan = dims[1]; 00605 if (getHDFattr(sd_id[0],"Extract Pixel Count","",(VOIDP)&itemp) == 0 && itemp > 0) { 00606 npix = itemp*4; 00607 getHDFattr(sd_id[0],"Extract Pixel Offset","",(VOIDP)&itemp); 00608 spix = itemp*4; 00609 getHDFattr(sd_id[0],"Extract Line Offset","",(VOIDP)&itemp); 00610 sscan = itemp*nd250; 00611 have_extract = 1; 00612 if(want_verbose) { 00613 printf("File was generated from L1A extract of size %d x %d\n",npix,nscan); 00614 printf(" Pixels %d - %d\n",spix+1,spix+npix); 00615 printf(" Lines %d - %d\n",sscan+1,sscan+nscan); 00616 } // want_verbose 00617 } 00618 ndet = nd250; 00619 np1000 = npix/4; 00620 sp1000 = spix/4; 00621 np500 = npix/2; 00622 sp500 = spix/2; 00623 get_modis_calfile(sd_id[0],file->calfile); /* read list of inputs to L1B */ 00624 break; 00625 case 500: 00626 sds_id = SDselect(sd_id[1],SDnametoindex(sd_id[1],("EV_500_RefSB"))); 00627 if (SDgetinfo(sds_id,NULL,&rank,dims,&type,&numattr) == -1) { 00628 fprintf(stderr,"-E- %s line %d: error getting dimension info.\n", 00629 __FILE__,__LINE__); 00630 return(HDF_FUNCTION_ERROR); 00631 } 00632 npix = dims[2]; 00633 nscan = dims[1]; 00634 if (getHDFattr(sd_id[1],"Extract Pixel Count","",(VOIDP)&itemp) == 0 && itemp > 0) { 00635 npix = itemp*2; 00636 getHDFattr(sd_id[1],"Extract Pixel Offset","",(VOIDP)&itemp); 00637 spix = itemp*2; 00638 getHDFattr(sd_id[1],"Extract Line Offset","",(VOIDP)&itemp); 00639 sscan = itemp*nd500; 00640 have_extract = 1; 00641 if(want_verbose) { 00642 printf("File was generated from L1A extract of size %d x %d\n",npix,nscan); 00643 printf(" Pixels %d - %d\n",spix+1,spix+npix); 00644 printf(" Lines %d - %d\n",sscan+1,sscan+nscan); 00645 } // want_verbose 00646 } 00647 ndet = nd500; 00648 np1000 = npix/2; 00649 sp1000 = spix/2; 00650 np500 = npix; 00651 sp500 = spix; 00652 get_modis_calfile(sd_id[1],file->calfile); /* read list of inputs to L1B */ 00653 break; 00654 default: 00655 sds_id = SDselect(sd_id[2],SDnametoindex(sd_id[2],("EV_1KM_RefSB"))); 00656 if (SDgetinfo(sds_id,NULL,&rank,dims,&type,&numattr) == -1) { 00657 fprintf(stderr,"-E- %s line %d: error getting dimension info.\n", 00658 __FILE__,__LINE__); 00659 return(HDF_FUNCTION_ERROR); 00660 } 00661 npix = dims[2]; 00662 nscan = dims[1]; 00663 if (getHDFattr(sd_id[2],"Extract Pixel Count","",(VOIDP)&itemp) == 0 && itemp > 0) { 00664 npix = itemp; 00665 getHDFattr(sd_id[2],"Extract Pixel Offset","",(VOIDP)&itemp); 00666 spix = itemp; 00667 getHDFattr(sd_id[2],"Extract Line Offset","",(VOIDP)&itemp); 00668 sscan = itemp*nd1000; 00669 have_extract = 1; 00670 if(want_verbose) { 00671 printf("File was generated from L1A extract of size %d x %d\n",npix,nscan); 00672 printf(" Pixels %d - %d\n",spix+1,spix+npix); 00673 printf(" Lines %d - %d\n",sscan+1,sscan+nscan); 00674 } // want_verbose 00675 } 00676 ndet = nd1000; 00677 np1000 = npix; 00678 sp1000 = spix; 00679 np500 = npix; 00680 sp500 = spix; 00681 get_modis_calfile(sd_id[2],file->calfile); /* read list of inputs to L1B */ 00682 break; 00683 } 00684 00685 file->npix = npix; 00686 file->ndets = ndet; 00687 file->nscan = nscan; 00688 file->sd_id = sd_id[2]; 00689 00690 00691 /* Open the HDF geolocation input file.*/ 00692 sd_id_g = SDstart(file->geofile, DFACC_RDONLY); 00693 if(sd_id_g == FAIL){ 00694 printf("Error opening geolocation file.\n"); 00695 fprintf(stderr,"-E- %s line %d: SDstart(%s, %d) failed.\n", 00696 __FILE__,__LINE__,file->geofile,DFACC_RDONLY); 00697 return(HDF_FUNCTION_ERROR); 00698 } 00699 00700 /* Has geolocation been corrected for terrain height? */ 00701 if (getHDFattr(sd_id_g,"Cumulated gflags","",(VOIDP)&gflags) == 0) { 00702 file->terrain_corrected = (int32_t)(gflags[5]==0); 00703 } else printf("-E- %s line %d: Error reading gflags.\n",__FILE__,__LINE__); 00704 00705 get_modis_orbit(sd_id_g, file); /* read orbit info */ 00706 00707 00708 if (file->sensorID == HMODISA && resolution == 250) destripe=1; 00709 00710 return(LIFE_IS_GOOD); 00711 } 00712 00713 00714 /* ----------------------------------------------------------------------------------- */ 00715 /* readl1_hmodis_hdf() - reads 1 line (scan) from a MODIS L1B file, loads l1rec. */ 00716 /* */ 00717 /* 00718 Pos Band Wavelength Resolution EV_Position 00719 0 8 412 1000 0 00720 1 9 443 1000 1 00721 2 3 469 500 0 00722 3 10 488 1000 2 00723 4 11 531 1000 3 00724 5 12 551 1000 4 00725 6 4 555 500 1 00726 7 1 645 250 0 00727 8 13 667lo 1000 5 00728 9 14 678lo 1000 7 00729 10 15 748 1000 9 00730 11 2 858 250 1 00731 12 16 869 1000 10 00732 13 5 1240 500 2 00733 14 6 1640 500 3 00734 15 7 2130 500 4 00735 00736 16 20 3750 1000 00737 17 22 3959 1000 00738 18 23 4050 1000 00739 19 27 6715 1000 00740 20 28 7325 1000 00741 21 29 8550 1000 00742 22 31 11030 1000 00743 23 32 12020 1000 00744 */ 00745 00746 /* */ 00747 /* B. Franz, SAIC, February 2006. */ 00748 /* ----------------------------------------------------------------------------------- */ 00749 int readl1_hmodis_hdf(filehandle *file, int32 scan, l1str *l1rec, int lonlat) 00750 { 00751 static int firstCall = 1; 00752 static int firstScan = 1; 00753 static float pi = PI; 00754 00755 static uint16 *idata = NULL; 00756 static float32 *fdata = NULL; 00757 00758 static int16 *hgt = NULL; 00759 static int16 *senz = NULL; 00760 static int16 *sena = NULL; 00761 static int16 *solz = NULL; 00762 static int16 *sola = NULL; 00763 00764 int32 rank; 00765 int32 dims[3]; 00766 int32 type; 00767 int32 numattr; 00768 int32 sds_id; 00769 00770 static float m [MAXBANDS]; 00771 static float b [MAXBANDS]; 00772 static float mrad[MAXBANDS]; 00773 static float brad[MAXBANDS]; 00774 static float m_cirrus; 00775 static float b_cirrus; 00776 00777 float floatbuf[15]; 00778 00779 static int nb250 = 2; 00780 static int ib250 [] = {7,11}; 00781 static int pb250 [] = {0,1}; 00782 00783 static int nb500 = 5; 00784 static int ib500 [] = {2,6,13,14,15}; 00785 static int pb500 [] = {0,1,2,3,4}; 00786 00787 static int nb1000 = 9; 00788 static int ib1000[] = {0,1,3,4,5,8,9,10,12}; 00789 static int pb1000[] = {0,1,2,3,4,5,7, 9,10}; 00790 00791 static int nbir = 8; 00792 static int ibir [] = {16,17,18,19,20,21,22,23}; 00793 static int pbir [] = {0,2,3,6,7,8,10,11}; 00794 00795 static int nbcir = 1; 00796 static int ibcir [] = {24}; 00797 static int pbcir [] = {0}; 00798 00799 static float64 lastTAIsec = -1; 00800 static float64 TAIsec, usec, dsec; 00801 static int16 year, day; 00802 static int32 lday; 00803 static double fsol; 00804 static int16 mside, badmside; 00805 static double mnorm[3]; 00806 static int32 lastframe = -1; 00807 static int32 last_work_scan; 00808 00809 int32 npix = (int32)file->npix; 00810 int32 nbands = (int32)file->nbands; 00811 int32 frame = (rscan+scan) / ndet; 00812 int32 detnum = (rscan+scan) % ndet; 00813 int32 i,ip,iw,ib,ipb,status; 00814 int32 idet, iline; 00815 00816 static float *Fobar; 00817 double hour; 00818 double esdist; 00819 00820 if (firstCall) { 00821 00822 firstCall = 0; 00823 00824 /* We need this for the reflectance to radiance conversion */ 00825 rdsensorinfo(file->sensorID,file->input->evalmask,"Fobar", (void **) &Fobar); 00826 00827 /* Buffer for one line, all bands */ 00828 if ((idata = (uint16 *) calloc(npix*MAXBANDS, sizeof(uint16))) == NULL) { 00829 printf("-E- %s line %d: Error allocating data space.\n", 00830 __FILE__,__LINE__); 00831 return(1); 00832 } 00833 if ((fdata = (float32 *) calloc(npix*MAXBANDS, sizeof(float32))) == NULL) { 00834 printf("-E- %s line %d: Error allocating data space.\n", 00835 __FILE__,__LINE__); 00836 return(1); 00837 } 00838 00839 /* Need to buffer a full (10-det) scan of 1km geolocation */ 00840 if ((lon1000 = alloc2d_float(np1000,nd1000)) == NULL) { 00841 printf("-E- %s line %d: Error allocating geolocation space.\n", 00842 __FILE__,__LINE__); 00843 return(1); 00844 } 00845 if ((lat1000 = alloc2d_float(np1000,nd1000)) == NULL) { 00846 printf("-E- %s line %d: Error allocating geolocation space.\n", 00847 __FILE__,__LINE__); 00848 return(1); 00849 } 00850 if ((hgt1000 = alloc2d_float(np1000,nd1000)) == NULL) { 00851 printf("-E- %s line %d: Error allocating geolocation space.\n", 00852 __FILE__,__LINE__); 00853 return(1); 00854 } 00855 if ((solz1000 = alloc2d_float(np1000,nd1000)) == NULL) { 00856 printf("-E- %s line %d: Error allocating geolocation space.\n", 00857 __FILE__,__LINE__); 00858 return(1); 00859 } 00860 if ((sola1000 = alloc2d_float(np1000,nd1000)) == NULL) { 00861 printf("-E- %s line %d: Error allocating geolocation space.\n", 00862 __FILE__,__LINE__); 00863 return(1); 00864 } 00865 if ((senz1000 = alloc2d_float(np1000,nd1000)) == NULL) { 00866 printf("-E- %s line %d: Error allocating geolocation space.\n", 00867 __FILE__,__LINE__); 00868 return(1); 00869 } 00870 if ((sena1000 = alloc2d_float(np1000,nd1000)) == NULL) { 00871 printf("-E- %s line %d: Error allocating geolocation space.\n", 00872 __FILE__,__LINE__); 00873 return(1); 00874 } 00875 00876 /* Need buffers for single scan of scaled values */ 00877 if ((hgt = (int16 *) calloc(np1000, sizeof(int16))) == NULL) { 00878 printf("-E- %s line %d: Error allocating hgt space.\n", 00879 __FILE__,__LINE__); 00880 return(1); 00881 } 00882 if ((solz = (int16 *) calloc(np1000, sizeof(int16))) == NULL) { 00883 printf("-E- %s line %d: Error allocating solz space.\n", 00884 __FILE__,__LINE__); 00885 return(1); 00886 } 00887 if ((sola = (int16 *) calloc(np1000, sizeof(int16))) == NULL) { 00888 printf("-E- %s line %d: Error allocating sola space.\n", 00889 __FILE__,__LINE__); 00890 return(1); 00891 } 00892 if ((senz = (int16 *) calloc(np1000, sizeof(int16))) == NULL) { 00893 printf("-E- %s line %d: Error allocating senz space.\n", 00894 __FILE__,__LINE__); 00895 return(1); 00896 } 00897 if ((sena = (int16 *) calloc(np1000, sizeof(int16))) == NULL) { 00898 printf("-E- %s line %d: Error allocating sena space.\n", 00899 __FILE__,__LINE__); 00900 return(1); 00901 } 00902 00903 if (resolution == 250) { 00904 00905 /* reflectance scale and offset factors */ 00906 00907 if (getHDFattr(sd_id[0],"reflectance_scales","EV_250_RefSB",(VOIDP)floatbuf) != 0) { 00908 printf("-E- %s line %d: Error reading reflectance scale attribute.\n", 00909 __FILE__,__LINE__); 00910 return(1); 00911 } 00912 for (ib=0; ib<nb250; ib++) 00913 m[ib250[ib]] = floatbuf[pb250[ib]]; 00914 if (getHDFattr(sd_id[0],"reflectance_offsets","EV_250_RefSB",(VOIDP)floatbuf) != 0) { 00915 printf("-E- %s line %d: Error reading reflectance offset attribute.\n", 00916 __FILE__,__LINE__); 00917 return(1); 00918 } 00919 for (ib=0; ib<nb250; ib++) 00920 b[ib250[ib]] = floatbuf[pb250[ib]]; 00921 00922 if (getHDFattr(sd_id[1],"reflectance_scales","EV_500_RefSB",(VOIDP)floatbuf) != 0) { 00923 printf("-E- %s line %d: Error reading reflectance scale attribute.\n", 00924 __FILE__,__LINE__); 00925 return(1); 00926 } 00927 for (ib=0; ib<nb500; ib++) 00928 m[ib500[ib]] = floatbuf[pb500[ib]]; 00929 00930 if (getHDFattr(sd_id[1],"reflectance_offsets","EV_500_RefSB",(VOIDP)floatbuf) != 0) { 00931 printf("-E- %s line %d: Error reading reflectance offset attribute.\n", 00932 __FILE__,__LINE__); 00933 return(1); 00934 } 00935 for (ib=0; ib<nb500; ib++) 00936 b[ib500[ib]] = floatbuf[pb500[ib]]; 00937 00938 if (getHDFattr(sd_id[2],"reflectance_scales","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 00939 printf("-E- %s line %d: Error reading reflectance scale attribute.\n", 00940 __FILE__,__LINE__); 00941 return(1); 00942 } 00943 for (ib=0; ib<nb1000; ib++) 00944 m[ib1000[ib]] = floatbuf[pb1000[ib]]; 00945 00946 if (getHDFattr(sd_id[2],"reflectance_offsets","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 00947 printf("-E- %s line %d: Error reading reflectance offset attribute.\n", 00948 __FILE__,__LINE__); 00949 return(1); 00950 } 00951 for (ib=0; ib<nb1000; ib++) 00952 b[ib1000[ib]] = floatbuf[pb1000[ib]]; 00953 00954 /* radiance scale and offset factors (for destriping subframes) */ 00955 00956 if (getHDFattr(sd_id[0],"radiance_scales","EV_250_RefSB",(VOIDP)floatbuf) != 0) { 00957 printf("-E- %s line %d: Error reading radiance scale attribute.\n", 00958 __FILE__,__LINE__); 00959 return(1); 00960 } 00961 for (ib=0; ib<nb250; ib++) 00962 mrad[ib250[ib]] = floatbuf[pb250[ib]]; 00963 00964 if (getHDFattr(sd_id[0],"radiance_offsets","EV_250_RefSB",(VOIDP)floatbuf) != 0) { 00965 printf("-E- %s line %d: Error reading radiance offset attribute.\n", 00966 __FILE__,__LINE__); 00967 return(1); 00968 } 00969 for (ib=0; ib<nb250; ib++) 00970 brad[ib250[ib]] = floatbuf[pb250[ib]]; 00971 00972 if (getHDFattr(sd_id[1],"radiance_scales","EV_500_RefSB",(VOIDP)floatbuf) != 0) { 00973 printf("-E- %s line %d: Error reading radiance scale attribute.\n", 00974 __FILE__,__LINE__); 00975 return(1); 00976 } 00977 for (ib=0; ib<nb500; ib++) 00978 mrad[ib500[ib]] = floatbuf[pb500[ib]]; 00979 00980 if (getHDFattr(sd_id[1],"radiance_offsets","EV_500_RefSB",(VOIDP)floatbuf) != 0) { 00981 printf("-E- %s line %d: Error reading radiance offset attribute.\n", 00982 __FILE__,__LINE__); 00983 return(1); 00984 } 00985 for (ib=0; ib<nb500; ib++) 00986 brad[ib500[ib]] = floatbuf[pb500[ib]]; 00987 00988 if (getHDFattr(sd_id[2],"radiance_scales","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 00989 printf("-E- %s line %d: Error reading radiance scale attribute.\n", 00990 __FILE__,__LINE__); 00991 return(1); 00992 } 00993 for (ib=0; ib<nb1000; ib++) 00994 mrad[ib1000[ib]] = floatbuf[pb1000[ib]]; 00995 00996 if (getHDFattr(sd_id[2],"radiance_offsets","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 00997 printf("-E- %s line %d: Error reading radiance offset attribute.\n", 00998 __FILE__,__LINE__); 00999 return(1); 01000 } 01001 for (ib=0; ib<nb1000; ib++) 01002 brad[ib1000[ib]] = floatbuf[pb1000[ib]]; 01003 01004 01005 01006 /* allocate data buffers for interpolation */ 01007 01008 for (ib=0; ib<nb1000; ib++) { 01009 if ((fdatabuf[ib1000[ib]] = (float32 **) alloc2d_float(np1000,nd1000)) == NULL) { 01010 printf("-E- %s line %d: Error allocating data buffer space for band %d\n", 01011 __FILE__,__LINE__,ib1000[ib]); 01012 exit(1); 01013 } 01014 } 01015 01016 for (ib=0; ib<nbir; ib++) { 01017 if ((fdatabuf[ibir[ib]] = (float32 **) alloc2d_float(np1000,nd1000)) == NULL) { 01018 printf("-E- %s line %d: Error allocating data buffer space for band %d\n", 01019 __FILE__,__LINE__,ibir[ib]); 01020 exit(1); 01021 } 01022 } 01023 01024 for (ib=0; ib<nbcir; ib++) { 01025 if ((fdatabuf[ibcir[ib]] = (float32 **) alloc2d_float(np1000,nd1000)) == NULL) { 01026 printf("-E- %s line %d: Error allocating data buffer space for band %d\n", 01027 __FILE__,__LINE__,ibcir[ib]); 01028 exit(1); 01029 } 01030 } 01031 01032 for (ib=0; ib<nb500; ib++) { 01033 if ((fdatabuf[ib500[ib]] = (float32 **) alloc2d_float(np500,nd500)) == NULL) { 01034 printf("-E- %s line %d: Error allocating data buffer space for band %d\n", 01035 __FILE__,__LINE__,ib500[ib]); 01036 exit(1); 01037 } 01038 } 01039 01040 } else if (resolution == 500) { 01041 01042 /* reflectance scale and offset factors */ 01043 01044 if (getHDFattr(sd_id[1],"reflectance_scales","EV_250_Aggr500_RefSB",(VOIDP)floatbuf) != 0) { 01045 printf("-E- %s line %d: Error reading reflectance scale attribute.\n", 01046 __FILE__,__LINE__); 01047 return(1); 01048 } 01049 for (ib=0; ib<nb250; ib++) 01050 m[ib250[ib]] = floatbuf[pb250[ib]]; 01051 if (getHDFattr(sd_id[1],"reflectance_offsets","EV_250_Aggr500_RefSB",(VOIDP)floatbuf) != 0) { 01052 printf("-E- %s line %d: Error reading reflectance offset attribute.\n", 01053 __FILE__,__LINE__); 01054 return(1); 01055 } 01056 for (ib=0; ib<nb250; ib++) 01057 b[ib250[ib]] = floatbuf[pb250[ib]]; 01058 01059 if (getHDFattr(sd_id[1],"reflectance_scales","EV_500_RefSB",(VOIDP)floatbuf) != 0) { 01060 printf("-E- %s line %d: Error reading reflectance scale attribute.\n", 01061 __FILE__,__LINE__); 01062 return(1); 01063 } 01064 for (ib=0; ib<nb500; ib++) 01065 m[ib500[ib]] = floatbuf[pb500[ib]]; 01066 01067 if (getHDFattr(sd_id[1],"reflectance_offsets","EV_500_RefSB",(VOIDP)floatbuf) != 0) { 01068 printf("-E- %s line %d: Error reading reflectance offset attribute.\n", 01069 __FILE__,__LINE__); 01070 return(1); 01071 } 01072 for (ib=0; ib<nb500; ib++) 01073 b[ib500[ib]] = floatbuf[pb500[ib]]; 01074 01075 if (getHDFattr(sd_id[2],"reflectance_scales","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 01076 printf("-E- %s line %d: Error reading reflectance scale attribute.\n", 01077 __FILE__,__LINE__); 01078 return(1); 01079 } 01080 for (ib=0; ib<nb1000; ib++) 01081 m[ib1000[ib]] = floatbuf[pb1000[ib]]; 01082 01083 if (getHDFattr(sd_id[2],"reflectance_offsets","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 01084 printf("-E- %s line %d: Error reading reflectance offset attribute.\n", 01085 __FILE__,__LINE__); 01086 return(1); 01087 } 01088 for (ib=0; ib<nb1000; ib++) 01089 b[ib1000[ib]] = floatbuf[pb1000[ib]]; 01090 01091 if (getHDFattr(sd_id[1],"reflectance_scales","EV_250_Aggr500_RefSB",(VOIDP)floatbuf) != 0) { 01092 printf("-E- %s line %d: Error reading reflectance scale attribute.\n", 01093 __FILE__,__LINE__); 01094 return(1); 01095 } 01096 01097 01098 /* radiance scale and offset factors (for destriping subframes) */ 01099 01100 for (ib=0; ib<nb250; ib++) 01101 mrad[ib250[ib]] = floatbuf[pb250[ib]]; 01102 if (getHDFattr(sd_id[1],"radiance_offsets","EV_250_Aggr500_RefSB",(VOIDP)floatbuf) != 0) { 01103 printf("-E- %s line %d: Error reading radiance offset attribute.\n", 01104 __FILE__,__LINE__); 01105 return(1); 01106 } 01107 for (ib=0; ib<nb250; ib++) 01108 brad[ib250[ib]] = floatbuf[pb250[ib]]; 01109 01110 if (getHDFattr(sd_id[1],"radiance_scales","EV_500_RefSB",(VOIDP)floatbuf) != 0) { 01111 printf("-E- %s line %d: Error reading radiance scale attribute.\n", 01112 __FILE__,__LINE__); 01113 return(1); 01114 } 01115 for (ib=0; ib<nb500; ib++) 01116 mrad[ib500[ib]] = floatbuf[pb500[ib]]; 01117 01118 if (getHDFattr(sd_id[1],"radiance_offsets","EV_500_RefSB",(VOIDP)floatbuf) != 0) { 01119 printf("-E- %s line %d: Error reading radiance offset attribute.\n", 01120 __FILE__,__LINE__); 01121 return(1); 01122 } 01123 for (ib=0; ib<nb500; ib++) 01124 brad[ib500[ib]] = floatbuf[pb500[ib]]; 01125 01126 if (getHDFattr(sd_id[2],"radiance_scales","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 01127 printf("-E- %s line %d: Error reading radiance scale attribute.\n", 01128 __FILE__,__LINE__); 01129 return(1); 01130 } 01131 for (ib=0; ib<nb1000; ib++) 01132 mrad[ib1000[ib]] = floatbuf[pb1000[ib]]; 01133 01134 if (getHDFattr(sd_id[2],"radiance_offsets","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 01135 printf("-E- %s line %d: Error reading radiance offset attribute.\n", 01136 __FILE__,__LINE__); 01137 return(1); 01138 } 01139 for (ib=0; ib<nb1000; ib++) 01140 brad[ib1000[ib]] = floatbuf[pb1000[ib]]; 01141 01142 /* allocate data buffers for interpolation */ 01143 01144 for (ib=0; ib<nb1000; ib++) { 01145 if ((fdatabuf[ib1000[ib]] = (float32 **) alloc2d_float(np1000,nd1000)) == NULL) { 01146 printf("-E- %s line %d: Error allocating data buffer space for band %d\n", 01147 __FILE__,__LINE__,ib1000[ib]); 01148 exit(1); 01149 } 01150 } 01151 01152 for (ib=0; ib<nbir; ib++) { 01153 if ((fdatabuf[ibir[ib]] = (float32 **) alloc2d_float(np1000,nd1000)) == NULL) { 01154 printf("-E- %s line %d: Error allocating data buffer space for band %d\n", 01155 __FILE__,__LINE__,ibcir[ib]); 01156 exit(1); 01157 } 01158 } 01159 01160 for (ib=0; ib<nbcir; ib++) { 01161 if ((fdatabuf[ibcir[ib]] = (float32 **) alloc2d_float(np1000,nd1000)) == NULL) { 01162 printf("-E- %s line %d: Error allocating data buffer space for band %d\n", 01163 __FILE__,__LINE__,ibcir[ib]); 01164 exit(1); 01165 } 01166 } 01167 01168 } else { /* 1KM */ 01169 01170 /* reflectance scale and offset factors */ 01171 01172 if (getHDFattr(sd_id[2],"reflectance_scales","EV_250_Aggr1km_RefSB",(VOIDP)floatbuf) != 0) { 01173 printf("-E- %s line %d: Error reading reflectance scale attribute.\n", 01174 __FILE__,__LINE__); 01175 return(1); 01176 } 01177 for (ib=0; ib<nb250; ib++) 01178 m[ib250[ib]] = floatbuf[pb250[ib]]; 01179 01180 if (getHDFattr(sd_id[2],"reflectance_offsets","EV_250_Aggr1km_RefSB",(VOIDP)floatbuf) != 0) { 01181 printf("-E- %s line %d: Error reading reflectance offset attribute.\n", 01182 __FILE__,__LINE__); 01183 return(1); 01184 } 01185 for (ib=0; ib<nb250; ib++) 01186 b[ib250[ib]] = floatbuf[pb250[ib]]; 01187 01188 if (getHDFattr(sd_id[2],"reflectance_scales","EV_500_Aggr1km_RefSB",(VOIDP)floatbuf) != 0) { 01189 printf("-E- %s line %d: Error reading reflectance scale attribute.\n", 01190 __FILE__,__LINE__); 01191 return(1); 01192 } 01193 for (ib=0; ib<nb500; ib++) 01194 m[ib500[ib]] = floatbuf[pb500[ib]]; 01195 01196 if (getHDFattr(sd_id[2],"reflectance_offsets","EV_500_Aggr1km_RefSB",(VOIDP)floatbuf) != 0) { 01197 printf("-E- %s line %d: Error reading reflectance offset attribute.\n", 01198 __FILE__,__LINE__); 01199 return(1); 01200 } 01201 for (ib=0; ib<nb500; ib++) 01202 b[ib500[ib]] = floatbuf[pb500[ib]]; 01203 01204 if (getHDFattr(sd_id[2],"reflectance_scales","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 01205 printf("-E- %s line %d: Error reading reflectance scale attribute.\n", 01206 __FILE__,__LINE__); 01207 return(1); 01208 } 01209 for (ib=0; ib<nb1000; ib++) 01210 m[ib1000[ib]] = floatbuf[pb1000[ib]]; 01211 01212 if (getHDFattr(sd_id[2],"reflectance_offsets","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 01213 printf("-E- %s line %d: Error reading reflectance offset attribute.\n", 01214 __FILE__,__LINE__); 01215 return(1); 01216 } 01217 for (ib=0; ib<nb1000; ib++) 01218 b[ib1000[ib]] = floatbuf[pb1000[ib]]; 01219 01220 /* radiance scale and offset factors (for destriping subframes) */ 01221 01222 if (getHDFattr(sd_id[2],"radiance_scales","EV_250_Aggr1km_RefSB",(VOIDP)floatbuf) != 0) { 01223 printf("-E- %s line %d: Error reading radiance scale attribute.\n", 01224 __FILE__,__LINE__); 01225 return(1); 01226 } 01227 for (ib=0; ib<nb250; ib++) 01228 mrad[ib250[ib]] = floatbuf[pb250[ib]]; 01229 01230 if (getHDFattr(sd_id[2],"radiance_offsets","EV_250_Aggr1km_RefSB",(VOIDP)floatbuf) != 0) { 01231 printf("-E- %s line %d: Error reading radiance offset attribute.\n", 01232 __FILE__,__LINE__); 01233 return(1); 01234 } 01235 for (ib=0; ib<nb250; ib++) 01236 brad[ib250[ib]] = floatbuf[pb250[ib]]; 01237 01238 if (getHDFattr(sd_id[2],"radiance_scales","EV_500_Aggr1km_RefSB",(VOIDP)floatbuf) != 0) { 01239 printf("-E- %s line %d: Error reading radiance scale attribute.\n", 01240 __FILE__,__LINE__); 01241 return(1); 01242 } 01243 for (ib=0; ib<nb500; ib++) 01244 mrad[ib500[ib]] = floatbuf[pb500[ib]]; 01245 01246 if (getHDFattr(sd_id[2],"radiance_offsets","EV_500_Aggr1km_RefSB",(VOIDP)floatbuf) != 0) { 01247 printf("-E- %s line %d: Error reading radiance offset attribute.\n", 01248 __FILE__,__LINE__); 01249 return(1); 01250 } 01251 for (ib=0; ib<nb500; ib++) 01252 brad[ib500[ib]] = floatbuf[pb500[ib]]; 01253 01254 if (getHDFattr(sd_id[2],"radiance_scales","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 01255 printf("-E- %s line %d: Error reading radiance scale attribute.\n", 01256 __FILE__,__LINE__); 01257 return(1); 01258 } 01259 for (ib=0; ib<nb1000; ib++) 01260 mrad[ib1000[ib]] = floatbuf[pb1000[ib]]; 01261 01262 if (getHDFattr(sd_id[2],"radiance_offsets","EV_1KM_RefSB",(VOIDP)floatbuf) != 0) { 01263 printf("-E- %s line %d: Error reading radiance offset attribute.\n", 01264 __FILE__,__LINE__); 01265 return(1); 01266 } 01267 for (ib=0; ib<nb1000; ib++) 01268 brad[ib1000[ib]] = floatbuf[pb1000[ib]]; 01269 } 01270 01271 /* radiance scale and offset factors for emmissive bands */ 01272 01273 if (getHDFattr(sd_id[2],"radiance_scales","EV_1KM_Emissive",(VOIDP)floatbuf) != 0) { 01274 printf("-E- %s line %d: Error reading emissive radiance scale attribute.\n", 01275 __FILE__,__LINE__); 01276 return(1); 01277 } 01278 for (ib=0; ib<nbir; ib++) 01279 mrad[ibir[ib]] = floatbuf[pbir[ib]]; 01280 01281 if (getHDFattr(sd_id[2],"radiance_offsets","EV_1KM_Emissive",(VOIDP)floatbuf) != 0) { 01282 printf("-E- %s line %d: Error reading emissive radiance offset attribute.\n", 01283 __FILE__,__LINE__); 01284 return(1); 01285 } 01286 for (ib=0; ib<nbir; ib++) 01287 brad[ibir[ib]] = floatbuf[pbir[ib]]; 01288 01289 /* rreflectance scale and offset factors for cirrus band */ 01290 01291 if (getHDFattr(sd_id[2],"reflectance_scales","EV_Band26",(VOIDP)&m_cirrus) != 0) { 01292 printf("-E- %s line %d: Error reading reflectance scale attribute.\n", 01293 __FILE__,__LINE__); 01294 return(1); 01295 } 01296 01297 if (getHDFattr(sd_id[2],"reflectance_offsets","EV_Band26",(VOIDP)&b_cirrus) != 0) { 01298 printf("-E- %s line %d: Error reading reflectance offset attribute.\n", 01299 __FILE__,__LINE__); 01300 return(1); 01301 } 01302 } 01303 01304 for (ip=0; ip<npix; ip++) 01305 l1rec->pixnum[ip] = spix + ip; 01306 01307 if (frame != lastframe) { 01308 01309 /* Get mirror side (0=A,1=B) and mirror normal (in ECR) */ 01310 READ_SDS_ID(sd_id_g,"Mirror side",&mside, frame,0,0,0, 1,1,1,1); 01311 READ_SDS_ID(sd_id_g,"T_inst2ECR" ,&mnorm, frame,0,0,0, 1,3,1,1); 01312 01313 /* Fix mirror-side if set to -1 in geolocation (used as array index) */ 01314 if (mside < 0) { 01315 mside = 0; 01316 badmside = 1; 01317 } else 01318 badmside = 0; 01319 01320 /* Read granule time and compute scan times */ 01321 READ_SDS_ID(sd_id_g,"EV start time",(float64 *)&TAIsec, frame,0,0,0, 1,1,1,1); 01322 01323 /* Work-around for single-scan time errors */ 01324 if (TAIsec < 0.0) { 01325 printf("-W- %s: bad time in geolocation file at frame %d, using previous.\n", 01326 __FILE__,frame); 01327 TAIsec = lastTAIsec; 01328 } else 01329 lastTAIsec = TAIsec; 01330 01331 usec = TAIsec + 725846400.0; 01332 unix2yds(usec,&year,&day,&dsec); 01333 01334 *(l1rec->year) = (int32_t)year; 01335 *(l1rec->day) = (int32_t)day; 01336 *(l1rec->msec) = (int32_t)(dsec*1000.0); 01337 01338 /* Get earth-sun distance correction for this frame */ 01339 esdist = esdist_(l1rec->year,l1rec->day,l1rec->msec); 01340 fsol = pow(1.0/esdist,2); 01341 01342 /* Load 1km geolocation buffer */ 01343 for (idet=0; idet<nd1000; idet++) { 01344 iline = frame * nd1000 + idet; 01345 READ_SDS_ID(sd_id_g,"Longitude",&lon1000[idet][0], iline,sp1000,0,0, 1,np1000,1,1); 01346 READ_SDS_ID(sd_id_g,"Latitude" ,&lat1000[idet][0], iline,sp1000,0,0, 1,np1000,1,1); 01347 if(!lonlat) { 01348 READ_SDS_ID(sd_id_g,"Height" ,hgt, iline,sp1000,0,0, 1,np1000,1,1); 01349 READ_SDS_ID(sd_id_g,"SolarZenith" ,solz, iline,sp1000,0,0, 1,np1000,1,1); 01350 READ_SDS_ID(sd_id_g,"SolarAzimuth" ,sola, iline,sp1000,0,0, 1,np1000,1,1); 01351 READ_SDS_ID(sd_id_g,"SensorZenith" ,senz, iline,sp1000,0,0, 1,np1000,1,1); 01352 READ_SDS_ID(sd_id_g,"SensorAzimuth" ,sena, iline,sp1000,0,0, 1,np1000,1,1); 01353 for (ip=0; ip<np1000; ip++) { 01354 hgt1000[idet][ip] = (float)hgt[ip]; 01355 solz1000[idet][ip] = (float)solz[ip] * 0.01; 01356 sola1000[idet][ip] = (float)sola[ip] * 0.01; 01357 senz1000[idet][ip] = (float)senz[ip] * 0.01; 01358 sena1000[idet][ip] = (float)sena[ip] * 0.01; 01359 } 01360 } // if not lonlat 01361 } 01362 01363 /* Load data buffers, if required */ 01364 01365 if(!lonlat) { 01366 if (resolution < 1000) { 01367 for (iw=0; iw<nb1000; iw++) { 01368 ib = ib1000[iw]; 01369 for (idet=0; idet<nd1000; idet++) { 01370 iline = frame * nd1000 + idet; 01371 READ_SDS_ID(sd_id[2],"EV_1KM_RefSB",&idata[ib*npix], pb1000[iw],iline,sp1000,0,1,1,np1000,1); 01372 for (ip=0; ip<np1000; ip++) 01373 fdatabuf[ib][idet][ip] = (float) idata[ib*npix+ip]; 01374 } 01375 } 01376 for (iw=0; iw<nbir; iw++) { 01377 ib = ibir[iw]; 01378 for (idet=0; idet<nd1000; idet++) { 01379 iline = frame * nd1000 + idet; 01380 READ_SDS_ID(sd_id[2],"EV_1KM_Emissive",&idata[ib*npix], pbir[iw],iline,sp1000,0,1,1,np1000,1); 01381 for (ip=0; ip<np1000; ip++) 01382 fdatabuf[ib][idet][ip] = (float) idata[ib*npix+ip]; 01383 } 01384 } 01385 for (iw=0; iw<nbcir; iw++) { 01386 ib = ibcir[iw]; 01387 for (idet=0; idet<nd1000; idet++) { 01388 iline = frame * nd1000 + idet; 01389 READ_SDS_ID(sd_id[2],"EV_Band26",&idata[ib*npix], iline,sp1000,0,0,1,np1000,1,1); 01390 for (ip=0; ip<np1000; ip++) 01391 fdatabuf[ib][idet][ip] = (float) idata[ib*npix+ip]; 01392 } 01393 } 01394 } 01395 01396 if (resolution < 500) { 01397 for (iw=0; iw<nb500; iw++) { 01398 ib = ib500[iw]; 01399 for (idet=0; idet<nd500; idet++) { 01400 iline = frame * nd500 + idet; 01401 READ_SDS_ID(sd_id[1],"EV_500_RefSB",&idata[ib*npix], pb500[iw],iline,sp500,0,1,1,np500,1); 01402 for (ip=0; ip<np500; ip++) 01403 fdatabuf[ib][idet][ip] = (float) idata[ib*npix+ip]; 01404 if (destripe) { 01405 subframe_calibration(l1rec->sensorID,ib,mrad[ib],brad[ib],sp500,np500,&fdatabuf[ib][idet][0]); 01406 } 01407 } 01408 } 01409 } 01410 } // if not lonlat 01411 01412 lastframe = frame; 01413 } 01414 01415 *(l1rec->year) = (int32_t)year; 01416 *(l1rec->day) = (int32_t)day; 01417 *(l1rec->msec) = (int32_t)(dsec*1000.0); 01418 01419 01420 /* Get position and path geometry */ 01421 01422 switch (resolution) { 01423 case 250: 01424 case 500: 01425 for (ip=0; ip<npix; ip++) 01426 modis_geo_interp(ip,detnum,resolution, 01427 &l1rec->lon [ip], &l1rec->lat [ip], &l1rec->height[ip], 01428 &l1rec->solz[ip], &l1rec->sola[ip], 01429 &l1rec->senz[ip], &l1rec->sena[ip], lonlat); 01430 break; 01431 default: 01432 memcpy(l1rec->lon, &lon1000 [detnum][0],npix*sizeof(float32)); 01433 memcpy(l1rec->lat, &lat1000 [detnum][0],npix*sizeof(float32)); 01434 if(!lonlat) { 01435 memcpy(l1rec->height,&hgt1000 [detnum][0],npix*sizeof(float32)); 01436 memcpy(l1rec->solz, &solz1000[detnum][0],npix*sizeof(float32)); 01437 memcpy(l1rec->sola, &sola1000[detnum][0],npix*sizeof(float32)); 01438 memcpy(l1rec->senz, &senz1000[detnum][0],npix*sizeof(float32)); 01439 memcpy(l1rec->sena, &sena1000[detnum][0],npix*sizeof(float32)); 01440 } // if not lonlat 01441 break; 01442 } 01443 01444 if(lonlat) { 01445 return(LIFE_IS_GOOD); 01446 } 01447 01448 01449 /* Check for nav errors */ 01450 01451 for (ip=0; ip<npix; ip++) { 01452 if (l1rec->lon[ip] < -181.0 || l1rec->lon[ip] > 181.0 || 01453 l1rec->lat[ip] < -91.0 || l1rec->lat[ip] > 91.0 || 01454 badmside) 01455 l1rec->navfail[ip] = 1; 01456 } 01457 01458 /* Compute polarization frame rotation angles */ 01459 compute_alpha(l1rec->lon, l1rec->lat, l1rec->senz, l1rec->sena, 01460 mnorm, npix, l1rec->alpha); 01461 01462 01463 /* Read L1B VNIR data, scale to radiance, and copy relevant bands to L1 record */ 01464 /* MODIS format is band sequential. We want band interlaced. We use imbedded */ 01465 /* scale factors to get reflectance, then convert to radiance through mean F0. */ 01466 /* Note that use of MCST radiance scale factor gives different result. */ 01467 01468 switch (resolution) { 01469 case 250: 01470 for (iw=0; iw<nb250; iw++) { 01471 ib = ib250[iw]; 01472 READ_SDS_ID(sd_id[0],"EV_250_RefSB",&idata[ib*npix], pb250[iw],scan,spix,0,1,1,npix,1); 01473 for (ip=0; ip<npix; ip++) 01474 fdata[ib*npix+ip] = (float) idata[ib*npix+ip]; 01475 if (destripe) { 01476 subframe_calibration(l1rec->sensorID,ib,mrad[ib],brad[ib],spix,npix,&fdata[ib*npix]); 01477 } 01478 } 01479 for (iw=0; iw<nb500; iw++) { 01480 ib = ib500[iw]; 01481 for (ip=0; ip<npix; ip++) { 01482 modis_data_interp(ib,ip,detnum,500,250,&fdata[ib*npix+ip]); 01483 idata[ib*npix+ip] = (uint16) rint(fdata[ib*npix+ip]); 01484 } 01485 } 01486 for (iw=0; iw<nb1000; iw++) { 01487 ib = ib1000[iw]; 01488 for (ip=0; ip<npix; ip++) { 01489 modis_data_interp(ib,ip,detnum,1000,250,&fdata[ib*npix+ip]); 01490 idata[ib*npix+ip] = (uint16) rint(fdata[ib*npix+ip]); 01491 } 01492 } 01493 for (iw=0; iw<nbir; iw++) { 01494 ib = ibir[iw]; 01495 for (ip=0; ip<npix; ip++) { 01496 modis_data_interp(ib,ip,detnum,1000,250,&fdata[ib*npix+ip]); 01497 idata[ib*npix+ip] = (uint16) rint(fdata[ib*npix+ip]); 01498 } 01499 } 01500 for (iw=0; iw<nbcir; iw++) { 01501 ib = ibcir[iw]; 01502 for (ip=0; ip<npix; ip++) { 01503 modis_data_interp(ib,ip,detnum,1000,250,&fdata[ib*npix+ip]); 01504 idata[ib*npix+ip] = (uint16) rint(fdata[ib*npix+ip]); 01505 } 01506 } 01507 break; 01508 01509 case 500: 01510 for (iw=0; iw<nb250; iw++) { 01511 ib = ib250[iw]; 01512 READ_SDS_ID(sd_id[1],"EV_250_Aggr500_RefSB",&idata[ib*npix], pb250[iw],scan,spix,0,1,1,npix,1); 01513 for (ip=0; ip<npix; ip++) 01514 fdata[ib*npix+ip] = (float) idata[ib*npix+ip]; 01515 } 01516 for (iw=0; iw<nb500; iw++) { 01517 ib = ib500[iw]; 01518 READ_SDS_ID(sd_id[1],"EV_500_RefSB",&idata[ib*npix], pb500[iw],scan,spix,0,1,1,npix,1); 01519 for (ip=0; ip<npix; ip++) 01520 fdata[ib*npix+ip] = (float) idata[ib*npix+ip]; 01521 if (destripe) { 01522 subframe_calibration(l1rec->sensorID,ib,mrad[ib],brad[ib],spix,npix,&fdata[ib*npix]); 01523 } 01524 } 01525 for (iw=0; iw<nb1000; iw++) { 01526 ib = ib1000[iw]; 01527 for (ip=0; ip<npix; ip++) { 01528 modis_data_interp(ib,ip,detnum,1000,500,&fdata[ib*npix+ip]); 01529 idata[ib*npix+ip] = (uint16) rint(fdata[ib*npix+ip]); 01530 } 01531 } 01532 for (iw=0; iw<nbir; iw++) { 01533 ib = ibir[iw]; 01534 for (ip=0; ip<npix; ip++) { 01535 modis_data_interp(ib,ip,detnum,1000,500,&fdata[ib*npix+ip]); 01536 idata[ib*npix+ip] = (uint16) rint(fdata[ib*npix+ip]); 01537 } 01538 } 01539 for (iw=0; iw<nbcir; iw++) { 01540 ib = ibcir[iw]; 01541 for (ip=0; ip<npix; ip++) { 01542 modis_data_interp(ib,ip,detnum,1000,500,&fdata[ib*npix+ip]); 01543 idata[ib*npix+ip] = (uint16) rint(fdata[ib*npix+ip]); 01544 } 01545 } 01546 break; 01547 01548 default: 01549 for (iw=0; iw<nb250; iw++) { 01550 ib = ib250[iw]; 01551 READ_SDS_ID(sd_id[2],"EV_250_Aggr1km_RefSB",&idata[ib*npix], pb250[iw],scan,spix,0,1,1,npix,1); 01552 for (ip=0; ip<npix; ip++) 01553 fdata[ib*npix+ip] = (float) idata[ib*npix+ip]; 01554 } 01555 for (iw=0; iw<nb500; iw++) { 01556 ib = ib500[iw]; 01557 READ_SDS_ID(sd_id[2],"EV_500_Aggr1km_RefSB",&idata[ib*npix], pb500[iw],scan,spix,0,1,1,npix,1); 01558 for (ip=0; ip<npix; ip++) 01559 fdata[ib*npix+ip] = (float) idata[ib*npix+ip]; 01560 } 01561 for (iw=0; iw<nb1000; iw++) { 01562 ib = ib1000[iw]; 01563 READ_SDS_ID(sd_id[2],"EV_1KM_RefSB",&idata[ib*npix], pb1000[iw],scan,spix,0,1,1,npix,1); 01564 for (ip=0; ip<npix; ip++) 01565 fdata[ib*npix+ip] = (float) idata[ib*npix+ip]; 01566 } 01567 for (iw=0; iw<nbir; iw++) { 01568 ib = ibir[iw]; 01569 READ_SDS_ID(sd_id[2],"EV_1KM_Emissive",&idata[ib*npix], pbir[iw],scan,spix,0,1,1,npix,1); 01570 for (ip=0; ip<npix; ip++) 01571 fdata[ib*npix+ip] = (float) idata[ib*npix+ip]; 01572 } 01573 for (iw=0; iw<nbcir; iw++) { 01574 ib = ibcir[iw]; 01575 READ_SDS_ID(sd_id[2],"EV_Band26",&idata[ib*npix], scan,spix,0,0, 1,npix,1,1); 01576 for (ip=0; ip<npix; ip++) 01577 fdata[ib*npix+ip] = (float) idata[ib*npix+ip]; 01578 } 01579 break; 01580 } 01581 01582 for (iw=0; iw<nbands; iw++) { 01583 01584 l1rec->Fo[iw] = Fobar[iw] * fsol; 01585 01586 for (ip=0; ip<npix; ip++) { 01587 01588 ipb = ip*NBANDS+iw; 01589 l1rec->Lt[ipb] = BAD_FLT; 01590 01591 /* skip night data for visible bands */ 01592 if (l1rec->solz[ip] >= SOLZNIGHT) 01593 continue; 01594 01595 /* check for sentinel values and flag as appropriate */ 01596 if (idata[iw*npix+ip] >= 65500) { 01597 switch (idata[iw*npix+ip]) { 01598 case 65527: /* not in earth view */ 01599 l1rec->navfail[ip] = 1; 01600 break; 01601 case 65529: /* saturation */ 01602 case 65533: 01603 l1rec->Lt[ipb] = 1000.0; 01604 if (iw <= 12) l1rec->hilt[ip] = 1; 01605 default: /* other stuff */ 01606 break; 01607 } 01608 01609 } else 01610 l1rec->Lt[ipb] = (fdata[iw*npix+ip] - b[iw]) * m[iw] * l1rec->Fo[iw]/pi; 01611 } 01612 } 01613 01614 for (iw=nbands; iw<nbands+nbir; iw++) { 01615 01616 ib = iw-nbands; 01617 01618 switch (resolution) { 01619 case 250: idet=detnum/4; break; 01620 case 500: idet=detnum/2; break; 01621 case 1000: idet=detnum/1; break; 01622 } 01623 01624 for (ip=0; ip<npix; ip++) { 01625 01626 ipb = ip*NBANDSIR+ib; 01627 l1rec->Ltir[ipb] = 0.0; 01628 01629 if (idata[iw*npix+ip] >= 65500) { 01630 ; 01631 } else { 01632 01633 l1rec->Ltir[ipb] = (fdata[iw*npix+ip] - brad[iw]) * mrad[iw]/10.0; 01634 01635 01636 if (l1rec->sensorID == MODIST || l1rec->sensorID == HMODIST) { 01637 l1rec->Ltir[ipb] *= (1.0-radoff[ib][idet]); 01638 if (mside == 1) l1rec->Ltir[ipb] *= (1.0 + mfact[ib]); 01639 } 01640 if (((file->input->evalmask & SSTMODS) != 0) && 01641 (l1rec->sensorID == MODISA || l1rec->sensorID == HMODISA)) { 01642 if (iw == 0 && detnum == 2-1) { 01643 l1rec->Ltir[ipb] *= ch20d2cor; 01644 } else if (iw == 0 && detnum == 3-1) { 01645 l1rec->Ltir[ipb] *= ch20d3cor; 01646 } else if (iw == 0 && detnum == 10-1) { 01647 l1rec->Ltir[ipb] *= ch20d10cor; 01648 } else if (iw == 1 && detnum == 1-1) { 01649 l1rec->Ltir[ipb] *= ch22d1cor; 01650 } else if (iw == 1 && detnum == 8-1) { 01651 l1rec->Ltir[ipb] *= ch22d8cor; 01652 } else if (iw == 1 && detnum == 9-1) { 01653 l1rec->Ltir[ipb] *= ch22d9cor; 01654 } else if (iw == 1 && detnum == 10-1) { 01655 l1rec->Ltir[ipb] *= ch22d10cor; 01656 } else if (iw == 2 && detnum == 1-1) { 01657 l1rec->Ltir[ipb] *= ch23d1cor; 01658 } 01659 } 01660 } 01661 } 01662 } 01663 01664 /* Read cirrus band */ 01665 01666 for (ip=0; ip<npix; ip++) { 01667 ib = ibcir[0]; 01668 l1rec->rho_cirrus[ip] = 0.0; 01669 if (idata[iw*npix+ip] < 65500) { 01670 l1rec->rho_cirrus[ip] = (fdata[ib*npix+ip] - b_cirrus) * m_cirrus / cos(l1rec->solz[ip]/RADEG); 01671 } 01672 } 01673 01674 firstScan = 0; 01675 01676 l1rec->sensorID = file->sensorID; 01677 l1rec->npix = file->npix; 01678 l1rec->detnum = (int32_t)detnum; 01679 l1rec->mside = (int32_t)mside; 01680 01681 /* Convert IR bands to brightness temperature */ 01682 radiance2bt(l1rec,resolution); 01683 01684 return(LIFE_IS_GOOD); 01685 } 01686 01687 01688 int closel1_hmodis_hdf(filehandle *file) 01689 { 01690 int i; 01691 01692 if ((sd_id_g != sd_id[0]) && SDend(sd_id_g)) { 01693 fprintf(stderr,"-E- %s line %d: SDend(%d) failed for file, %s.\n", 01694 __FILE__,__LINE__,sd_id[0],file->geofile); 01695 return(HDF_FUNCTION_ERROR); 01696 } 01697 for (i=0; i<3; i++) 01698 if (sd_id[i] > 0 && SDend(sd_id[i])) { 01699 fprintf(stderr,"-E- %s line %d: SDend(%d) failed for file, %s.\n", 01700 __FILE__,__LINE__,i,file->name); 01701 return(HDF_FUNCTION_ERROR); 01702 } 01703 01704 return(LIFE_IS_GOOD); 01705 }
1.7.6.1