/* * Copyright(c) 2015, Space Science and Engineering Center, UW-Madison * Refer to "McIDAS Software Acquisition and Distribution Policies" * in the file mcidas/data/license.txt */ /**** $Id: wariaget.c,v 1.22 2015/06/18 18:12:43 russd Tst $ ****/ /* ** Name: ** wariaget : ADDE data server for Himawari-8 image format. ** ** Interface: ** int ** main(int argc, char *argv[]) ** ** Input: ** argc : argument count ** argv : argument vector ** ** Input and Output: ** none ** ** Output: ** none ** ** Return values: ** 0 : success ** < 0 : failure, enumerated by text message returned to client ** ** Remarks: ** This program is reads Himawari-8 segmented data format ** ** Categories: ** ADDE server */ #include #include #include #if defined(__APPLE__) #include #endif #include #include "mcidas.h" #include "mcidasp.h" #include "mcncdf.h" #include "m0arg.h" #include "netcdf.h" #include "AreaDir.h" #include "SDIUtil.h" #include "m0glue.h" #include "math.h" #include "hisd.h" int **malloc2dint(int, int); int main ( int argc, char *argv[] ) /* CONSTANTS */ #define TRUE 1 /* TRUE flag */ #define FALSE 0 /* FALSE flag */ #define DIR_SIZE 64 /* size of the diectory block */ #define IMAGE_FILE_VERSION 4 /* AREA ID type element */ #define MAX_STR_LEN 80 /* card length */ #define WARI_MAX_BANDS 16 /* maximum number of bands supported */ #define WARI_MAX_TIMES 25000 /* size of file arrays */ #define WARI_MAX_SEGS 10 /* size of file arrays */ #define WARI_MASK 30 /* size of the Himawari no data region */ #define WARI_SSS 86 /* SSEC satellite number */ #define WARI_NAV_OFFSET 256 /* NAV block offset */ #define WARI_AUX_ENTRIES 0 /* AUX block entries */ #define SSEC_INVALID 0 /* SSEC INVALID_VALUE */ #define NORMAL_END 0 #define ERROR_ALLOCATE 72 #define ERROR_READ 79 #define ERROR_DATA 80 #define M_PI 3.14159265358979323846 #define DEGTORAD (M_PI/180.0) #define RADTODEG (180.0/M_PI) #define SCLUNIT 1.525878906250000e-05 #define NORMAL_END 0 #define ERROR_END 100 { /* internal Variables */ char llchar[4]={"LL "}; char substr[5]; /* segment substring buffer */ char level_map[WARI_MAX_BANDS]; /* prefix level map */ char *file_name; /* name of NetCDF file */ char *image_stype; /* image source type */ char *last_stype; /* previous source type */ char *dataset; /* dataset name */ char *def_unit; /* default UNIT */ char *dir_unit; /* directory UNIT */ char *comment; /* ADDE comment field */ char *format; /* ADDE dataset format */ char *group; /* dataset group */ char *info; /* ADDE dataset info */ char *mask; /* dataset file mask */ char *place; /* coordinate system */ char *punit; /* parameter UNITS */ char *request; /* request string */ char *type; /* ADDE dataset type */ char *unit; /* UNITS keyword */ char *cBand; /* file name BAND */ char *coverage; /* file name COVERAGE */ char *resolution; /* file name RESOLUTION */ char *segment; /* file name SEGMENT */ char **flist = NULL; /* hdf file list */ /* image tracking arrays */ int **Nsegs; int **lines; int **elems; int **bands; int **iress; int *temp1; int **temp2; int ***positions; const char *dum; /* dummy variable for arg routines */ const char server[] = {"WARIAGET"}; /* server name */ double request_center_latitude; /* requested center latitude */ double request_center_longitude; /* requested center longitude */ double dparm; /* generic double parameter */ float WARI_RAW_SCALE = 1.0; /* scaling for RAW data */ float WARI_RAD_SCALE = 1000.0; /* scaling for RAD data */ float WARI_TEMP_SCALE = 100.0; /* scaling for TEMP data */ float WARI_ALB_SCALE = 100.0; /* scaling for ALB data */ float WARI_BRIT_SCALE = 1.0; /* scaling for BRIT data */ float *float_array; /* generic float array */ float temp; /* temperature */ float lat_diff; /* latitude difference */ float lon_diff; /* longitude difference */ float minimum; /* minimum difference */ float image_geo_latitude; /* CP latitude */ float image_geo_longitude; /* CP longitude */ float image_geo_line_resolution; /* image LALO line resolution */ float image_geo_element_resolution; /* image LALO element resolution */ float sector_min_lat; /* request sector minimum latitude */ float sector_min_lon; /* request sector minimum longitude */ float sector_max_lat; /* request sector maximum latitude */ float sector_max_lon; /* request sector maximum longitude */ float image_line_resolution; /* image line resolution */ float image_element_resolution; /* image element resolution */ float fResolution; /* file resolution */ float save_sspLon; float save_cfac; float save_lfac; float save_coff; float save_loff; float rvalue; /* Calibration block */ Fint4 calb_BlockLen; Fint4 calb_bandNo; Fint4 calb_bitPix; Fint4 calb_errorCount; Fint4 calb_outCount; float calb_waveLen; float calb_gainCnt2rad; float calb_cnstCnt2rad; float calb_rad2btpC0; float calb_rad2btpC1; float calb_rad2btpC2; float calb_btp2radC0; float calb_btp2radC1; float calb_btp2radC2; float calb_lightSpeed; float calb_planckConst; float calb_bolzConst; float calb_rad2albedo; float calb_invalidValue; Fint4 aux_head[100]; /* AUX block header */ Fint4 directory[DIR_SIZE]; /* image directory block */ Fint4 nav_block[128]; /* navigation block */ Fint4 cal_block[128]; /* calibration block */ Fint4 cal_opt[2]; /* calibration options*/ Freal Flat; Freal Flon; Freal Fline; Freal Felem; Freal Fdummy; double Cal_Table[65536]; int file_handle; /* handle to parsed file name */ int length; /* length of parsed file name */ int num_parm; /* number of parameters from parsed file name */ int i_parm; /* file name parameter index */ int YMD, Y, M, D, HMS; /* Year,Month,Day varaibles */ int aux_block_size; /* size of the AUX block */ int bdate; /* geginning CCYYDDD */ int btime; /* beginning time */ int curDate; /* current CCYYDDD */ int curTime; /* current HHMMSS */ int dat_block_size; /* data block size */ int def_space; /* default SPACE */ int edate; /* ending range date */ int etime; /* ending time */ int hiBound; /* ending dataset position */ int i, ii; /* loop index */ int j, jj; /* loop index */ int k, kk; /* array index */ int dim3; int dim2; int xx; int indx; /* index pointer */ int lalo_array_size; /* lat/lon array size */ int lalo_array_line_size; /* lat/lon array line size */ int lalo_array_elem_size; /* lat/lon array elem size */ int line; /* image line number */ int loBound; /* starting dataset position */ int ncdf_id; /* NetCDF id number */ int dim_id; /* NetCDF dimenson number */ int var_id; /* NetCDF variable number */ int num_flist; /* number of files */ int ounit_flag; /* output unit flag */ int pscale; /* parameter scale */ int psize; /* prefix size */ int rc; /* function return code */ int rt_flag; /* ADDE realtime flag */ int space; /* output spacing */ int totalSize; /* transfer byte count */ int trace; /* trace flag */ int transaction = 0; /* transaction flag */ int user = 0; /* user id */ int nbad; /* number of bad lines */ int status; /* function status */ int band; /* band number */ int Pos; /* ADDE position number */ int bPos; /* Beginning ADDE position number */ int ePos; /* Ending ADDE position number */ int iPos; /* ADDE position increment */ int elem; /* image lement */ int image_stype_flag; /* image source type flag */ int image_hms; /* image time */ int image_cyd; /* image date CCYYDDD format */ int image_iyd; /* image date AREA format */ int image_sec; /* image time in seconds */ int image_Nsegs; /* max segments expected */ int image_Isegs; /* this segment number */ int image_data_size; /* size of a data element in bytes */ int image_band; /* file image band */ int adde_pos; /* number of times */ /* int dates[WARI_MAX_TIMES]; dates array */ /* int times[WARI_MAX_TIMES]; times array */ int seconds[WARI_MAX_TIMES]; /* times (second) array */ int cover[WARI_MAX_TIMES]; /* image coverage array */ /* int Nsegs[WARI_MAX_TIMES][WARI_MAX_BANDS]; total number of segments expected */ /* int lines[WARI_MAX_TIMES][WARI_MAX_BANDS]; image line dimension array */ /* int elems[WARI_MAX_TIMES][WARI_MAX_BANDS]; image element dimension array */ /* int bands[WARI_MAX_TIMES][WARI_MAX_BANDS]; image band array */ /* int iress[WARI_MAX_TIMES][WARI_MAX_BANDS]; image resolution array */ /* int positions[WARI_MAX_TIMES][WARI_MAX_BANDS][WARI_MAX_SEGS]; array of computed ADDE position numbers */ int image_n_bands; /* number of image bands */ int iBand; /* band from file name */ int iCoverage; /* coverage from file name */ int iSegment; /* segment from file name */ int iResolution; /* SSEC resolution */ int image_bands[WARI_MAX_BANDS]; /* array of image bands */ int var_ids[WARI_MAX_BANDS]; /* NetCDF variable id array */ int image_starting_line; /* image starting line */ int image_starting_element; /* image starting element */ int image_dataset_position; /* image position */ int ADDE_dataset_position; /* ADDE position */ int image_segment_position; /* image segment sub-position */ int image_Nline_size; /* segment image line size */ int image_Nseg_size; /* segment size */ int image_line_size; /* image line size */ int image_element_size; /* image element size */ int image_line_base_resolution; /* image base line resolution */ int image_element_base_resolution; /* image base element resolution */ int bRes; int image_geo_line_size; /* image LALO line resolution */ int image_geo_element_size; /* image LALO element resolution */ int image_geo_line_magnification; /* image LALO line magnification */ int image_geo_element_magnification; /* image LALO element magnification */ int image_geo_data_size; /* image LALO byte size */ int image_geo_latitude_resolution; /* image LALO latitude resolution */ int image_geo_longitude_resolution; /* image LALO longitude resolution */ int image_geo_nav_offset; /* LALO navigation element offset */ int image_lines, rows; /* image lines/rows */ int image_elems, cells; /* image elemnets/cells */ int CAL_BLOCK_SIZE = 128; /* CAL block word size */ int NAV_BLOCK_SIZE = 128; /* NAV block word size */ int WARI_CAL_OFFSET = 768; /* CAL block offset */ int WARI_AUX_OFFSET = 1280; /* AUX block offset */ int nav_offset; /* naviagtion element offset */ int request_data_size; /* request element byte size */ int request_n_bands; /* request number of bands */ int request_bands[WARI_MAX_BANDS]; /* request band array */ int NJ,NK; int request_lalo_flag; /* request LALO flag */ int request_image_flag; /* request image flag */ int request_area_flag; /* request area flag */ int request_line_size; /* request line size */ int request_element_size; /* request element size */ int request_line_resolution; /* request line resoluton */ int request_element_resolution; /* request element resolution */ int request_line_magnification; /* request line magnification */ int request_element_magnification; /* request element magnification */ int request_starting_line; /* request starting line */ int request_starting_element; /* request starting element */ int found_n_bands; /* image number of bands */ int found_bands[WARI_MAX_BANDS]; /* array of image bands */ int read_line_size; /* image read line size */ int read_element_size; /* image read element size */ int read_starting_line; /* image read startng line */ int read_starting_element; /* image read starting eleemnt */ int read_ending_line; /* image read ending line */ int read_ending_element; /* image read ending element */ int read_geo_line_size; /* image LALO line size */ int read_geo_element_size; /* image LALO element size */ int read_geo_starting_line; /* image LALO starting line */ int read_geo_starting_element; /* image LALO starting element */ int read_geo_ending_line; /* image LALO ending line */ int read_geo_ending_element; /* image LALO ending element */ int test; /* test for value */ int extend; /* test for multiple of 4 */ int size; /* size test variable */ int OUTSIZE; /* output line size */ int SAMPLE_size; /* input data sample size */ int LALO_size; /* input LALO sample size */ int Lat_size; /* input Latitude array size */ int Lon_size; /* input Longitude array size */ int DATA_size; /* input data array size */ int init; int TimeBaseSearch; int TimeBasePos; int ONE=1; int TWO=2; int LL; int currentSeg; int missingSeg; int seg_i; int thisSeg; int thisPos; int startPos; int inc_time; /* bump time needed to calculate target image start time */ servacct requestBlock; /* request block text */ short *SAMPLE_array; /* input sample array */ short *DATA_array; /* output sample array */ static char trace_string[500]; /* traceing text */ unsigned char *uchar_array; /* 1 byte output array */ unsigned int iexp; /* band map varaible */ Fint4 *ulong_array; /* 4 byte output array */ unsigned long int value; /* 4 byte output variable */ unsigned long long llparm; /* LongLong output variable */ unsigned short int *ushort_array; /* 2 byte output array */ /* string contains parsing mask */ McArgSyntax filesyntax = {"/_", "=", ";", "{\"", "}", "'", "'", NULL, "X", " "}; /* for JMA source code support */ HisdHeader *header; FILE *fp; double pixLon; double pixLat; float pixLin; float pixEle; double *Lon_array; double *Lat_array; double *Data_array; unsigned char *Vis_array; unsigned short int *Cnt_array; float Pix[5]; float Lin[5]; float Rii; float Rjj; double Rad; double Tbb; double Alb; unsigned short Cnt; int NCnt; int NCnt_valid; int Cnt_valid; int Cnt_error; double Cnt_min, Cnt_max; double Rad_min, Rad_max; double Alb_min, Alb_max; int Cnt_num; int Rad_num; int Alb_num; int CP_line; int CP_elem; /* Planck constant (Joule second) */ float PC = 6.6260755e-34; /* Speed of light in vacuum (meters per second) */ float SL = 2.9979246e+8; /* Boltzmann constant (Joules per Kelvin) */ float BC = 1.380658e-23; /* for JMA source code support */ /* ************************************************************************* Step 1: process the transaction request Use standard function calls to retrieve the parameters in the comm block. ************************************************************************* */ /* initialize local server */ rc = M0InitLocalServer(server, argc, argv, &requestBlock, &request); /* see if logging should be done, determined by hidden TRACE= flag */ trace = M0IsTraceSet(request); (void *)sprintf(trace_string,"%s TEST Himawari-8 Data Server V2.3 ",server); M0sxtrce( trace_string ); /* move the request back into the comm block */ (void) strncpy(requestBlock.text, request, sizeof(requestBlock.text)); /* fill in the data type */ transaction = atoi(argv[5]); (void) memcpy(requestBlock.transaction, &transaction, sizeof(requestBlock.transaction)); /* fill in the user */ user = atoi(argv[3]); (void) memcpy(requestBlock.user, &user, sizeof(requestBlock.user)); /* initialize accounting */ m0vserv_(&requestBlock); /* get dataset info */ rc = M0sxdatasetinfo(request, &group, &dataset, &type, &format, &mask, &info, &comment, &loBound, &hiBound, &rt_flag); /* Dataset Name */ sprintf(trace_string, "%s Dataset = %s", server,group); M0sxtrce( trace_string ); /* ************************************************************************* Step 2: Get a list of file names matching the ADDE file mask ************************************************************************* */ /* produce a list of files matching the file mask */ sprintf(trace_string, "%s File Mask = %s", server,mask); M0sxtrce( trace_string ); flist = M0GetMaskFileList(mask); if( flist == (char **)NULL ) { (void)strcpy(requestBlock.errormsg, "Failed to generate file mask list"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } else { num_flist = VecLen(flist); if( num_flist <= 0 ) { (void)strcpy(requestBlock.errormsg, "No files found matching file mask"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } sprintf(trace_string, "%s Number of files found = %d", server,num_flist); M0sxtrce( trace_string ); } /* malloc the 2D image tracking arrays */ Nsegs = malloc2dint( num_flist, WARI_MAX_BANDS ); lines = malloc2dint( num_flist, WARI_MAX_BANDS ); elems = malloc2dint( num_flist, WARI_MAX_BANDS ); bands = malloc2dint( num_flist, WARI_MAX_BANDS ); iress = malloc2dint( num_flist, WARI_MAX_BANDS ); /* malloc the 3D image tracking arrays */ dim3 = num_flist * WARI_MAX_BANDS * WARI_MAX_SEGS; temp1 = (int*) malloc(dim3 * sizeof(int)); dim2 = num_flist * WARI_MAX_BANDS; temp2 = (int**) malloc(dim2 * sizeof(int*)); positions = (int***) malloc( num_flist * sizeof(int**)); for( i=0; i< num_flist; i++) { positions[i] = temp2 + (WARI_MAX_BANDS*i); for( j=0; jEday TIME - image Btime<->Etime PLACE - request sector reference positioning A=AREA, I-IMAGE, E=EARTH U=UPPER, C= CENTER example: EC = EARTH CENTER LMAG - line magnification EMAG - element magnification SIZE - line and element request size BAND - image band STYPE - Source type (RAW, BRIT) SPACE - bytes per element (1,2 or 4) UNIT - output units (RAW, BRIT) ************************************************************************* */ /* This is where the POS, DATE and TIME parameter handelers used to be. They were moved because they are based on ADDE image position. The Himawari images are segemented which means there may be multiple files for each complete image so the files must be examined to determine which file(s) belong to a specific DATE and TIME. After that is done we will know the number of images we truely have */ /* PLACE */ rc = Mcargstr(0, "", 3, "AU", (const char **)&place); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid Coordinate format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } request_lalo_flag = FALSE; request_image_flag = FALSE; request_area_flag = FALSE; if( strncmp(place,"EC",2) == 0 || strncmp(place,"EU",2) == 0 ) request_lalo_flag = TRUE; if( strncmp(place,"IC",2) == 0 || strncmp(place,"IU",2) == 0 ) request_image_flag = TRUE; if( strncmp(place,"AC",2) == 0 || strncmp(place,"AU",2) == 0 ) request_area_flag = TRUE; /* Process the request: Latitude/Longitude*/ if( request_lalo_flag == TRUE ) { rc = Mcargdll(0, " ", 4, (double)999.0, (double)-90.0, (double)90.0, &request_center_latitude, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid LATitude format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } rc = Mcargdll(0, " ", 5, (double)999.0, (double)-180.0, (double)180.0, &request_center_longitude, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid LONgitude format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } request_center_longitude = -request_center_longitude; (void *)sprintf(trace_string,"%s Center Earth COORD: Lat=%f Lon=%f" ,server,request_center_latitude,request_center_longitude); M0sxtrce( trace_string ); } /* Image LINELE */ else if( request_image_flag == TRUE ) { rc = Mcargint(0, "", 4, 1, 1, 0, &request_starting_line, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid LINe format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } rc = Mcargint(0, "", 5, 1, 1, 0, &request_starting_element, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid ELEment format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } (void *)sprintf(trace_string,"%s Starting Image COORD: Lin=%d Ele=%d" ,server,request_starting_line,request_starting_element); M0sxtrce( trace_string ); } /* Area LINELE */ else if( request_area_flag == TRUE ) { rc = Mcargint(0, "", 4, 0, 1, 0, &read_starting_line, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid LINe format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } rc = Mcargint(0, "", 5, 0, 1, 0, &read_starting_element, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid ELEment format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } (void *)sprintf(trace_string,"%s Starting Area COORD: Lin=%d Ele=%d" ,server,read_starting_line,read_starting_element); M0sxtrce( trace_string ); } /* Process the request: Line Magnification */ rc = Mcargint(0, "LMA.G", 1, 1, -50, 1, &request_line_magnification, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid Line MAG specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } if( request_line_magnification == 0 ) request_line_magnification = 1; if( request_line_magnification < 0 ) request_line_magnification = -(request_line_magnification);; /* Process the request: Element Magnification */ rc = Mcargint(0, "EMA.G", 1, 1, -50, 1, &request_element_magnification, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid Element MAG specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } if( request_element_magnification == 0 ) request_element_magnification = 1; if( request_element_magnification < 0 ) request_element_magnification = -(request_element_magnification);; (void *)sprintf(trace_string,"%s MAG: Lin=%d Ele=%d" ,server,request_line_magnification,request_element_magnification); M0sxtrce( trace_string ); /* SIZE */ rc = Mcargint(0, "", 7, 480, 1, 99999, &request_line_size, &dum ); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid Line SIZE specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } rc = Mcargint(0, "", 8, 640, 1, 99999, &request_element_size, &dum ); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid Element SIZE specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } (void *)sprintf(trace_string,"%s SIZE: Lin=%d Ele=%d" ,server,request_line_size,request_element_size); M0sxtrce( trace_string ); /* BAND */ request_n_bands = Mcargnum(0, "BAN.D" ); if( request_n_bands == 0 ) { (void)strcpy(requestBlock.errormsg, "No BAND(s) specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } else { char *cALL; /* BAND=ALL*/ rc = Mcargstr(0, "BAN.D", 1, " ", (const char **)&cALL); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid BAND specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } /* Transfer ALL bands */ if( strncmp( cALL, "ALL",3 ) == 0 ) { request_n_bands = WARI_MAX_BANDS - 2 ; for( i=0; i 1 ) { (void)strcpy(requestBlock.errormsg, "Multi-band transfer is not permitted for this data type"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } /* STYPE */ image_stype_flag = Mcargnum(0, "STY.PE" ); if( image_stype_flag != 0 ) { char *cstype; /* STYPE keyword */ rc = Mcargstr(0, "STY.PE", 1, "VISR", (const char **)&cstype); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid STYPE specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } if( strncmp( cstype, "VISR",4 ) != 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid STYPE specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } if( request_n_bands > 1 ) { (void)strcpy(requestBlock.errormsg, "STYPE with Multi-band not allowed"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } image_stype = stralloc( "VISR", NULL); def_unit = stralloc( "BRIT", NULL); last_stype = stralloc( "WARI", NULL); def_space = 1; CAL_BLOCK_SIZE = 0; WARI_CAL_OFFSET = 0; WARI_AUX_OFFSET = 768; } else { image_stype = stralloc( "WARI", NULL); def_unit = stralloc( "RAW",NULL); last_stype = stralloc( " ", NULL); def_space = 2; } /* SPACE */ rc = Mcargint(0, "SPA.CE", 1, def_space, 0, 4, &space, &dum); if( rc < 0 || space == 3 ) { (void)strcpy(requestBlock.errormsg, "Invalid SPACE specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } /* UNIT */ rc = Mcargstr(0, "UNI.T", 1, def_unit, (const char **)&unit); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid UNIT specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } if( strncmp( unit, "RAW",3 ) == 0 ) { if( space == 0 ) space = 2; if( space < 2 ) { (void)strcpy(requestBlock.errormsg, "Invalid Spacing for RAW specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } dir_unit = stralloc( "RAW ",NULL); request_data_size = space; pscale = WARI_RAW_SCALE; punit = stralloc(" ",NULL); ounit_flag = 0; } else if( strncmp( unit, "RAD" ,3) == 0 ) { if( space == 0 ) space = 4; if( request_n_bands > 1 ) { (void)strcpy(requestBlock.errormsg, "Multi-band transfer UNITS must be RAW"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } dir_unit = stralloc( "RAD ",NULL); request_data_size = space; pscale = WARI_RAD_SCALE; punit = stralloc("WM**",NULL); ounit_flag = 1; } else if( strncmp( unit, "TEMP" ,4) == 0 ) { if( space == 0 ) space = 2; if( request_n_bands > 1 ) { (void)strcpy(requestBlock.errormsg, "Multi-band transfer UNITS must be RAW"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } dir_unit = stralloc( "TEMP",NULL); request_data_size = space; pscale = WARI_TEMP_SCALE; punit = stralloc("K ",NULL); ounit_flag = 2; } else if( strncmp( unit, "ALB" ,3) == 0 ) { if( space == 0 ) space = 2; if( request_n_bands > 1 ) { (void)strcpy(requestBlock.errormsg, "Multi-band transfer UNITS must be RAW"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } dir_unit = stralloc( "ALB ",NULL); request_data_size = space; pscale = WARI_ALB_SCALE; punit = stralloc("% ",NULL); ounit_flag = 3; } else if( strncmp( unit, "BRIT" ,4) == 0 ) { if( space == 0 ) space = 1; if( request_n_bands > 1 ) { (void)strcpy(requestBlock.errormsg, "Multi-band transfer UNITS must be RAW"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } dir_unit = stralloc( "BRIT",NULL); request_data_size = space; pscale = WARI_BRIT_SCALE; punit = stralloc(" ",NULL); ounit_flag = 4; } else { (void)strcpy(requestBlock.errormsg, "Invalid UNIT specified"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } (void *)sprintf(trace_string,"%s SPACE: %d",server,space); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"%s UNIT: %s",server,unit); M0sxtrce( trace_string ); /* ************************************************************************* Step 4: Analyze the file names Use the file name to get the image date and time. ************************************************************************* */ /* open each file to get the day and time of the data */ adde_pos = 0; for( Pos=0; Pos= 0 ) { Y = YMD / 10000; M = (YMD % 10000) / 100; D = (YMD % 100 ); rc = Mcdmytocyd( D,M,Y, &image_cyd); rc = Mccydtoiyd( image_cyd, &image_iyd); } rc = Mcargint( file_handle, " ", i_parm+3, 0,0,-1, &HMS, &dum); if( rc >= 0 ) { image_hms = HMS*100; status = Mcdaytimetosec( image_cyd, image_hms, &image_sec ); } /* sprintf(trace_string, "%s Date=%d", server, image_cyd); M0sxtrce( trace_string ); sprintf(trace_string, "%s Time=%d", server, image_hms); M0sxtrce( trace_string ); */ /* break the band from the file name */ rc = Mcargstr( file_handle, " ", i_parm+4, " ", &cBand); iBand = 0; if (strncmp(cBand, "B01", 3) == 0) iBand = 1; if (strncmp(cBand, "B02", 3) == 0) iBand = 2; if (strncmp(cBand, "B03", 3) == 0) iBand = 3; if (strncmp(cBand, "B04", 3) == 0) iBand = 4; if (strncmp(cBand, "B05", 3) == 0) iBand = 5; if (strncmp(cBand, "B06", 3) == 0) iBand = 6; if (strncmp(cBand, "B07", 3) == 0) iBand = 7; if (strncmp(cBand, "B08", 3) == 0) iBand = 8; if (strncmp(cBand, "B09", 3) == 0) iBand = 9; if (strncmp(cBand, "B10", 3) == 0) iBand = 10; if (strncmp(cBand, "B11", 3) == 0) iBand = 11; if (strncmp(cBand, "B12", 3) == 0) iBand = 12; if (strncmp(cBand, "B13", 3) == 0) iBand = 13; if (strncmp(cBand, "B14", 3) == 0) iBand = 14; if (strncmp(cBand, "B15", 3) == 0) iBand = 15; if (strncmp(cBand, "B16", 3) == 0) iBand = 16; if( iBand == 0 ) { sprintf(trace_string, "%s Could not read the Band Number", server); M0sxtrce( trace_string ); } /* sprintf(trace_string, "%s Band=%s iBand=%d", server,cBand,iBand); M0sxtrce( trace_string ); */ /* break the coverage from the file name */ rc = Mcargstr( file_handle, " ", i_parm+5, " ", &coverage); iCoverage = 0; if (strncmp(coverage, "FLDK", 4) == 0) iCoverage = 1; if (strncmp(coverage, "HNDK", 4) == 0) iCoverage = 2; if (strncmp(coverage, "HSDK", 4) == 0) iCoverage = 3; if (strncmp(coverage, "JP01", 4) == 0) iCoverage = 41; if (strncmp(coverage, "JP02", 4) == 0) iCoverage = 42; if (strncmp(coverage, "JP03", 4) == 0) iCoverage = 43; if (strncmp(coverage, "JP04", 4) == 0) iCoverage = 44; if (strncmp(coverage, "R301", 4) == 0) iCoverage = 51; if (strncmp(coverage, "R302", 4) == 0) iCoverage = 52; if (strncmp(coverage, "R303", 4) == 0) iCoverage = 53; if (strncmp(coverage, "R304", 4) == 0) iCoverage = 54; if (strncmp(coverage, "R401", 4) == 0) iCoverage = 601; if (strncmp(coverage, "R402", 4) == 0) iCoverage = 602; if (strncmp(coverage, "R403", 4) == 0) iCoverage = 603; if (strncmp(coverage, "R404", 4) == 0) iCoverage = 604; if (strncmp(coverage, "R405", 4) == 0) iCoverage = 605; if (strncmp(coverage, "R406", 4) == 0) iCoverage = 606; if (strncmp(coverage, "R407", 4) == 0) iCoverage = 607; if (strncmp(coverage, "R408", 4) == 0) iCoverage = 608; if (strncmp(coverage, "R409", 4) == 0) iCoverage = 609; if (strncmp(coverage, "R410", 4) == 0) iCoverage = 610; if (strncmp(coverage, "R411", 4) == 0) iCoverage = 611; if (strncmp(coverage, "R412", 4) == 0) iCoverage = 612; if (strncmp(coverage, "R413", 4) == 0) iCoverage = 613; if (strncmp(coverage, "R414", 4) == 0) iCoverage = 614; if (strncmp(coverage, "R415", 4) == 0) iCoverage = 615; if (strncmp(coverage, "R416", 4) == 0) iCoverage = 616; if (strncmp(coverage, "R417", 4) == 0) iCoverage = 617; if (strncmp(coverage, "R418", 4) == 0) iCoverage = 618; if (strncmp(coverage, "R419", 4) == 0) iCoverage = 619; if (strncmp(coverage, "R420", 4) == 0) iCoverage = 620; if (strncmp(coverage, "R501", 4) == 0) iCoverage = 701; if (strncmp(coverage, "R502", 4) == 0) iCoverage = 702; if (strncmp(coverage, "R503", 4) == 0) iCoverage = 703; if (strncmp(coverage, "R504", 4) == 0) iCoverage = 704; if (strncmp(coverage, "R505", 4) == 0) iCoverage = 705; if (strncmp(coverage, "R506", 4) == 0) iCoverage = 706; if (strncmp(coverage, "R507", 4) == 0) iCoverage = 707; if (strncmp(coverage, "R508", 4) == 0) iCoverage = 708; if (strncmp(coverage, "R509", 4) == 0) iCoverage = 709; if (strncmp(coverage, "R510", 4) == 0) iCoverage = 710; if (strncmp(coverage, "R511", 4) == 0) iCoverage = 711; if (strncmp(coverage, "R512", 4) == 0) iCoverage = 712; if (strncmp(coverage, "R513", 4) == 0) iCoverage = 713; if (strncmp(coverage, "R514", 4) == 0) iCoverage = 714; if (strncmp(coverage, "R515", 4) == 0) iCoverage = 715; if (strncmp(coverage, "R516", 4) == 0) iCoverage = 716; if (strncmp(coverage, "R517", 4) == 0) iCoverage = 717; if (strncmp(coverage, "R518", 4) == 0) iCoverage = 718; if (strncmp(coverage, "R519", 4) == 0) iCoverage = 719; if (strncmp(coverage, "R520", 4) == 0) iCoverage = 720; if( iCoverage == 0 ) { sprintf(trace_string, "%s Could not read the Coverge Code", server); M0sxtrce( trace_string ); } /* patch for Target sector image time */ if( iCoverage == 51 ) { sprintf(trace_string, "%s TARGET Coverage=%s", server,coverage); M0sxtrce( trace_string ); inc_time = 215; status = Mcinctime( image_cyd, image_hms, inc_time, &image_cyd, &image_hms ); status = Mcdaytimetosec( image_cyd, image_hms, &image_sec ); sprintf(trace_string, "%s NEW IMAGE TIME=%d", server,image_hms); M0sxtrce( trace_string ); } else if( iCoverage == 52 ) { sprintf(trace_string, "%s TARGET Coverage=%s", server,coverage); M0sxtrce( trace_string ); inc_time = 545; status = Mcinctime( image_cyd, image_hms, inc_time, &image_cyd, &image_hms ); status = Mcdaytimetosec( image_cyd, image_hms, &image_sec ); sprintf(trace_string, "%s NEW IMAGE TIME=%d", server,image_hms); M0sxtrce( trace_string ); } else if( iCoverage == 53 ) { sprintf(trace_string, "%s TARGET Coverage=%s", server,coverage); M0sxtrce( trace_string ); inc_time = 715; status = Mcinctime( image_cyd, image_hms, inc_time, &image_cyd, &image_hms ); status = Mcdaytimetosec( image_cyd, image_hms, &image_sec ); sprintf(trace_string, "%s NEW IMAGE TIME=%d", server,image_hms); M0sxtrce( trace_string ); } else if( iCoverage == 54 ) { sprintf(trace_string, "%s TARGET Coverage=%s", server,coverage); M0sxtrce( trace_string ); inc_time = 945; status = Mcinctime( image_cyd, image_hms, inc_time, &image_cyd, &image_hms ); status = Mcdaytimetosec( image_cyd, image_hms, &image_sec ); sprintf(trace_string, "%s NEW IMAGE TIME=%d", server,image_hms); M0sxtrce( trace_string ); } /* sprintf(trace_string, "%s Coverage=%s", server,coverage); M0sxtrce( trace_string ); */ /* break the resolution from the file name */ rc = Mcargstr( file_handle, " ", i_parm+6, " ", &resolution); fResolution = 0.0; iResolution = 0; if (strncmp(resolution, "R05", 3) == 0) { fResolution = 0.5; iResolution = 1; } else if (strncmp(resolution, "R10", 3) == 0) { fResolution = 1.0; iResolution = 2; } else if (strncmp(resolution, "R20", 3) == 0) { fResolution = 2.0; iResolution = 4; } else if (strncmp(resolution, "R40", 3) == 0) { fResolution = 4.0; iResolution = 8; } if( iResolution == 0 ) { sprintf(trace_string, "%s Could not read the SSP Resolution", server); M0sxtrce( trace_string ); } /* sprintf(trace_string, "%s iResolution=%d", server,iResolution); M0sxtrce( trace_string ); */ /* break the segment from the file name */ rc = Mcargstr( file_handle, " ", i_parm+7, " ", &segment); (void) strncpy(substr, segment+1, 4); substr[4] = '\0'; rc = Mcstrtoint(substr,&iSegment); image_Isegs = iSegment / 100; image_Nsegs = iSegment % 100; /* sprintf(trace_string, "%s Segment=%s iSegment=%d", server,segment,iSegment); M0sxtrce( trace_string ); sprintf(trace_string, "%s Isegs=%d Nsegs=%d", server,image_Isegs,image_Nsegs); M0sxtrce( trace_string ); */ /* line/elem dimesions based on coverage and band */ if( iCoverage >= 1 && iCoverage <= 3 ) { if( iBand == 1 ) { image_lines = 11000; image_elems = 11000; } else if( iBand == 2) { image_lines = 11000; image_elems = 11000; } else if( iBand == 3) { image_lines = 22000; image_elems = 22000; } else if( iBand == 4) { image_lines = 11000; image_elems = 11000; } else { image_lines = 5500; image_elems = 5500; } } else if( iCoverage >= 41 && iCoverage <= 44) { if( iBand == 1 ) { image_lines = 2400; image_elems = 3000; } else if( iBand == 2) { image_lines = 2400; image_elems = 3000; } else if( iBand == 3) { image_lines = 4800; image_elems = 6000; } else if( iBand == 4) { image_lines = 2400; image_elems = 3000; } else { image_lines = 1200; image_elems = 1500; } } else if( iCoverage >= 51 && iCoverage <= 54) { if( iBand == 1 ) { image_lines = 1000; image_elems = 1000; } else if( iBand == 2) { image_lines = 1000; image_elems = 1000; } else if( iBand == 3) { image_lines = 2000; image_elems = 2000; } else if( iBand == 4) { image_lines = 1000; image_elems = 1000; } else { image_lines = 500; image_elems = 500; } } else if( iCoverage >= 601 && iCoverage <= 620) { if( iBand == 1 ) { image_lines = 500; image_elems = 1000; } else if( iBand == 2) { image_lines = 500; image_elems = 1000; } else if( iBand == 3) { image_lines = 1000; image_elems = 2000; } else if( iBand == 4) { image_lines = 500; image_elems = 1000; } else { image_lines = 500; image_elems = 500; } } else if( iCoverage >= 701 && iCoverage <= 720) { if( iBand == 1 ) { image_lines = 500; image_elems = 1000; } else if( iBand == 2) { image_lines = 500; image_elems = 1000; } else if( iBand == 3) { image_lines = 1000; image_elems = 2000; } else if( iBand == 4) { image_lines = 500; image_elems = 1000; } else { image_lines = 250; image_elems = 500; } } /* sprintf(trace_string, "%s Lines=%d", server,image_lines); M0sxtrce( trace_string ); sprintf(trace_string, "%s Elems=%d", server,image_elems); M0sxtrce( trace_string ); */ /* SSEC image line/elem resolution factor based on band number sprintf(trace_string, "%s iBand =%d", server,iBand); M0sxtrce( trace_string ); if( iBand == 1 ) { image_line_resolution = 2; image_element_resolution = 2; image_geo_line_resolution = 1.0; image_geo_element_resolution = 1.0; } else if( iBand == 2 ) { image_line_resolution = 2; image_element_resolution = 2; image_geo_line_resolution = 1.0; image_geo_element_resolution = 1.0; } else if( iBand == 3 ) { image_line_resolution = 1; image_element_resolution = 1; image_geo_line_resolution = 0.5; image_geo_element_resolution = 0.5; } else if( iBand == 4 ) { image_line_resolution = 2; image_element_resolution = 2; image_geo_line_resolution = 1.0; image_geo_element_resolution = 1.0; } else { image_line_resolution = 4; image_element_resolution = 4; image_geo_line_resolution = 2.0; image_geo_element_resolution = 2.0; } sprintf(trace_string, "%s Image Line Res =%f", server,image_line_resolution); M0sxtrce( trace_string ); sprintf(trace_string, "%s Image Elem Res =%f", server,image_element_resolution); M0sxtrce( trace_string ); */ /* ************************************************************************* Step 6: Arrays of image date, time and dimensions into time order ************************************************************************* */ /* make a list of the files in ascending time order */ if( image_sec > 0 || image_lines > 0 || image_elems > 0 ) { /* first entry into date/time tables */ if( adde_pos == 0 ) { /* zero the file position array */ for( i=0; ii; j-- ) { seconds[j] = seconds[j-1]; /* dates[j] = dates[j-1]; */ /* times[j] = times[j-1]; */ cover[j] = cover[j-1]; for( ii=0; ii FALSE */ TimeBaseSearch = 0; /* DAY */ if( Mcargnum(0, "DAY") != 0 ) { /* set TIME Base switch to TRUE */ /* get the specified DATE */ rc = Mcargiyd(0, "DAY", 1, curDate, 1972001, curDate, &bdate, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid DAY format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } if( iPos <= 0 && rc != 300) { bPos = 1; ePos = adde_pos; (void *)sprintf(trace_string,"%s ********> DAY RANGE",server); M0sxtrce( trace_string ); TimeBaseSearch = 1; TimeBasePos = iPos; } rc = Mcargiyd(0, "DAY", 2, bdate, bdate, curDate, &edate, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid DAY format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } } else { bdate = 1972001; edate = curDate; } (void *)sprintf(trace_string,"%s DAY: range= %d to %d",server,bdate,edate); M0sxtrce( trace_string ); /* TIME */ if( Mcargnum(0, "TIM.E" ) != 0 ) { rc = Mcargihr(0, "TIM.E", 1, 0, 0, 235959, &btime, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid TIME format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } if( iPos <= 0 && rc != 400) { (void *)sprintf(trace_string,"%s ********> TIME RANGE",server); M0sxtrce( trace_string ); bPos = 1; ePos = adde_pos; TimeBaseSearch = 1; TimeBasePos = iPos; } rc = Mcargihr(0, "TIM.E", 2, btime, 0, 235959, &etime, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid TIME format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } if( rc == 400 ) etime = 235959; } else { btime = 0; etime = 235959; } (void *)sprintf(trace_string,"%s TIME: range= %d to %d",server,btime,etime); M0sxtrce( trace_string ); if( TimeBaseSearch == 1 ) { (void *)sprintf(trace_string,"%s ***** > TIME BASED SEARCH ",server); M0sxtrce( trace_string ); } /* Process the request: Dataset Position */ loBound = (adde_pos * -1) +1; hiBound = adde_pos; rc = Mcargint(0, "", 2, 0, loBound, hiBound, &image_dataset_position, &dum); if( rc < 0 ) { (void)strcpy(requestBlock.errormsg, "Invalid Dataset Position format"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } /* ************************************************************************* Step 7: MAIN LOOP At this point we have a list of files that have passed the date, time and dimension tests. The next step is to scan this list to verify that the dates and times meet any transaction conditions. If an image file passes the condition tests we read the file datasets to locate all of the "parts" we need to satisfy a data transfer. The necessary parts are: - navigation info - calibration info - Data array ************************************************************************* */ sprintf(trace_string, "%s BEFORE MAIN LOOP --- bPos=%d, ePos=%d", server, bPos, ePos); M0sxtrce( trace_string ); for( ADDE_dataset_position=bPos; ADDE_dataset_position<=ePos; ADDE_dataset_position++) { found_n_bands = 0; /* for Time Based Search invert the position */ if( TimeBaseSearch == 0 ) { image_dataset_position = ADDE_dataset_position; } else { image_dataset_position = adde_pos - (ADDE_dataset_position-1); } /* convert seconds to date and time */ status = Mcsectodaytime( seconds[image_dataset_position-1], &image_cyd, &image_hms ); /* convert julian date to image date */ /* image_cyd = dates[image_dataset_position-1]; */ status = Mccydtoiyd( image_cyd, &image_iyd ); /* Check client search conditions: DAY=bDAY eDAY */ if( (image_cyd < bdate) || (image_cyd > edate) ) continue; /* Check client search conditions: TIME=bTIME eTIME */ /* image_hms = times[image_dataset_position-1]; */ if( (image_hms < btime) || (image_hms > etime) ) continue; /* Check client search conditions: BAND=band */ image_band = bands[image_dataset_position-1][request_bands[0]-1]; if( image_band == -1 ) continue; if( image_band != request_bands[0] ) continue; /* Check for TimeBaseSearch */ if( TimeBaseSearch == 1 ) { if( TimeBasePos < 0 ) { TimeBasePos++; continue; } } sprintf(trace_string, "%s looking for band=%d", server, image_band); M0sxtrce( trace_string ); /* get the resolution of the image band */ image_line_base_resolution = iress[image_dataset_position-1][image_band-1]; image_element_base_resolution = iress[image_dataset_position-1][image_band-1]; sprintf(trace_string, "%s IMAGE_LINE_BASE_RES=%d",server, image_line_base_resolution); M0sxtrce( trace_string ); /* loop through the image segments */ Cnt_valid = 0; Cnt_error = 0; image_Nline_size = -9999; for( image_segment_position=0; image_segment_positionnav->sspLat; image_geo_longitude = header->proj->subLon; image_element_size = header->data->nPix; image_Nline_size = header->data->nLin; image_Nseg_size = header->seg->strLineNo; image_line_size = image_Nline_size * Nsegs[image_dataset_position-1][image_band-1]; image_data_size = header->data->bitPix/8; image_n_bands = 1; sprintf(trace_string, "%s Lines=%d Elems=%d",server, image_Nline_size, image_element_size); M0sxtrce( trace_string ); sprintf(trace_string, "%s First line number of the image = %d ",server, image_Nseg_size); M0sxtrce( trace_string ); sprintf(trace_string, "%s PROJ SSP: %f NADIR: %f",server, header->proj->subLon, header->nav->nadirLon); M0sxtrce( trace_string ); /* convert center lon/lat to elem/line */ if( request_lalo_flag == TRUE ) { pixLon = (double)request_center_longitude; pixLat = (double)request_center_latitude; lonlat_to_pixlin(header,pixLon,pixLat,&pixEle,&pixLin); CP_line = (int)(pixLin+0.5); CP_elem = (int)(pixEle+0.5); sprintf(trace_string, "%s CP Lat=%f Lon=%f Line=%d Elem=%d",server, pixLat, pixLon, CP_line, CP_elem); M0sxtrce( trace_string ); } /* image starting coordinates */ image_starting_line = 1; image_starting_element = 1; /* reserve memory for the image arrays */ SAMPLE_size = image_element_size; if( image_segment_position == 0 ) { sprintf(trace_string, "%s Malloc Data=%d Count=%d",server, SAMPLE_size, image_element_size); M0sxtrce( trace_string ); Data_array = malloc( image_element_size * sizeof(double) ); if( Data_array == NULL ) { (void)strcpy(requestBlock.errormsg, "FAILED -- malloc Data array"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } Cnt_array = malloc( image_element_size * sizeof(unsigned short) ); if( Cnt_array == NULL ) { (void)strcpy(requestBlock.errormsg, "FAILED -- malloc Cnt array"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } /* initialize the 1D array index */ kk = 0; /* initialize the min/max */ Cnt_min = 9999.0000; Cnt_max = -9999.0000; Cnt_num = 0.0000; Rad_min = 9999.0000; Rad_max = -9999.0000; Rad_num = 0.0000; Alb_min = 9999.0000; Alb_max = -9999.0000; Alb_num = 0.0000; /* naviagtion parameters */ /* save_sspLon = header->proj->subLon; */ save_sspLon = 140.70; save_cfac = header->proj->cfac; save_lfac = header->proj->lfac; save_coff = header->proj->coff; save_loff = header->proj->loff; /* calibration parameters */ calb_BlockLen = (int)header->calib->BlockLen; calb_bandNo = (int)header->calib->bandNo; calb_bitPix = (int)header->calib->bitPix; calb_errorCount = (int)header->calib->errorCount; calb_outCount = (int)header->calib->outCount; calb_invalidValue = SSEC_INVALID; calb_waveLen = (float)header->calib->waveLen; calb_gainCnt2rad = (float)header->calib->gain_cnt2rad; calb_cnstCnt2rad = (float)header->calib->cnst_cnt2rad; if( (header->calib->bandNo>=7 && strstr(header->basic->satName,"Himawari")!=NULL ) || (header->calib->bandNo>=2 && strstr(header->basic->satName,"MTSAT-2") !=NULL ) ) { calb_rad2btpC0 = (float)header->calib->rad2btp_c0; calb_rad2btpC1 = (float)header->calib->rad2btp_c1; calb_rad2btpC2 = (float)header->calib->rad2btp_c2; calb_btp2radC0 = (float)header->calib->btp2rad_c0; calb_btp2radC1 = (float)header->calib->btp2rad_c1; calb_btp2radC2 = (float)header->calib->btp2rad_c2; } else { calb_rad2albedo = (float)header->calib->rad2albedo; } sprintf(trace_string, " GAINCNT=%e",calb_gainCnt2rad); M0sxtrce( trace_string ); sprintf(trace_string, " CNSTCNT=%e",calb_cnstCnt2rad); M0sxtrce( trace_string ); sprintf(trace_string, " RAD2ALB=%e",calb_rad2albedo); M0sxtrce( trace_string ); sprintf(trace_string, " RAD2btpC0=%e",calb_rad2btpC0); M0sxtrce( trace_string ); sprintf(trace_string, " RAD2btpC1=%e",calb_rad2btpC1); M0sxtrce( trace_string ); sprintf(trace_string, " RAD2btpC2=%e",calb_rad2btpC2); M0sxtrce( trace_string ); sprintf(trace_string, " INVALID=%d",calb_invalidValue); M0sxtrce( trace_string ); sprintf(trace_string, " outCount=%d",header->calib->outCount); M0sxtrce( trace_string ); sprintf(trace_string, " errorCount=%d",header->calib->errorCount); M0sxtrce( trace_string ); /* build the Cal_table */ for( jj=0; jj<65536; jj++ ) { /* Calibration is RAW */ if( ounit_flag == 0 ) { if( jj == header->calib->outCount || jj == header->calib->errorCount ) { Cal_Table[jj] = INVALID_VALUE; } else { Cal_Table[jj] = jj; } /* proceed with RAD calibration */ } else { if( jj == header->calib->outCount || jj == header->calib->errorCount ) { Rad = INVALID_VALUE; } else { Rad = (double)jj * header->calib->gain_cnt2rad + header->calib->cnst_cnt2rad; } /* check if request is for RAD */ if( ounit_flag == 1 ) { Cal_Table[jj] = Rad; /* proceed with TEMP or ALB calibration */ } else { /* compute the TEMP if radiance band */ if( (header->calib->bandNo>=7 && strstr(header->basic->satName,"Himawari")!=NULL ) || (header->calib->bandNo>=2 && strstr(header->basic->satName,"MTSAT-2") !=NULL )) { if( Rad == INVALID_VALUE ) { Tbb = INVALID_VALUE; } else { if( Rad >= 2.000 && Rad <= 2.010 ) { sprintf(trace_string, " ******** Rad=%f",Rad); M0sxtrce( trace_string ); } hisd_radiance_to_tbb(header,Rad,&Tbb ); if( Rad >= 2.000 && Rad <= 2.010 ) { sprintf(trace_string, " ******** Tbb=%f",Tbb); M0sxtrce( trace_string ); } } /* check if request is for TEMP */ if( ounit_flag == 2 ) { if( Tbb == INVALID_VALUE ) { Cal_Table[jj] = 0; } else { Cal_Table[jj] = Tbb; } /* proceed with BRIT calibration */ } else { if( Tbb == INVALID_VALUE ) { value = 255; } else { rvalue = (float)Tbb; (void)gryscl_(&rvalue,&value); } if( value < 0 ) value = 0; if( value > 255 ) value = 255; Cal_Table[jj] = value; } /* compute the ALB */ } else { if(Rad == INVALID_VALUE){ Alb = INVALID_VALUE; } else { Alb = header->calib->rad2albedo * Rad; if( Alb < 0 ) Alb = INVALID_VALUE; } /* check if request is for ALB */ if( ounit_flag == 3 ) { if( Alb == INVALID_VALUE ) { Cal_Table[jj] = 0; } else { Cal_Table[jj] = Alb; } /* proceed with BRIT calibration */ } else { if( Alb == INVALID_VALUE ) { value = 0; } else { rvalue = 255.0*sqrt(Alb)+.5; value = (int)(rvalue); } if( value < 0 ) value = 0; if( value > 255 ) value = 255; Cal_Table[jj] = value; } } } } } } /* free the header */ hisd_free(header); /* close the file */ fclose( fp ); /* get out of this loop */ break; } /* increment the number of bands found */ found_bands[found_n_bands] = request_bands[0]; found_n_bands=1; sprintf(trace_string, "%s FOUND THIS BAND = %d",server,found_n_bands); M0sxtrce( trace_string ); /* ************************************************************************* Step 8: Sectorize the input image We can test if all the parts of the request we found in this image file. If so, "sector-izing" of the lat, lon and data is performed here. ************************************************************************* */ if( found_n_bands == request_n_bands ) { /* compute the SSEC resolution factor of the requested data */ sprintf(trace_string, "%s LINE: Base Res=%d Mag=%d",server,image_line_base_resolution, request_line_magnification); M0sxtrce( trace_string ); request_line_resolution = image_line_base_resolution * request_line_magnification; request_element_resolution = image_element_base_resolution * request_element_magnification; /* set request lines if 99999 was passed */ if( request_line_size == 99999 ) request_line_size = image_line_size; if( request_element_size == 99999 ) request_element_size = image_element_size; /* make sure that request is a multiple of 4 extend = request_element_size; test = request_element_size / 4 * 4; if (test != request_element_size ) { request_element_size = ((request_element_size/4)+1) * 4; } */ (void *)sprintf(trace_string,"%s COORDINATES: place=%s", server,place); M0sxtrce( trace_string ); /* check for Earth Centered/Upper load */ if( strncmp( place, "EC",2 ) == 0 || strncmp( place, "EU",2 ) == 0 ) { (void *)sprintf(trace_string,"%s COORDINATES: EARTH", server); M0sxtrce( trace_string ); /* initialize the navigation block */ for( i=0; i image_line_size || read_ending_line < 0 || read_starting_element > image_element_size || read_ending_element < 0 ) { (void)strcpy(requestBlock.errormsg, "The portion of the image requested does not exist"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } GEO_REQEST: (void *)sprintf(trace_string,"%s GEO REQUEST",server); M0sxtrce( trace_string ); /* determine the size of the geo request */ read_geo_line_size = image_line_size; read_geo_element_size = image_element_size; /* determine the coordinates of the geo request */ read_geo_starting_line = 0; read_geo_starting_element = 0; read_geo_ending_line = read_geo_starting_line + read_geo_line_size - 1; read_geo_ending_element = read_geo_starting_element + read_geo_element_size - 1; (void *)sprintf(trace_string,"%s READ GEO: line size=%d",server,read_geo_line_size); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"%s READ GEO: element size=%d",server,read_geo_element_size); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"%s READ GEO: starting line=%d",server,read_geo_starting_line); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"%s READ GEO: starting element=%d",server,read_geo_starting_element); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"%s READ GEO: ending line=%d",server,read_geo_ending_line); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"%s READ GEO: ending element=%d",server,read_geo_ending_element); M0sxtrce( trace_string ); /* ************************************************************************* Step 9: Transfer starts here Construct and send the following blocks: DIRECTORY NAVIGATION CALIBRATION AUXILLARY DATA ************************************************************************* */ (void *)sprintf(trace_string,"%s TRANSFER",server); M0sxtrce( trace_string ); /* size of the AUX block */ aux_block_size = 0; /* size of the line prefix */ psize = 0; /* size of the DAT block */ dat_block_size = (read_line_size / request_line_magnification) * ( (read_element_size / request_element_magnification) * request_n_bands * request_data_size + psize ); (void *)sprintf(trace_string,"%s dat_block_size=%d",server,dat_block_size); M0sxtrce( trace_string ); /* send the size of the total transmission */ totalSize = (DIR_SIZE*4) + (NAV_BLOCK_SIZE*4) + (CAL_BLOCK_SIZE*4) + (aux_block_size) + (dat_block_size); (void *)sprintf(trace_string,"%s total_block_size=%d",server,totalSize); M0sxtrce( trace_string ); requestBlock.reply_length += totalSize; M0swbyt4(&totalSize, 1); rc = M0sxsend(4, &totalSize ); /* ***************** DIRECTORY BLOCK ****************** */ /* initialize the directory block */ for( i=0; i= 41 && cover[image_dataset_position-1] <= 44) { memcpy( (void *)&directory[24], (void *)"HIMAWARI-8 Japan",16); } else if( cover[image_dataset_position-1] >= 51 && cover[image_dataset_position-1] <= 54) { memcpy( (void *)&directory[24], (void *)"HIMAWARI-8 Target",16); } else { memcpy( (void *)&directory[24], (void *)"HIMAWARI-8 Landmark",18); } /* band map */ for( i=0; i=0 && i=0 && j0 ) NCnt_valid++; } else { value = 0; } /* WRITE the McIDAS data array */ if( space == 1 ) { uchar_array[jj] = (unsigned char)value; } else if( space == 2 ) { ushort_array[jj] = (unsigned short)value; } else if( space == 4 ) { ulong_array[jj] = (unsigned long)value; } jj++; } /* correct byte order */ if( space == 4 ) { size = request_n_bands*request_element_size; M0swbyt4( &ulong_array[0], size ); } if( space == 2 ) { size = request_n_bands*request_element_size; (void)swbyt2_( ushort_array+(psize/2), &size ); size = size + (psize/2); } if( space == 1 ) { size = psize + (request_n_bands*request_element_size); } /* send out the line array */ if( space == 1 ) { rc = M0sxsend( size, uchar_array ); } else if( space == 2 ) { rc = M0sxsend( size*2, ushort_array ); } else if( space == 4 ) { rc = M0sxsend( size*4, ulong_array ); } } /* FREE the data array */ free( Data_array ); /* FREE the transmission array */ if( space == 1 ) free ( uchar_array) ; if( space == 2 ) free ( ushort_array ); if( space == 4 ) free ( ulong_array ); } } sprintf(trace_string, "%s Done NCnt=%d NCnt_valid=%d", server, NCnt, NCnt_valid); M0sxtrce( trace_string ); /* did we find anything??? */ sprintf(trace_string, "%s found_n_bands=%d", server, found_n_bands); M0sxtrce( trace_string ); if( found_n_bands <= 0 ) { sprintf(trace_string, "%s IN THE NO BANDS FOUND TEST", server); M0sxtrce( trace_string ); (void)strcpy(requestBlock.errormsg, "No images satisfy the selection criteria"); requestBlock.returncode = -30; M0sxdone(&requestBlock); return(0); } /* terminate the request */ M0sxdone(&requestBlock); return (0); } /* ---------------------------------------------------------------------------- Sample source code for reading Himawari Standard Data Copyright (C) 2014 MSC (Meteorological Satellite Center) of JMA Disclaimer : MSC does not guarantee regarding the correctness, accuracy, reliability, or any other aspect regarding use of these sample codes. Detail of Himawari Standard Format : For data structure of Himawari Standard Format, prelese refer to MSC Website and Himawari Standard Data User's Guid. MSC Website http://mscweb.kishou.go.jp/index.htm Himawari Standard Data User's Guid http://mscweb.kishou.go.jp/himawari89/space_segment/hsd_sample/HS_D_users_guide_en.pdf History April, 2014 First release ---------------------------------------------------------------------------- */ static int hisd_comp_table(HisdHeader *header); /* ---------------------------------------------------------------------------- hisd_read_header() ----------------------------------------------------------------------------*/ int hisd_read_header(HisdHeader *header ,FILE *fp){ char byteFlag; int ii; if(fp == NULL)return(ERROR_READ); /* allocate--------------------------------------------------------------*/ header->basic = (Basic_Info *)calloc(1,sizeof(Basic_Info)); header->data = (Data_Info *)calloc(1,sizeof(Data_Info)); header->proj = (Proj_Info *)calloc(1,sizeof(Proj_Info)); header->nav = (Navi_Info *)calloc(1,sizeof(Navi_Info)); header->calib = (Calib_Info *)calloc(1,sizeof(Calib_Info)); header->interCalib = (InterCalib_Info *)calloc(1,sizeof(InterCalib_Info)); header->seg = (Segm_Info *)calloc(1,sizeof(Segm_Info)); header->navcorr = (NaviCorr_Info *)calloc(1,sizeof(NaviCorr_Info)); header->obstime = (ObsTime_Info *)calloc(1,sizeof(ObsTime_Info)); header->error = (Error_Info *)calloc(1,sizeof(Error_Info)); header->spare = (Spare *)calloc(1,sizeof(Spare)); header->correct_table=(Correct_Table *)calloc(1,sizeof(Correct_Table)); if( header->basic == NULL || header->data == NULL || header->proj == NULL || header->nav == NULL || header->calib == NULL || header->interCalib == NULL || header->seg == NULL || header->navcorr == NULL || header->obstime == NULL || header->error == NULL || header->spare == NULL || header->correct_table == NULL ){ fprintf(stderr,"calloc error\n"); return(ERROR_ALLOCATE); } /* #1--------------------------------------------------------------------*/ if( (1>fread(&header->basic->HeaderNum,1,1,fp)) || (1>fread(&header->basic->BlockLen,2,1,fp)) || (1>fread(&header->basic->headerNum,2,1,fp)) || (1>fread(&header->basic->byteOrder,1,1,fp)) || (1>fread(&header->basic->satName,16,1,fp)) || (1>fread(&header->basic->proName,16,1,fp)) || (1>fread(&header->basic->ObsType1,4,1,fp)) || (1>fread(&header->basic->ObsType2,2,1,fp)) || (1>fread(&header->basic->TimeLine,2,1,fp)) || (1>fread(&header->basic->ObsStartTime,8,1,fp)) || (1>fread(&header->basic->ObsEndTime,8,1,fp)) || (1>fread(&header->basic->fileCreationMjd,8,1,fp)) || (1>fread(&header->basic->totalHeaderLen,4,1,fp)) || (1>fread(&header->basic->dataLen,4,1,fp)) || (1>fread(&header->basic->qflag1,1,1,fp)) || (1>fread(&header->basic->qflag2,1,1,fp)) || (1>fread(&header->basic->qflag3,1,1,fp)) || (1>fread(&header->basic->qflag4,1,1,fp)) || (1>fread(&header->basic->verName,32,1,fp)) || (1>fread(&header->basic->fileName,128,1,fp)) || (1>fread(&header->basic->spare,40,1,fp)) ){ fprintf(stderr,"header #1 read error\n"); return(ERROR_READ); } byteFlag=0; if(header->basic->BlockLen != 282){ /* byte swap */ byteFlag=1; swapBytes(&header->basic->BlockLen,2,1); swapBytes(&header->basic->headerNum,2,1); swapBytes(&header->basic->ObsStartTime,8,1); swapBytes(&header->basic->ObsEndTime,8,1); swapBytes(&header->basic->TimeLine,2,1); swapBytes(&header->basic->fileCreationMjd,8,1); swapBytes(&header->basic->totalHeaderLen,4,1); swapBytes(&header->basic->dataLen,4,1); } if( header->basic->HeaderNum != 1 || header->basic->BlockLen !=282 ){ fprintf(stderr,"header #1 read error\n"); fprintf(stderr,"HeaderNum=%d\n",header->basic->HeaderNum); fprintf(stderr,"BlockLen=%d\n" ,header->basic->BlockLen); return(ERROR_READ); } /* #2--------------------------------------------------------------------*/ if( (1>fread(&header->data->HeaderNum,1,1,fp)) || (1>fread(&header->data->BlockLen,2,1,fp)) || (1>fread(&header->data->bitPix,2,1,fp)) || (1>fread(&header->data->nPix,2,1,fp)) || (1>fread(&header->data->nLin,2,1,fp)) || (1>fread(&header->data->comp,1,1,fp)) || (1>fread(&header->data->spare,40,1,fp)) ){ fprintf(stderr,"header #2 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->data->BlockLen,2,1); swapBytes(&header->data->bitPix,2,1); swapBytes(&header->data->nPix,2,1); swapBytes(&header->data->nLin,2,1); } if( header->data->HeaderNum != 2 || header->data->BlockLen != 50 ){ fprintf(stderr,"header #2 read error\n"); fprintf(stderr,"HeaderNum=%d\n",header->data->HeaderNum); fprintf(stderr,"BlockLen=%d\n" ,header->data->BlockLen); return(ERROR_READ); } /* #3--------------------------------------------------------------------*/ if( (1>fread(&header->proj->HeaderNum,1,1,fp)) || (1>fread(&header->proj->BlockLen,2,1,fp)) || (1>fread(&header->proj->subLon,8,1,fp)) || (1>fread(&header->proj->cfac,4,1,fp)) || (1>fread(&header->proj->lfac,4,1,fp)) || (1>fread(&header->proj->coff,4,1,fp)) || (1>fread(&header->proj->loff,4,1,fp)) || (1>fread(&header->proj->satDis,8,1,fp)) || (1>fread(&header->proj->eqtrRadius,8,1,fp)) || (1>fread(&header->proj->polrRadius,8,1,fp)) || (1>fread(&header->proj->projParam1,8,1,fp)) || (1>fread(&header->proj->projParam2,8,1,fp)) || (1>fread(&header->proj->projParam3,8,1,fp)) || (1>fread(&header->proj->projParamSd,8,1,fp)) || (1>fread(&header->proj->resampleKind,2,1,fp)) || (1>fread(&header->proj->resampleSize,2,1,fp)) || (1>fread(&header->proj->spare,40,1,fp)) ){ fprintf(stderr,"header #3 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->proj->BlockLen,2,1); swapBytes(&header->proj->subLon,8,1); swapBytes(&header->proj->cfac,4,1); swapBytes(&header->proj->lfac,4,1); swapBytes(&header->proj->coff,4,1); swapBytes(&header->proj->loff,4,1); swapBytes(&header->proj->satDis,8,1); swapBytes(&header->proj->eqtrRadius,8,1); swapBytes(&header->proj->polrRadius,8,1); swapBytes(&header->proj->projParam1,8,1); swapBytes(&header->proj->projParam2,8,1); swapBytes(&header->proj->projParam3,8,1); swapBytes(&header->proj->projParamSd,8,1); swapBytes(&header->proj->resampleKind,2,1); swapBytes(&header->proj->resampleSize,2,1); } if( header->proj->HeaderNum != 3 || header->proj->BlockLen !=127 ){ fprintf(stderr,"header #3 read error\n"); fprintf(stderr,"HeaderNum=%d\n",header->proj->HeaderNum); fprintf(stderr,"BlockLen=%d\n" ,header->proj->BlockLen); return(ERROR_READ); } /* #4--------------------------------------------------------------------*/ if( (1>fread(&header->nav->HeaderNum,1,1,fp)) || (1>fread(&header->nav->BlockLen,2,1,fp)) || (1>fread(&header->nav->navMjd,8,1,fp)) || (1>fread(&header->nav->sspLon,8,1,fp)) || (1>fread(&header->nav->sspLat,8,1,fp)) || (1>fread(&header->nav->satDis,8,1,fp)) || (1>fread(&header->nav->nadirLon,8,1,fp)) || (1>fread(&header->nav->nadirLat,8,1,fp)) || (1>fread(&header->nav->sunPos_x,8,1,fp)) || (1>fread(&header->nav->sunPos_y,8,1,fp)) || (1>fread(&header->nav->sunPos_z,8,1,fp)) || (1>fread(&header->nav->moonPos_x,8,1,fp)) || (1>fread(&header->nav->moonPos_y,8,1,fp)) || (1>fread(&header->nav->moonPos_z,8,1,fp)) || (1>fread(&header->nav->spare,40,1,fp)) ){ fprintf(stderr,"header #4 read error\n"); return(1); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->nav->BlockLen,2,1); swapBytes(&header->nav->navMjd,8,1); swapBytes(&header->nav->sspLon,8,1); swapBytes(&header->nav->sspLat,8,1); swapBytes(&header->nav->satDis,8,1); swapBytes(&header->nav->nadirLon,8,1); swapBytes(&header->nav->nadirLat,8,1); swapBytes(&header->nav->sunPos_x,8,1); swapBytes(&header->nav->sunPos_y,8,1); swapBytes(&header->nav->sunPos_z,8,1); swapBytes(&header->nav->moonPos_x,8,1); swapBytes(&header->nav->moonPos_y,8,1); swapBytes(&header->nav->moonPos_z,8,1); } if( header->nav->HeaderNum != 4 || header->nav->BlockLen !=139 ){ fprintf(stderr,"header #4 read error\n"); fprintf(stderr,"HeaderNum=%d\n",header->nav->HeaderNum); fprintf(stderr,"BlockLen=%d\n" ,header->nav->BlockLen); return(ERROR_READ); } /* #5--------------------------------------------------------------------*/ if( (1>fread(&header->calib->HeaderNum,1,1,fp)) || (1>fread(&header->calib->BlockLen,2,1,fp)) || (1>fread(&header->calib->bandNo,2,1,fp)) || (1>fread(&header->calib->waveLen,8,1,fp)) || (1>fread(&header->calib->bitPix,2,1,fp)) || (1>fread(&header->calib->errorCount,2,1,fp)) || (1>fread(&header->calib->outCount,2,1,fp)) || (1>fread(&header->calib->gain_cnt2rad,8,1,fp)) || (1>fread(&header->calib->cnst_cnt2rad,8,1,fp)) ){ fprintf(stderr,"header #5 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->calib->BlockLen,2,1); swapBytes(&header->calib->bandNo,2,1); swapBytes(&header->calib->waveLen,8,1); swapBytes(&header->calib->bitPix,2,1); swapBytes(&header->calib->errorCount,2,1); swapBytes(&header->calib->outCount,2,1); swapBytes(&header->calib->gain_cnt2rad,8,1); swapBytes(&header->calib->cnst_cnt2rad,8,1); } if( (header->calib->bandNo>=7 && strstr(header->basic->satName,"Himawari")!=NULL ) || (header->calib->bandNo>=2 && strstr(header->basic->satName,"MTSAT-2") !=NULL ) ){ if( (1>fread(&header->calib->rad2btp_c0,8,1,fp)) || (1>fread(&header->calib->rad2btp_c1,8,1,fp)) || (1>fread(&header->calib->rad2btp_c2,8,1,fp)) || (1>fread(&header->calib->btp2rad_c0,8,1,fp)) || (1>fread(&header->calib->btp2rad_c1,8,1,fp)) || (1>fread(&header->calib->btp2rad_c2,8,1,fp)) || (1>fread(&header->calib->lightSpeed,8,1,fp)) || (1>fread(&header->calib->planckConst,8,1,fp)) || (1>fread(&header->calib->bolzConst,8,1,fp)) || (1>fread(&header->calib->spare,40,1,fp)) ){ fprintf(stderr,"header #5 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->calib->rad2btp_c0,8,1); swapBytes(&header->calib->rad2btp_c1,8,1); swapBytes(&header->calib->rad2btp_c2,8,1); swapBytes(&header->calib->btp2rad_c0,8,1); swapBytes(&header->calib->btp2rad_c1,8,1); swapBytes(&header->calib->btp2rad_c2,8,1); swapBytes(&header->calib->lightSpeed,8,1); swapBytes(&header->calib->planckConst,8,1); swapBytes(&header->calib->bolzConst,8,1); } }else{ if( (1>fread(&header->calib->rad2albedo,8,1,fp)) || (1>fread(&header->calib->spareV,104,1,fp))){ fprintf(stderr,"header #5 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->calib->rad2albedo,8,1); } } if( header->calib->HeaderNum != 5 || header->calib->BlockLen !=147 ){ fprintf(stderr,"header #5 read error\n"); fprintf(stderr,"HeaderNum=%d\n",header->calib->HeaderNum); fprintf(stderr,"BlockLen=%d\n" ,header->calib->BlockLen); return(ERROR_READ); } /* #6--------------------------------------------------------------------*/ if( (1>fread(&header->interCalib->HeaderNum,1,1,fp)) || (1>fread(&header->interCalib->BlockLen,2,1,fp)) || (1>fread(&header->interCalib->gsicsCorr_C,8,1,fp)) || (1>fread(&header->interCalib->gsicsCorr_C_er,8,1,fp)) || (1>fread(&header->interCalib->gsicsCorr_1,8,1,fp)) || (1>fread(&header->interCalib->gsicsCorr_1_er,8,1,fp)) || (1>fread(&header->interCalib->gsicsCorr_2,8,1,fp)) || (1>fread(&header->interCalib->gsicsCorr_2_er,8,1,fp)) || (1>fread(&header->interCalib->gsicsCorr_StrMJD,8,1,fp)) || (1>fread(&header->interCalib->gsicsCorr_EndMJD,8,1,fp)) || (1>fread(&header->interCalib->gsicsCorrInfo,1,64,fp)) || (1>fread(&header->interCalib->spare,128,1,fp)) ){ fprintf(stderr,"header #6 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->interCalib->BlockLen,2,1); swapBytes(&header->interCalib->gsicsCorr_2,8,1); swapBytes(&header->interCalib->gsicsCorr_2_er,8,1); swapBytes(&header->interCalib->gsicsCorr_1,8,1); swapBytes(&header->interCalib->gsicsCorr_1_er,8,1); swapBytes(&header->interCalib->gsicsCorr_C,8,1); swapBytes(&header->interCalib->gsicsCorr_C_er,8,1); swapBytes(&header->interCalib->gsicsCorr_StrMJD,8,1); swapBytes(&header->interCalib->gsicsCorr_EndMJD,8,1); } if( header->interCalib->HeaderNum != 6 || header->interCalib->BlockLen !=259 ){ fprintf(stderr,"header #6 read error\n"); fprintf(stderr,"HeaderNum=%d\n",header->interCalib->HeaderNum); fprintf(stderr,"BlockLen=%d\n" ,header->interCalib->BlockLen); return(ERROR_READ); } /* #7--------------------------------------------------------------------*/ if( (1>fread(&header->seg->HeaderNum,1,1,fp)) || (1>fread(&header->seg->BlockLen,2,1,fp)) || (1>fread(&header->seg->totalSegNum,1,1,fp)) || (1>fread(&header->seg->segSeqNo,1,1,fp)) || (1>fread(&header->seg->strLineNo,2,1,fp)) || (1>fread(&header->seg->spare,40,1,fp)) ){ fprintf(stderr,"header #7 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->seg->BlockLen,2,1); swapBytes(&header->seg->strLineNo,2,1); } if( header->seg->HeaderNum != 7 || header->seg->BlockLen !=47 ){ fprintf(stderr,"header #7 read error\n"); fprintf(stderr,"HeaderNum=%d\n",header->seg->HeaderNum); fprintf(stderr,"BlockLen=%d\n" ,header->seg->BlockLen); return(ERROR_READ); } /* #8--------------------------------------------------------------------*/ if( (1>fread(&header->navcorr->HeaderNum,1,1,fp)) || (1>fread(&header->navcorr->BlockLen,2,1,fp)) || (1>fread(&header->navcorr->RoCenterColumn,4,1,fp)) || (1>fread(&header->navcorr->RoCenterLine,4,1,fp)) || (1>fread(&header->navcorr->RoCorrection,8,1,fp)) || (1>fread(&header->navcorr->correctNum,2,1,fp)) ){ fprintf(stderr,"header #8 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->navcorr->BlockLen,2,1); swapBytes(&header->navcorr->RoCenterColumn,4,1); swapBytes(&header->navcorr->RoCenterLine,4,1); swapBytes(&header->navcorr->RoCorrection,8,1); swapBytes(&header->navcorr->correctNum,2,1); } header->navcorr->lineNo=(unsigned short *) calloc(header->navcorr->correctNum,sizeof(unsigned short)); header->navcorr->columnShift=(float *) calloc(header->navcorr->correctNum,sizeof(float)); header->navcorr->lineShift=(float *) calloc(header->navcorr->correctNum,sizeof(float)); if(header->navcorr->lineNo == NULL || header->navcorr->columnShift == NULL || header->navcorr->lineShift == NULL ){ fprintf(stderr,"calloc error\n"); return(ERROR_READ); } for(ii=0;iinavcorr->correctNum;ii++){ if( (1>fread(&header->navcorr->lineNo[ii],2,1,fp)) || (1>fread(&header->navcorr->columnShift[ii],4,1,fp)) || (1>fread(&header->navcorr->lineShift[ii],4,1,fp)) ){ fprintf(stderr,"header #8 read error\n"); return(ERROR_READ); } } if( (1>fread(&header->navcorr->spare,40,1,fp))){ fprintf(stderr,"header #8 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->navcorr->lineNo[0],2,header->navcorr->correctNum); swapBytes(&header->navcorr->columnShift[0],4, header->navcorr->correctNum); swapBytes(&header->navcorr->lineShift[0],4,header->navcorr->correctNum); } if( header->navcorr->HeaderNum != 8 ){ fprintf(stderr,"header #8 read error\n"); fprintf(stderr,"HeaderNum=%d\n",header->navcorr->HeaderNum); fprintf(stderr,"BlockLen=%d\n" ,header->navcorr->BlockLen); return(ERROR_READ); } /* #9--------------------------------------------------------------------*/ if( (1>fread(&header->obstime->HeaderNum,1,1,fp)) || (1>fread(&header->obstime->BlockLen,2,1,fp)) || (1>fread(&header->obstime->obsNum,2,1,fp)) ){ fprintf(stderr,"header #9 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->obstime->BlockLen,2,1); swapBytes(&header->obstime->obsNum,2,1); } header->obstime->lineNo=(unsigned short *) calloc(header->obstime->obsNum,sizeof(unsigned short)); header->obstime->obsMJD=(double *) calloc(header->obstime->obsNum,sizeof(double)); if(header->obstime->lineNo == NULL || header->obstime->obsMJD == NULL ){ fprintf(stderr,"calloc error\n"); return(ERROR_READ); } for(ii=0;iiobstime->obsNum;ii++){ if( (1>fread(&header->obstime->lineNo[ii],2,1,fp)) || (1>fread(&header->obstime->obsMJD[ii],8,1,fp)) ){ fprintf(stderr,"header #9 read error\n"); return(ERROR_READ); } } if( (1>fread(&header->obstime->spare,40,1,fp)) ){ fprintf(stderr,"header #9 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->obstime->lineNo[0],2,header->obstime->obsNum); swapBytes(&header->obstime->obsMJD[0],8,header->obstime->obsNum); } if( header->obstime->HeaderNum != 9 ){ fprintf(stderr,"header #9 read error\n"); fprintf(stderr,"HeaderNum=%d\n",header->obstime->HeaderNum); fprintf(stderr,"BlockLen=%d\n" ,header->obstime->BlockLen); return(ERROR_READ); } /* #10-------------------------------------------------------------------*/ if( (1>fread(&header->error->HeaderNum,1,1,fp)) || (1>fread(&header->error->BlockLen,4,1,fp)) || (1>fread(&header->error->errorNum,2,1,fp)) ){ fprintf(stderr,"header #10 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->error->BlockLen,4,1); swapBytes(&header->error->errorNum,2,1); } if(header->error->errorNum>=1){ header->error->lineNo=(unsigned short *) calloc(header->error->errorNum,sizeof(unsigned short)); header->error->errPixNum=(unsigned short *) calloc(header->error->errorNum,sizeof(unsigned short)); if(header->error->lineNo == NULL || header->error->errPixNum == NULL ){ fprintf(stderr,"calloc error\n"); return(ERROR_READ); } for(ii=0;iierror->errorNum;ii++){ if( (1>fread(&header->error->lineNo[ii],2,1,fp)) || (1>fread(&header->error->errPixNum[ii],2,1,fp)) ){ fprintf(stderr,"header #10 read error\n"); return(ERROR_READ); } } } if( (1>fread(&header->error->spare,40,1,fp)) ){ fprintf(stderr,"header #10 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->error->lineNo[0],2,header->error->errorNum); swapBytes(&header->error->errPixNum[0],2,header->error->errorNum); } if( header->error->HeaderNum !=10 ){ fprintf(stderr,"header #10 read error\n"); fprintf(stderr,"HeaderNum=%d\n",header->error->HeaderNum); fprintf(stderr,"BlockLen=%d\n" ,header->error->BlockLen); return(ERROR_READ); } /* #11-------------------------------------------------------------------*/ if( (1>fread(&header->spare->HeaderNum,1,1,fp)) || (1>fread(&header->spare->BlockLen,2,1,fp)) || (1>fread(&header->spare->spare,256,1,fp)) ){ fprintf(stderr,"header #11 read error\n"); return(ERROR_READ); } if(byteFlag==1){ swapBytes(&header->spare->BlockLen,2,1); } if( header->spare->HeaderNum !=11 || header->spare->BlockLen !=259 ){ fprintf(stderr,"header #11 read error\n"); fprintf(stderr,"HeaderNum=%d\n",header->spare->HeaderNum); fprintf(stderr,"BlockLen=%d\n" ,header->spare->BlockLen); return(ERROR_READ); } /* ----------------------------------------------------------------------*/ return(NORMAL_END); } /* ---------------------------------------------------------------------------- hisd_free() ----------------------------------------------------------------------------*/ int hisd_free(HisdHeader *header){ if(header->navcorr->lineNo != NULL){ free(header->navcorr->lineNo); } if(header->navcorr->columnShift != NULL){ free(header->navcorr->columnShift);} if(header->navcorr->lineShift != NULL){ free(header->navcorr->lineShift); } if(header->obstime->lineNo != NULL){ free(header->obstime->lineNo); } if(header->obstime->obsMJD != NULL){ free(header->obstime->obsMJD); } if(header->error->lineNo != NULL){ free(header->error->lineNo); } if(header->error->errPixNum != NULL){ free(header->error->errPixNum); } if(header->correct_table->cmpCoff != NULL){ free(header->correct_table->cmpCoff);} if(header->correct_table->cmpLoff != NULL){ free(header->correct_table->cmpLoff);} if(header->correct_table != NULL){ free(header->correct_table);} if(header->basic != NULL) { free(header->basic); } if(header->data != NULL) { free(header->data); } if(header->proj != NULL) { free(header->proj); } if(header->nav != NULL) { free(header->nav); } if(header->calib != NULL) { free(header->calib); } if(header->interCalib != NULL) { free(header->interCalib); } if(header->seg != NULL) { free(header->seg); } if(header->navcorr != NULL) { free(header->navcorr); } if(header->obstime != NULL) { free(header->obstime); } if(header->error != NULL) { free(header->error); } if(header->spare != NULL) { free(header->spare); } return(NORMAL_END); } /* ---------------------------------------------------------------------------- hisd_getdata_by_pixlin() ----------------------------------------------------------------------------*/ int hisd_getdata_by_pixlin( HisdHeader *header, FILE *fp, int iLin, int iPix, unsigned short *sout ) { float pShift; float lShift; int ii; float hPix; float hLin; static double cosR_A; static double sinR_A; float R_C = header->navcorr->RoCenterColumn; float R_L = header->navcorr->RoCenterLine; double R_A = header->navcorr->RoCorrection / 1000. / 1000.; /* make line/pixel float */ hPix = 0.0; /* hLin = (float)iLin; */ hLin = (float)iLin; /* init */ /* *sout=header->calib->outCount; count value of out of scan area */ /* shift correction table */ if(header->correct_table->flag == 0){ hisd_comp_table(header); cosR_A = cos(R_A); sinR_A = sin(R_A); } /* rotational correction */ float hLin2= (hLin - R_L) * cosR_A - (hPix - R_C) * sinR_A + R_L; float hPix2= (hPix - R_C) * cosR_A + (hLin - R_L) * sinR_A + R_C; hLin = hLin2; hPix = hPix2; if(header->correct_table->flag == 1){ /* shift amount for column and line direction */ ii = (short)(hLin+0.5) - header->correct_table->startLineNo; if(ii < 0 ){ ii = 0; }else if(header->correct_table->lineNum -1 < ii){ ii = header->correct_table->lineNum -1; } pShift = header->correct_table->cmpCoff[ii]; lShift = header->correct_table->cmpLoff[ii]; /* shift correction */ hPix = hPix - pShift; hLin = hLin - lShift; } /* Nearest Neighbor int ll = (int)(hLin + 0.5) - header->seg->strLineNo ; int pp = (int)(hPix + 0.5); */ int ll = (int)(hLin) - header->seg->strLineNo ; int pp = (int)(hPix); /* check ll and pp */ if(ll < 0 || header->data->nLin < ll){ return(ERROR_DATA); } if(pp < 0 || header->data->nPix < pp){ return(ERROR_DATA); } /* read data */ unsigned long seek_pt; seek_pt = header->basic->totalHeaderLen + (ll * header->data->nPix + pp ) * sizeof(unsigned short); if(seek_pt > header->basic->totalHeaderLen + header->basic->dataLen){ return(ERROR_DATA); } fseek(fp,seek_pt,SEEK_SET); fread(sout,sizeof(unsigned short), iPix,fp); /* byte swap */ if(byteOrder() != header->basic->byteOrder ){ swapBytes(sout,sizeof(unsigned short),iPix); } return(NORMAL_END); } /* ---------------------------------------------------------------------------- hisd_radiance_to_tbb() ----------------------------------------------------------------------------*/ void hisd_radiance_to_tbb(HisdHeader *header,double radiance,double *tbb){ double effective_temperature; static char trace_string[500]; int init; /* central wave length */ double lambda = header->calib->waveLen / 1000000.0; init = 1; /*if( radiance >= 5.850 && radiance <= 5.860 ) init = 0;*/ /* radiance */ radiance = radiance * 1000000.0; /* planck_c1 = (2 * h * c^2 / lambda^5) */ double planck_c1= 2.0 * header->calib->planckConst * pow(header->calib->lightSpeed,2.0) / pow(lambda,5.0) ; /* planck_c2 = (h * c / k / lambda ) */ double planck_c2= header->calib->planckConst * header->calib->lightSpeed / header->calib->bolzConst / lambda ; if(radiance > 0 ){ effective_temperature = planck_c2 / log( (planck_c1 / radiance ) + 1.0 ); *tbb = header->calib->rad2btp_c0 + header->calib->rad2btp_c1 * effective_temperature + header->calib->rad2btp_c2 * pow(effective_temperature,2.0); /* if( init == 0 ) { (void *)sprintf(trace_string,"RAD=%f",radiance); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"RAD2BTP:C0=%f",rad2btp_c0); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"RAD2BTP:C1=%f",rad2btp_c1); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"RAD2BTP:C2=%f",rad2btp_c2); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"PLANCK: C1=%f",planck_c1); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"PLANCK: C2=%f",planck_c2); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"PLANCK: EffT=%f",effective_temperature); M0sxtrce( trace_string ); (void *)sprintf(trace_string,"TBB=%f ",*tbb); M0sxtrce( trace_string ); } */ }else{ *tbb = INVALID_VALUE; } } /* ---------------------------------------------------------------------------- hisd_comp_table() ----------------------------------------------------------------------------*/ static int hisd_comp_table(HisdHeader *header){ int ii,jj; int lineNo,lineNo1,lineNo2; float coff1,coff2,loff1,loff2; /* allocate */ if(header->correct_table==NULL){ header->correct_table = (Correct_Table *)calloc(1,sizeof(Correct_Table)); if(header->correct_table == NULL){ fprintf(stderr,"calloc error\n"); return(ERROR_ALLOCATE); } } /* make table */ if(header->navcorr->correctNum <2){ return(ERROR_DATA); }else{ header->correct_table->startLineNo = header->navcorr->lineNo[0]; header->correct_table->lineNum = header->navcorr->lineNo[header->navcorr->correctNum-1] - header->navcorr->lineNo[0] +1; if(header->correct_table->lineNum < 2){ return(ERROR_DATA); } /* allocate */ header->correct_table->cmpCoff = (float *)calloc(header->correct_table->lineNum,sizeof(float)); header->correct_table->cmpLoff = (float *)calloc(header->correct_table->lineNum,sizeof(float)); if(header->correct_table->cmpCoff == NULL || header->correct_table->cmpLoff == NULL ){ fprintf(stderr,"calloc error\n"); return(ERROR_ALLOCATE); } /* init */ for(ii=0;iicorrect_table->lineNum;ii++){ header->correct_table->cmpCoff[ii] = 0.0; header->correct_table->cmpLoff[ii] = 0.0; } /* navigation correction table */ for(ii=0;iinavcorr->correctNum-1;ii++){ lineNo1 = (int)header->navcorr->lineNo[ii]; lineNo2 = (int)header->navcorr->lineNo[ii+1]; coff1 = header->navcorr->columnShift[ii]; coff2 = header->navcorr->columnShift[ii+1]; loff1 = header->navcorr->lineShift[ii]; loff2 = header->navcorr->lineShift[ii+1]; for(lineNo=lineNo1;lineNo<=lineNo2;lineNo++){ jj = lineNo - header->correct_table->startLineNo; if( 0<= jj && jj <= header->correct_table->lineNum){ header->correct_table->cmpCoff[jj] = coff1 + (coff2 - coff1 ) / (float)(lineNo2 - lineNo1 ) * (float)(lineNo - lineNo1); header->correct_table->cmpLoff[jj] = loff1 + (loff2 - loff1 ) / (float)(lineNo2 - lineNo1 ) * (float)(lineNo - lineNo1); } } } } header->correct_table->flag = 1; return(NORMAL_END); } int swapBytes( void *buf, /**< inout */ int size, /**< in */ int nn ) /**< in */ { char *ba, *bb, *buf2 = buf; while( nn-- ) { bb = ( ba = buf2 ) + size -1; do { char a; a = *ba; *ba = *bb; *bb = a; } while( ++ba < --bb ); buf2 += size; } return 0; } int byteOrder(void){ int i=1; if( *((char *)&i) ) return 0; /* Little */ else if( *( (char *)&i + (sizeof(int)-1) ) ) return 1; /* Big */ else return -1; } int lonlat_to_pixlin(HisdHeader *head,double lon,double lat, float *pix,float *lin){ /* (1) init */ *pix = -9999.0; *lin = -9999.0; /* (2) check latitude */ if(lat < -90.0 || 90.0 < lat ){ return(ERROR_END); } /* (3) check longitude */ while(lon > 180.0){ lon-=360.0; } while(lon <-180.0){ lon+=360.0; } /* (4) degree to radian */ lon = lon * DEGTORAD; lat = lat * DEGTORAD; /* (5) geocentric latitude */ /* Global Specification 4.4.3.2 */ /* phi = arctan( (Rpol^2)/(Req^2) * tan(lat) ) */ /* */ /* (Rpol^2)/(Req^2) = head->proj->projParam2 */ double phi = atan( head->proj->projParam2 * tan(lat) ); /* (6) The length of Re */ /* Re = (Rpol) / sqrt( 1 - (Req^2 - Rpol^2) / Req^2 * cos^2(phi) ) */ /* */ /* Rpol = head->proj->polrRadius */ /* (Req^2 - Rpol^2) / Req^2 = head->proj->projParam1 */ double Re = (head->proj->polrRadius) / sqrt(1 - head->proj->projParam1 * cos(phi) * cos(phi)) ; /* (7) The cartesian components of the vector rs result as follows: */ /* r1 = h - Re * cos(phi) * cos(Le-Ld) */ /* r2 = -Re * cos(phi) * sin(Le-Ld) */ /* r3 = Re * sin(phi) */ /* */ /* Le : longitude */ /* Ld : sub_lon = head->proj->subLon */ /* h : distance from Earth's center to satellite (=head->proj->satDis) */ double r1 = head->proj->satDis - Re * cos(phi) * cos( lon - head->proj->subLon * DEGTORAD ); double r2 = - Re * cos(phi) * sin( lon - head->proj->subLon * DEGTORAD ); double r3 = Re * sin(phi); /* (8) check seeablibity */ double vx = Re * cos(phi) * cos( lon - head->proj->subLon * DEGTORAD ); if(0 < -r1 * vx - r2 * r2 + r3 * r3){ return(ERROR_END); } /* (9) The projection function is as follows: */ /* x = arctan(-r2/r1) */ /* y = arcsin(r3/rn) */ /* rn = sqrt(r1^2 + r2^2 + r3^2) */ double rn = sqrt(r1*r1 + r2*r2 + r3*r3); double x = atan2(-r2,r1) * RADTODEG; double y = asin(-r3/rn) * RADTODEG; /* (10) */ /* Global Specification 4.4.4 */ /* c = COFF + nint(x * 2^-16 * CFAC) */ /* l = LOFF + nint(y * 2^-16 * LFAC) */ float c = head->proj->coff + x * SCLUNIT * head->proj->cfac; float l = head->proj->loff + y * SCLUNIT * head->proj->lfac; *pix = c; *lin = l; return(NORMAL_END); } int pixlin_to_lonlat(HisdHeader *head,float pix,float lin, double *lon,double *lat){ double Sd,Sn,S1,S2,S3,Sxy; /* (0) init */ *lon = -9999; *lat = -9999; /* (1) pix,lin ==> c,l */ float c = pix; float l = lin; /* (2) the intermediate coordinates (x,y) */ /* Global Specification 4.4.4 Scaling Function */ /* c = COFF + nint(x * 2^-16 * CFAC) */ /* l = LOFF + nint(y * 2^-16 * LFAC) */ /* The intermediate coordinates (x,y) are as follows : */ /* x = (c -COFF) / (2^-16 * CFAC) */ /* y = (l -LOFF) / (2^-16 * LFAC) */ /* SCLUNIT = 2^-16 */ double x = DEGTORAD * ( c - head->proj->coff) / ( SCLUNIT * head->proj->cfac); double y = DEGTORAD * ( l - head->proj->loff) / ( SCLUNIT * head->proj->lfac); /* (3) longtitude,latitude */ /* Global Specification 4.4.3.2 */ /* The invers projection function is as follows : */ /* lon = arctan(S2/S1) + sub_lon */ /* lat = arctan( (Req^2/Rpol^2) * S3 / Sxy ) */ /* */ /* Thererin the variables S1,S2,S3,Sxy are as follows : */ /* S1 = Rs - Sn * cos(x) * cos(y) */ /* S2 = Sn * sin(x) * cos(y) */ /* S3 =-Sn * sin(y) */ /* Sxy = sqrt(S1^2 + S2^2) */ /* Sn =(Rs * cos(x) * cos(y) - Sd ) / */ /* (cos(y) * cos(y) + (Req^2/Rpol^2) * sin(y) * sin(y)) */ /* Sd =sqrt( (Rs * cos(x) * cos(y))^2 */ /* - ( cos(y) * cos(y) + (Req^2/Rpol^2) * sin(y) * sin(y) ) */ /* * (Rs^2 - Req^2) */ /* The variables Rs,Rpol,Req,(Req^2/Rpol^2),(Rs^2 - Req^2) are as follows : */ /* Rs : distance from Earth center to satellite= head->proj->satDis */ /* Rpol: polar radius of the Earth = head->proj->polrRadius */ /* Req : equator raidus of the Earth = head->proj->eqtrRadius */ /* (Req^2/Rpol^2) = head->proj->projParam3 */ /* (Rs^2 - Req^2) = head->proj->projParamSd */ Sd = (head->proj->satDis * cos(x) * cos(y)) * (head->proj->satDis * cos(x) * cos(y)) - (cos(y) * cos(y) + head->proj->projParam3 * sin(y) * sin(y)) * head->proj->projParamSd; if(Sd < 0 ) { return(ERROR_END);} else { Sd = sqrt(Sd);} Sn = (head->proj->satDis * cos(x) * cos(y) -Sd) / (cos(y) * cos(y) + head->proj->projParam3 * sin(y) * sin(y)); S1 = head->proj->satDis - (Sn * cos(x) * cos(y)); S2 = Sn * sin(x) * cos(y); S3 =-Sn * sin(y); Sxy=sqrt( S1 * S1 + S2 * S2); *lon = RADTODEG * atan2(S2,S1) + head->proj->subLon; *lat = RADTODEG * atan(head->proj->projParam3 * S3 / Sxy); /* (4) check longtitude */ while(*lon > 180.0 ){ *lon = *lon-360.0;} while(*lon <-180.0 ){ *lon = *lon+360.0;} return(NORMAL_END); } int **malloc2dint(int nr, int nc) { /* Allocate integer 2d array */ int dim2, i; int **buf; int *temp1; dim2 = nr * nc; temp1 = (int*) malloc(dim2 * sizeof(int)); buf = (int**) malloc(nr * sizeof(int*)); for ( i = 0; i < nr; i++) { buf[i] = temp1 + nc*i; } return buf; } float **malloc2dfloat(int nr, int nc) { /* Allocate float 2d array */ int dim2, i; float **buf; float *temp1; dim2 = nr * nc; temp1 = (float*) malloc(dim2 * sizeof(float)); buf = (float**) malloc(nr * sizeof(float*)); for ( i = 0; i < nr; i++) { buf[i] = temp1 + nc*i; } return buf; }