/***************************************************************************************** slice_response_code.h Author: Michael P. Owen 8/07/2007 MERS Laboratory, Brigham Young University Description: C subroutines to obtain and correctly interpolate the QuikSCAT slice spatial response. A bi-linear interpolation and transformation algorithm is implemented, interpolating first in orbit time and then in antenna azimuth from the tabulated files. The spatial response data files must be loaded prior to usage and the strings Vslicefile and Hslicefile must point to the correct locations of each file. Notes: -This code loads the entire response functions into memory (can be larger than 2GB). The code can be modified to load only the desired sections. -All global variables that are allocated need to be freed by the caller after usage. The necessary calls are: free(Vdel_lat); free(Vdel_lon); free(Vresponse); free(Vpoints); free(Hdel_lat); free(Hdel_lon); free(Hresponse); free(Hpoints); free(pulse_del_lat); free(pulse_del_lon); free(pulse_points); -The arrays to hold all of the loaded data are allocated dynamically. -Global variables are used to simplify data organization and pointer usage. -interpolate_slice performs all needed interpolation and indexing to update the pulse pointers. -To obtain the response functions for a specific slice the pulse_response pointer can be indexed as follows: pulse_response[jjjj + mmm*nnn*ii]; where ii is the slice number from 0 to 7, and jjjj is an index variable that is always < mmm*nnn. ********************************************************************************************/ #include #include #include /*--------------------- Global vars -------------------------*/ /* Filenames and paths of the spatial response data files */ char Vslicefile[200]="/auto/users/stuart/src/linux/Xshape/slice.outer.backup"; char Hslicefile[200]="/auto/users/stuart/src/linux/Xshape/slice.inner.backup"; /* File pointer */ FILE *resp_file; /* File size indicators */ int mmm, nnn; /* array(m by n) dimensions of one response */ int num_looks, num_positions; /* number of azimuth looks and orbit positions stored in the data */ /*Vpol arrays for the entire response file*/ int *Vpoints, *Vcen_points; float *Vdel_lat, *Vdel_lon, *Vresponse; /*Hpol arrays for the entire response file*/ int *Hpoints, *Hcen_points; float *Hdel_lat, *Hdel_lon, *Hresponse; /* arrays to hold the interpolated response for a given pulse*/ int *pulse_points, *pulse_cen_points; float *pulse_del_lat, *pulse_del_lon; /* pointers to the interpolated lat and lon fields for a pulse */ float *pulse_response; /* never allocated since it is merely a pointer to a place in the H or V response arrays */ /*---------------- function prototypes --------------------*/ /* load all of the global response arrays*/ void loadresponse(); /* performs interpolation from the data files to obtain the slice response, lat and lon fields for the given orbit time and azimuth position */ /* c_lat and c_lon are the pulse centers as reported in the L1B data and are required for proper interpolation */ void interpolate_slice(float antenna_azimuth, float orbit_position, int polarization, float c_lat, float c_lon); /* get the significant points for a pulse to use in lat/lon field interpolation */ void getpoints(float *x1, float *x2, float *x3, float *y1, float *y2, float *y3, int look, int position, int polarization, float c_lat); /* update the pulse_response pointer to the correct place in the response arrays*/ void getresponse(int look, int position , int polarization); /* free dynamically allocated arrays*/ //void free_response(); int nint(float r) /* nearest integer */ { int ret_val = r; if (ret_val - r > 0.5) ret_val--; if (r - ret_val > 0.5) ret_val++; return(ret_val); } /*-------------- Functions -----------------*/ /* Coastal: Subroutines used in Coastal Processing using the Actual antenna response */ void interpolate_slice(float antenna_azimuth, float orbit_time , int polarization, float c_lat, float c_lon ){ /* polarization == 0 (Hpol), == 1(Vpol) */ int lower_deg, upper_deg; float *upper_lat, *upper_lon, *lower_lat, *lower_lon; int i; float x,y,r_x, r_y, t_x, t_y; float ux1, ux2, ux3, uy1, uy2, uy3;/* upper look points*/ float lx1, lx2, lx3, ly1, ly2, ly3;/* lower look points*/ float stretch_x, stretch_y, stretch_x1, stretch_y1, stretch_x2, stretch_y2, rotate_deg; float temp, distux, distuy, distlx, distly, distx1, disty1, distx2, disty2; float Ux1, Uy1, Ux2, Uy2, Lx1, Ly1, Lx2, Ly2; float Urx1, Ury1, Urx2, Ury2, Lrx1, Lry1, Lrx2, Lry2; float orbit_scale, azimuth_scale; float lowx1, lowx2, lowy1, lowy2; float upx1, upx2, upy1, upy2; float a, b, c, d; float deg2rad = 3.14159 /180; float rad2deg = 180/(3.14159); float *look_response; float orbit_time_table[36]= { 0., 515., 745., 945., 1115., 1255., 1365., 1445., 1495., 1515., 1535., 1585., 1665., 1775., 1915., 2085., 2285., 2515., 3030., 3545., 3775., 3975., 4145., 4285., 4395., 4475., 4525., 4545., 4565., 4615., 4695., 4805., 4945., 5115., 5315., 5545. } ; float table_orbit_period = 6061.0; int orbit_position = 0; int low_position, up_position; int bilinear = 1; float Re, Rphi, Kflat = 1/298.257; float Ra = 6378.1363; /* find the correct look angle to interpolate from */ lower_deg =( nint(antenna_azimuth/10.0)); if(antenna_azimuth > lower_deg*10) upper_deg = lower_deg + 1; else{ upper_deg = lower_deg; lower_deg = upper_deg - 1; } /* check boundary cases */ if(lower_deg == -1) lower_deg = 35; if(upper_deg == 36) upper_deg = 0; if(upper_deg == 37){ upper_deg = 1; lower_deg = 0; } /* find the correct orbit positions */ if (orbit_time < 0.0) orbit_time=orbit_time+table_orbit_period; if (orbit_time >= table_orbit_period){ orbit_time = fmod(orbit_time,table_orbit_period);} for (i=0; i= orbit_time_table[i] & orbit_time < orbit_time_table[i+1]){ low_position = i; up_position = i+1; goto label_5100; } } for(i=0;i<36;i++){ printf("Orbit_time_table[%d]: %f\n", i, orbit_time_table[i]); } printf("num_positions: %d\n", num_positions); printf("Orbit_time: %f\n", orbit_time); printf("Orbit_time_error: exiting...\n"); exit(-1); label_5100: //printf("orbit_time = %f\n", orbit_time); //printf("low_position = %f\n", orbit_time_table[low_position]); //printf("up_position = %f\n", orbit_time_table[up_position]); // printf("low_position = %d, up_position = %d, lower_deg = %d, upper_deg = %d\n", low_position, up_position, lower_deg, upper_deg); /* find the correct interpolated vectors between the same upper looks of the two orbit positions */ getpoints(&ux1, &ux2, &ux3, &uy1, &uy2, &uy3, upper_deg, up_position, polarization, c_lat); getpoints(&lx1, &lx2, &lx3, &ly1, &ly2, &ly3, upper_deg, low_position, polarization, c_lat); /* rotate the upper_deg reference points to the LEFT by 10deg to be on the same axis */ temp = (ux1 * cos( (10.0*deg2rad)) + (uy1 *-sin(10.0*deg2rad) ) ); uy1 = ( ux1 * sin( (10.0*deg2rad)) + (uy1 * cos(10.0*deg2rad) ) ); ux1 = temp; temp = (ux2 * cos( (10.0*deg2rad)) + (uy2 *-sin(10.0*deg2rad))); uy2 = ( ux2 * sin( (10.0*deg2rad)) + (uy2 * cos(10.0*deg2rad))); ux2 = temp; temp = (ux3 * cos( (10.0*deg2rad)) + (uy3 *-sin(10.0*deg2rad))); uy3 = ( ux3 * sin( (10.0*deg2rad)) + (uy3 * cos(10.0*deg2rad))); ux3 = temp; temp = (lx1 * cos( (10.0*deg2rad)) + (ly1 *-sin(10.0*deg2rad))); ly1 = ( lx1 * sin( (10.0*deg2rad)) + (ly1 * cos(10.0*deg2rad))); lx1 = temp; temp = (lx2 * cos( (10.0*deg2rad)) + (ly2 *-sin(10.0*deg2rad))); ly2 = ( lx2 * sin( (10.0*deg2rad)) + (ly2 * cos(10.0*deg2rad))); lx2 = temp; temp = (lx3 * cos( (10.0*deg2rad)) + (ly3 *-sin(10.0*deg2rad))); ly3 = ( lx3 * sin( (10.0*deg2rad)) + (ly3 * cos(10.0*deg2rad))); lx3 = temp; /* make set of points into two vectors for each position */ Lrx1 = lx2 - lx1; Lry1 = ly2 - ly1; Lrx2 = lx3 - lx1; Lry2 = ly3 - ly1; Urx1 = ux2 - ux1; Ury1 = uy2 - uy1; Urx2 = ux3 - ux1; Ury2 = uy3 - uy1; /* Scale factor between positions */ if(up_position == 0) orbit_scale = (orbit_time - orbit_time_table[low_position])/(table_orbit_period - orbit_time_table[low_position]); else orbit_scale = (orbit_time - orbit_time_table[low_position])/(orbit_time_table[up_position] - orbit_time_table[low_position]); /* calculate the interpolated vector components between orbit_positions for the lower_deg*/ upx1 = (Urx1-Lrx1)*orbit_scale + Lrx1; upy1 = (Ury1-Lry1)*orbit_scale + Lry1; upx2 = (Urx2-Lrx2)*orbit_scale + Lrx2; upy2 = (Ury2-Lry2)*orbit_scale + Lry2; /* find the correct interpolated vectors between the same looks of the two orbit positions */ getpoints(&ux1, &ux2, &ux3, &uy1, &uy2, &uy3, lower_deg, up_position, polarization, c_lat); getpoints(&lx1, &lx2, &lx3, &ly1, &ly2, &ly3, lower_deg, low_position, polarization, c_lat); /* DO NOT rotate the lower_deg reference points to the LEFT to be on the same axis */ /* make set of points into two vectors for each position */ Lx1 = lx2 - lx1; Ly1 = ly2 - ly1; Lx2 = lx3 - lx1; Ly2 = ly3 - ly1; Ux1 = ux2 - ux1; Uy1 = uy2 - uy1; Ux2 = ux3 - ux1; Uy2 = uy3 - uy1; /* calculate the interpolated vector components between orbit_positions for the lower_deg*/ lowx1 = (Ux1-Lx1)*orbit_scale + Lx1; lowy1 = (Uy1-Ly1)*orbit_scale + Ly1; lowx2 = (Ux2-Lx2)*orbit_scale + Lx2; lowy2 = (Uy2-Ly2)*orbit_scale + Ly2; /* Interpolate between the azimuth angles of the two orbit positions */ /* calculate the azimuth scale factor */ azimuth_scale = (antenna_azimuth - (float)(lower_deg*10))/10.0 ; /* linearly interpolate for antenna azimuth between the two interpolated positions */ upx1 = (upx1 - lowx1)*azimuth_scale + lowx1; upy1 = (upy1 - lowy1)*azimuth_scale + lowy1; upx2 = (upx2 - lowx2)*azimuth_scale + lowx2; upy2 = (upy2 - lowy2)*azimuth_scale + lowy2; /* determine which table input is closest in both orbit position and azimuth angle*/ if(orbit_scale < 0.5){ /*use the lower orbit position*/ if(azimuth_scale < 0.5){ /*use the lower azimuth angle*/ /* calculate the transform matrix coefficients between the lower_deg lower-orbit position vectors and the interpolated vectors */ d = (Lx2*upy1 - Lx1*upy2)/(Lx2*Ly1 - Lx1*Ly2); b = (Lx2*upx1 - Lx1*upx2)/(Lx2*Ly1 - Lx1*Ly2); c = (upy1 - d*Ly1)/Lx1; a = (upx1 - b*Ly1)/Lx1; /* calculate the needed rotation ammount */ rotate_deg = (antenna_azimuth - (float)(lower_deg*10)); /* get the lower set of del_lat and del_lon */ if(polarization == 0/* Hpol */){ lower_lat = &Hdel_lat[mmm*nnn*lower_deg + mmm*nnn*num_looks*low_position]; lower_lon = &Hdel_lon[mmm*nnn*lower_deg + mmm*nnn*num_looks*low_position]; } else{ lower_lat = &Vdel_lat[mmm*nnn*lower_deg + mmm*nnn*num_looks*low_position]; lower_lon = &Vdel_lon[mmm*nnn*lower_deg + mmm*nnn*num_looks*low_position]; } } else{ /* use the upper azimuth position*/ /* calculate the transform matrix coefficients between the lower_deg lower-orbit position vectors and the interpolated vectors */ d = (Lrx2*upy1 - Lrx1*upy2)/(Lrx2*Lry1 - Lrx1*Lry2); b = (Lrx2*upx1 - Lrx1*upx2)/(Lrx2*Lry1 - Lrx1*Lry2); c = (upy1 - d*Lry1)/Lrx1; a = (upx1 - b*Lry1)/Lrx1; /* calculate the needed rotation ammount */ rotate_deg = (antenna_azimuth - (float)(upper_deg*10)); /* get the lower set of del_lat and del_lon */ if(polarization == 0/* Hpol */){ lower_lat = &Hdel_lat[mmm*nnn*upper_deg + mmm*nnn*num_looks*low_position]; lower_lon = &Hdel_lon[mmm*nnn*upper_deg + mmm*nnn*num_looks*low_position]; } else{ lower_lat = &Vdel_lat[mmm*nnn*upper_deg + mmm*nnn*num_looks*low_position]; lower_lon = &Vdel_lon[mmm*nnn*upper_deg + mmm*nnn*num_looks*low_position]; } } } else{/* use the upper orbit position */ if(azimuth_scale < 0.5){ /*use the lower azimuth angle*/ /* calculate the transform matrix coefficients between the lower_deg lower-orbit position vectors and the interpolated vectors */ d = (Ux2*upy1 - Ux1*upy2)/(Ux2*Uy1 - Ux1*Uy2); b = (Ux2*upx1 - Ux1*upx2)/(Ux2*Uy1 - Ux1*Uy2); c = (upy1 - d*Uy1)/Ux1; a = (upx1 - b*Uy1)/Ux1; /* calculate the needed rotation ammount */ rotate_deg = (antenna_azimuth - (float)(lower_deg*10)); /* get the lower set of del_lat and del_lon */ if(polarization == 0/* Hpol */){ lower_lat = &Hdel_lat[mmm*nnn*lower_deg + mmm*nnn*num_looks*up_position]; lower_lon = &Hdel_lon[mmm*nnn*lower_deg + mmm*nnn*num_looks*up_position]; } else{ lower_lat = &Vdel_lat[mmm*nnn*lower_deg + mmm*nnn*num_looks*up_position]; lower_lon = &Vdel_lon[mmm*nnn*lower_deg + mmm*nnn*num_looks*up_position]; } } else{ /* use the upper azimuth position*/ /* calculate the transform matrix coefficients between the lower_deg lower-orbit position vectors and the interpolated vectors */ d = (Urx2*upy1 - Urx1*upy2)/(Urx2*Ury1 - Urx1*Ury2); b = (Urx2*upx1 - Urx1*upx2)/(Urx2*Ury1 - Urx1*Ury2); c = (upy1 - d*Ury1)/Urx1; a = (upx1 - b*Ury1)/Urx1; /* calculate the needed rotation ammount */ rotate_deg = (antenna_azimuth - (float)(upper_deg*10)); /* get the lower set of del_lat and del_lon */ if(polarization == 0/* Hpol */){ lower_lat = &Hdel_lat[mmm*nnn*upper_deg + mmm*nnn*num_looks*up_position]; lower_lon = &Hdel_lon[mmm*nnn*upper_deg + mmm*nnn*num_looks*up_position]; } else{ lower_lat = &Vdel_lat[mmm*nnn*upper_deg + mmm*nnn*num_looks*up_position]; lower_lon = &Vdel_lon[mmm*nnn*upper_deg + mmm*nnn*num_looks*up_position]; } } } /* rotate and transform the lat lon fields */ Re = (1 - Kflat*sin(c_lat*deg2rad)*sin(c_lat*deg2rad))*Ra; for( i = 0; i < mmm*nnn; i++){ /* change from a lat lon to a rectilinear grid */ Rphi = Re*cos((c_lat + lower_lat[i] )*deg2rad); x = Rphi*sin(lower_lon[i]*deg2rad); y = Re * sin(lower_lat[i]*deg2rad) + Rphi*(1-cos(lower_lon[i]*deg2rad))*sin(c_lat*deg2rad); // Transform points t_x = a*x + b*y; t_y = c*x + d*y; // Rotate right r_x = t_x * cos( rotate_deg*deg2rad) + t_y * sin(rotate_deg*deg2rad); r_y = t_x *-sin( rotate_deg*deg2rad) + t_y * cos(rotate_deg*deg2rad); // convert back to lat pulse_del_lon[i] = asin(r_x/Rphi)*rad2deg; pulse_del_lat[i] = asin( (r_y-(1-cos(pulse_del_lon[i]*deg2rad))*sin(c_lat*deg2rad)*Rphi)/Re)*rad2deg; } /* determine which table input is closest in both orbit position and azimuth angle and get the response function*/ if(orbit_scale < 0.5){ /*use the lower orbit position*/ if(azimuth_scale < 0.5){ /*use the lower azimuth angle*/ getresponse( lower_deg, low_position, polarization); } else{ /* use the upper azimuth position*/ getresponse( upper_deg, low_position, polarization); } } else{/* use the upper orbit position */ if(azimuth_scale < 0.5){ /*use the lower azimuth angle*/ getresponse( lower_deg, up_position, polarization); } else{ /* use the upper azimuth position*/ getresponse( upper_deg, up_position, polarization); } } } /* Coastal: get equivalent points from response data and return the points in rectilinear space */ void getpoints(float *x1, float *x2, float *x3, float *y1, float *y2, float *y3, int look, int position , int polarization, float c_lat){ float lat1, lon1, lat2, lon2, lat3, lon3; float Re, Rphi, Kflat = 1/298.257; float Ra = 6378.1363; float deg2rad = 3.14159 /180; float rad2deg = 180/(3.14159); /* polarization == 0 (Hpol), == 1(Vpol) */ /* get the equivalent points to be used to get the scale and rotate factors*/ if(polarization == 0/* Hpol */){ /* get the significant points in latitude and longitude */ lat1 = Hdel_lat[mmm*nnn*look + mmm*nnn*num_looks*position + (Hpoints[6*num_looks*position + 6*look +1]-1)*mmm + (Hpoints[6*num_looks*position + 6*look +0]-1) ] ; lon1 = Hdel_lon[mmm*nnn*look + mmm*nnn*num_looks*position + (Hpoints[6*num_looks*position + 6*look +1]-1)*mmm + (Hpoints[6*num_looks*position + 6*look +0]-1) ] ; lat2 = Hdel_lat[mmm*nnn*look + mmm*nnn*num_looks*position + (Hpoints[6*num_looks*position + 6*look +3]-1)*mmm + (Hpoints[6*num_looks*position + 6*look +2]-1) ] ; lon2 = Hdel_lon[mmm*nnn*look + mmm*nnn*num_looks*position + (Hpoints[6*num_looks*position + 6*look +3]-1)*mmm + (Hpoints[6*num_looks*position + 6*look +2]-1) ] ; lat3 = Hdel_lat[mmm*nnn*look + mmm*nnn*num_looks*position + (Hpoints[6*num_looks*position + 6*look +5]-1)*mmm + (Hpoints[6*num_looks*position + 6*look +4]-1) ] ; lon3 = Hdel_lon[mmm*nnn*look + mmm*nnn*num_looks*position + (Hpoints[6*num_looks*position + 6*look +5]-1)*mmm + (Hpoints[6*num_looks*position + 6*look +4]-1) ] ; } else{ lat1 = Vdel_lat[mmm*nnn*look + mmm*nnn*num_looks*position + (Vpoints[6*num_looks*position + 6*look +1]-1)*mmm + (Vpoints[6*num_looks*position + 6*look +0]-1) ] ; lon1 = Vdel_lon[mmm*nnn*look + mmm*nnn*num_looks*position + (Vpoints[6*num_looks*position + 6*look +1]-1)*mmm + (Vpoints[6*num_looks*position + 6*look +0]-1) ] ; lat2 = Vdel_lat[mmm*nnn*look + mmm*nnn*num_looks*position + (Vpoints[6*num_looks*position + 6*look +3]-1)*mmm + (Vpoints[6*num_looks*position + 6*look +2]-1) ] ; lon2 = Vdel_lon[mmm*nnn*look + mmm*nnn*num_looks*position + (Vpoints[6*num_looks*position + 6*look +3]-1)*mmm + (Vpoints[6*num_looks*position + 6*look +2]-1) ] ; lat3 = Vdel_lat[mmm*nnn*look + mmm*nnn*num_looks*position + (Vpoints[6*num_looks*position + 6*look +5]-1)*mmm + (Vpoints[6*num_looks*position + 6*look +4]-1) ] ; lon3 = Vdel_lon[mmm*nnn*look + mmm*nnn*num_looks*position + (Vpoints[6*num_looks*position + 6*look +5]-1)*mmm + (Vpoints[6*num_looks*position + 6*look +4]-1) ] ; } /* convert the points to x,y in km */ Re = (1 - Kflat*sin(c_lat*deg2rad)*sin(c_lat*deg2rad))*Ra; Rphi = Re*cos((c_lat + lat1 )*deg2rad); *x1 = Rphi*sin(lon1*deg2rad); *y1 = Re * sin(lat1*deg2rad) + Rphi*(1-cos(lon1*deg2rad))*sin(c_lat*deg2rad); Rphi = Re*cos((c_lat + lat2 )*deg2rad); *x2 = Rphi*sin(lon2*deg2rad); *y2 = Re * sin(lat2*deg2rad) + Rphi*(1-cos(lon2*deg2rad))*sin(c_lat*deg2rad); Rphi = Re*cos((c_lat + lat3 )*deg2rad); *x3 = Rphi*sin(lon3*deg2rad); *y3 = Re * sin(lat3*deg2rad) + Rphi*(1-cos(lon3*deg2rad))*sin(c_lat*deg2rad); } /* Coastal: set the global pointer to the response */ void getresponse( int look, int position , int polarization){ float *resp; /* polarization == 0 (Hpol), == 1(Vpol) */ /* pulse response will point to the correct segment of memory now */ if(polarization == 0 /* Hpol*/ ){ pulse_response = &Hresponse[mmm*nnn*8*look + mmm*nnn*8*num_looks*position]; pulse_cen_points = &Hcen_points[16*look + 16*num_looks*position]; } else{ pulse_response = &Vresponse[mmm*nnn*8*look + mmm*nnn*8*num_looks*position]; pulse_cen_points = &Vcen_points[16*look + 16*num_looks*position]; } } /* int temp, mm, nn,i,j,k; int num_looks, num_positions; int inbuffer [6]; int *points; float *del_lat, *del_lon, *response;*/ /* Coastal: Load the arrays of the response functions */ void loadresponse(){ /* note that the pointers that are passed to the function are the addresses of the global pointers...*/ /* to correctly use them in the subroutine you need to dereference each array pointer */ int inbuffer[2], temp, i,j,k; int polarization; for(polarization = 0; polarization < 2; polarization++){ if(polarization == 0){ resp_file = fopen(Hslicefile, "rb"); if (resp_file==NULL) { printf("*** ERROR: could not open response table file %s\n",Hslicefile); return; } } else{ resp_file = fopen(Vslicefile, "rb"); if (resp_file==NULL) { printf("*** ERROR: could not open response table file %s\n",Vslicefile); return; } } if(polarization == 0){/* Hpol */ fread( &temp , sizeof(int), 1, resp_file); fread( inbuffer , sizeof(int), temp/4, resp_file); num_positions = inbuffer[0]; num_looks = inbuffer[1]; fread( &temp , sizeof(int), 1, resp_file); /* for looks after the very first */ for(k = 0; k < num_positions; k++){ for(j = 0; j < num_looks; j++){ /* */ /* read dimensions*/ fread( &temp , sizeof(int), 1, resp_file); fread( inbuffer , sizeof(int), temp/4, resp_file); mmm = inbuffer[0]; nnn = inbuffer[1]; fread( &temp , sizeof(int), 1, resp_file); if(k == 0 && j == 0){/* first load allocate the memory */ /* allocate all needed memory */ Hdel_lat = (float *)malloc(mmm*nnn*sizeof(float)*num_looks*num_positions); Hdel_lon = (float *)malloc(mmm*nnn*sizeof(float)*num_looks*num_positions); Hresponse = (float *)malloc(8*mmm*nnn*sizeof(float)*num_looks*num_positions); Hpoints = (int *)malloc(6*sizeof(int)*num_looks*num_positions); Hcen_points = (int *)malloc(16*sizeof(int)*num_looks*num_positions); } /* read latitude data*/ fread( &temp , sizeof(int), 1, resp_file); /* load lat lon into later sections of del_lat and del_lon arrays */ fread( &Hdel_lat[mmm*nnn*j + mmm*nnn*num_looks*k] , sizeof(float), temp/4, resp_file); fread( &temp , sizeof(int), 1, resp_file); /* read longitude data*/ fread( &temp , sizeof(int), 1, resp_file); fread( &Hdel_lon[mmm*nnn*j + mmm*nnn*num_looks*k] , sizeof(float), temp/4, resp_file); fread( &temp , sizeof(int), 1, resp_file); /* read response data*/ fread( &temp , sizeof(int), 1, resp_file); fread( &Hresponse[mmm*nnn*8*j + mmm*nnn*8*num_looks*k] , sizeof(float), temp/4, resp_file); fread( &temp , sizeof(int), 1, resp_file); /* read significant points*/ fread( &temp , sizeof(int), 1, resp_file); fread( &Hpoints[6*j + 6*num_looks*k] , sizeof(int), 6, resp_file); fread( &temp , sizeof(int), 1, resp_file); /* read the center points */ fread( &temp , sizeof(int), 1, resp_file); fread( &Hcen_points[16*j + 16*num_looks*k] , sizeof(int), 16, resp_file); fread( &temp , sizeof(int), 1, resp_file); } } } else{/* Vertical polarization */ fread( &temp , sizeof(int), 1, resp_file); fread( inbuffer , sizeof(int), temp/4, resp_file); num_positions = inbuffer[0]; num_looks = inbuffer[1]; fread( &temp , sizeof(int), 1, resp_file); /* read all data */ /* for looks after the very first */ for(k = 0; k < num_positions; k++){ for(j = 0; j < num_looks; j++){ /* */ /* read dimensions*/ fread( &temp , sizeof(int), 1, resp_file); fread( inbuffer , sizeof(int), temp/4, resp_file); mmm = inbuffer[0]; nnn = inbuffer[1]; fread( &temp , sizeof(int), 1, resp_file); if(k == 0 && j == 0){ /* allocate all needed memory */ Vdel_lat = (float *)malloc(mmm*nnn*sizeof(float)*num_looks*num_positions); Vdel_lon = (float *)malloc(mmm*nnn*sizeof(float)*num_looks*num_positions); Vresponse = (float *)malloc(8*mmm*nnn*sizeof(float)*num_looks*num_positions); Vpoints = (int *)malloc(6*sizeof(int)*num_looks*num_positions); Vcen_points = (int *)malloc(16*sizeof(int)*num_looks*num_positions); } /* read latitude data*/ fread( &temp , sizeof(int), 1, resp_file); /* load lat lon into later sections of del_lat and del_lon arrays */ fread( &Vdel_lat[mmm*nnn*j + mmm*nnn*num_looks*k] , sizeof(float), temp/4, resp_file); fread( &temp , sizeof(int), 1, resp_file); /* read longitude data*/ fread( &temp , sizeof(int), 1, resp_file); fread( &Vdel_lon[mmm*nnn*j + mmm*nnn*num_looks*k] , sizeof(float), temp/4, resp_file); fread( &temp , sizeof(int), 1, resp_file); /* read response data*/ fread( &temp , sizeof(int), 1, resp_file); fread( &Vresponse[mmm*nnn*8*j + mmm*nnn*8*num_looks*k] , sizeof(float), temp/4, resp_file); fread( &temp , sizeof(int), 1, resp_file); /* read significant points*/ fread( &temp , sizeof(int), 1, resp_file); fread( &Vpoints[6*j + 6*num_looks*k] , sizeof(int), 6, resp_file); fread( &temp , sizeof(int), 1, resp_file); /* read the center points */ fread( &temp , sizeof(int), 1, resp_file); fread( &Vcen_points[16*j + 16*num_looks*k] , sizeof(int), 16, resp_file); fread( &temp , sizeof(int), 1, resp_file); } } } fclose(resp_file); } /* Allocate memory needed to hold the corrected response for any pulse*/ pulse_del_lat =(float *)malloc(mmm*nnn*sizeof(float)); pulse_del_lon =(float *)malloc(mmm*nnn*sizeof(float)); pulse_points = (int *)malloc(6*sizeof(int)); pulse_cen_points = (int *)malloc(16*sizeof(int)); } /* free dynamically allocated arrays*/ /* void free_response(){ */ /* free(Vdel_lat); */ /* free(Vdel_lon); */ /* free(Vresponse); */ /* free(Vpoints); */ /* free(Hdel_lat); */ /* free(Hdel_lon); */ /* free(Hresponse); */ /* free(Hpoints); */ /* free(pulse_del_lat); */ /* free(pulse_del_lon); */ /* free(pulse_points); */ /* } */