#include <cstdio> #include <cstdlib> #include <cmath> #include "icessmis.h" #include "icegrids.h" #include "grid_math.h" #include "team2.C" #include "pole_fill.C" /*C$$$ MAIN PROGRAM DOCUMENTATION BLOCK */ /*C */ /*C MAIN PROGRAM: SSMI-S CONSTRUCT SEA ICE CONCENTRATION GRIDS */ /*C PRGMMR: ROBERT GRUMBINE ORG: W/NP21 DATE: 97-06-25 */ /*C */ /*C ABSTRACT: READ IN SSMI-S DATA AND CONSTRUCT SEA ICE CONCENTRATION */ /*C GRIDS USING THE NASA TEAM ALGORITHM. */ /*C */ /*C PROGRAM HISTORY LOG: */ /*C 97-06-25 ROBERT GRUMBINE */ /*C */ /*C USAGE: */ /*C FILES ARE HANDLED BY ORDER IN THE INVOKING COMMAND LINE */ /*C INPUT FILES: */ /*C FILE 1 - File with names of SSMI-S data files. UNFORMATTED BINARY GENERATED */ /*C - BY SSMI-S.bufr.X */ /*C FILE 2 - NORTHERN HEMISPHERE LAND MASK */ /*C FILE 3 - SOUTHERN HEMISPHERE LAND MASK */ /*C OUTPUT FILES: */ /*C FILE 4 - NORTHERN HEMISPHERE GRIDDED SSMI-S data */ /*C FILE 5 - SOUTHERN HEMISPHERE GRIDDED SSMI-S data */ /*C FILE 6 - UNMASKED NORTHERN HEMISPHERE ice CONCENTRATIONS */ /*C FILE 7 - UNMASKED SOUTHERN HEMISPHERE ice CONCENTRATIONS */ /* Additional Arguments: */ /* 8 - Julian Day */ /* 9 - Year */ /* 10 - Satellite number (F11 = 244 .. F16 = 249 */ /*C */ /*C SUBPROGRAMS CALLED: */ /*C UNIQUE: process_bufr, process_short_bufr, check_bufr, */ /*C zero_bufr, newfilt, pole_fill */ /*C getfld, ice_zero, ice_add_data, ice_add_bufr, */ /*C ice_avg_data, nasa_team, nasa_team2 */ /*C LIBRARY: */ /*C OMB: mapll, mapxy */ /*C */ /*C EXIT STATES: */ /*C COND = 0 - SUCCESSFUL RUN */ /*C = -1 - ERROR OCCURRED. MESSAGE WRITTED TO STDOUT */ /*C */ /*C REMARKS: */ /*C */ /*C ATTRIBUTES: */ /*C LANGUAGE: STANDARD C++ */ /*C MACHINE: ANY */ /*C */ /*C$$$ */ /* Main program for decoding the ssmis orbital records. */ /* Derived from ssmi program of 1997-2006 */ /* ************************************************ */ /* The following are for data qc and error tracking */ /* for process_bufr */ int err_19h_range = 0; int err_19v_range = 0; int err_22v_range = 0; int err_37v_range = 0; int err_37h_range = 0; int err_92v_range = 0; int err_92h_range = 0; int err_150h_range = 0; int err_19_polar = 0; int err_37_polar = 0; int err_92_polar = 0; int err_lat = 0; int err_lon = 0; /* For the algorithm */ int bad_low = 0; int bad_high = 0; int crop_low = 0; int filt37 = 0; int filt22 = 0; /* For filtration */ int efilt_37_n = 0; int efilt_37_s = 0; int efilt_22_n = 0; int efilt_22_s = 0; // Generic file for conducting other training/testing: FILE *tester; /* ************************************************ */ //float fmax(float x, float y); //float fmax(float x, float y) { // if (x > y) { // return x; // } // return y; //} extern float hires(ssmis *map, int nx, int ny, int polei, int polej, unsigned char *mask) ; /* ************************************************ */ int main(int argc, char *argv[]) { ssmisupt a_bufr_line; FILE *fin, *input, *outn, *outs, *nland, *sland; FILE *nraw, *sraw; char fname[80]; ssmis *south, *north; ssmis_tmp *north_tmp, *south_tmp; float rfld_n[NY_NORTH][NX_NORTH]; float rfld_s[NY_SOUTH][NX_SOUTH]; unsigned char fld_n[NY_NORTH][NX_NORTH]; unsigned char fld_s[NY_SOUTH][NX_SOUTH]; int north_pts, south_pts; int i, j, nerrs, nproc; int jday, satno, ref_year; int lines_used = 0, err_spots = 0; // For team2 ssmis_team2_tables arctic, antarctic; // Team2 initialization: arctic.tbmfy.resize(n_atm, n_tb); arctic.tbmow.resize(n_atm, n_tb); arctic.tbmcc.resize(n_atm, n_tb); arctic.tbmthin.resize(n_atm, n_tb); arctic.pole = 'n'; arctic_tables(arctic); lookuptable(arctic); antarctic.tbmfy.resize(n_atm, n_tb); antarctic.tbmow.resize(n_atm, n_tb); antarctic.tbmcc.resize(n_atm, n_tb); antarctic.tbmthin.resize(n_atm, n_tb); antarctic.pole = 's'; antarctic_tables(antarctic); lookuptable(antarctic); // End Team2 initialization /* New the ssmis arrays */ north_tmp = new ssmis_tmp[NX_NORTH*NY_NORTH]; south_tmp = new ssmis_tmp[NX_SOUTH*NY_SOUTH]; north = new ssmis[NX_NORTH*NY_NORTH]; south = new ssmis[NX_SOUTH*NY_SOUTH]; if (north == (ssmis*) NULL || south == (ssmis*) NULL || south_tmp == (ssmis_tmp*) NULL || north_tmp == (ssmis_tmp*) NULL) { printf("Failed to new one of the ssmis fields!\n"); return -1; } /* Open data files */ fin = fopen(argv[1], "r"); nland = fopen(argv[2], "r"); sland = fopen(argv[3], "r"); outn = fopen(argv[4], "w"); outs = fopen(argv[5], "w"); nraw = fopen(argv[6], "w"); sraw = fopen(argv[7], "w"); if (fin == (FILE *) NULL || outn == (FILE *) NULL || outs == (FILE *) NULL || nland == (FILE *) NULL || sland == (FILE *) NULL || nraw == (FILE *) NULL || sraw == (FILE *) NULL) { printf("Failed to open one of the nine files!\n"); if (fin == (FILE *) NULL ) {printf(" - input file %s\n",argv[1]);} if (nland == (FILE *) NULL ) {printf(" - north land in %s\n", argv[2]);} if (sland == (FILE *) NULL ) {printf(" - south land in %s\n", argv[3]);} if (outn == (FILE *) NULL ) {printf(" - north out %s\n", argv[4]);} if (outs == (FILE *) NULL ) {printf(" - south out %s\n", argv[5]);} if (nraw == (FILE *) NULL ) {printf(" - north raw out %s\n", argv[6]);} if (sraw == (FILE *) NULL ) {printf(" - south raw out %s\n", argv[7]);} return -1; } #ifdef TESTER tester = fopen("tester","w"); #endif jday = atoi(argv[8]); ref_year = atoi(argv[9]); satno = atoi(argv[10]); printf("jday = %d, year = %d, satellite number %d\n",jday, ref_year, satno); printf("size of ssmisupt = %d\n",(int) sizeof(ssmisupt) ); fflush(stdout); /* Set internal constants */ nproc = NORBITS; north_pts = (NX_NORTH)*(NY_NORTH); south_pts = (NX_SOUTH)*(NY_SOUTH); ice_zero(north_tmp, south_tmp, north_pts, south_pts); /* Read in the land masks fld_n, fld_s */ fread(fld_n, sizeof(unsigned char), north_pts, nland); fread(fld_s, sizeof(unsigned char), south_pts, sland); while ( !feof(fin) ) { /* Read in name of next file to process */ fscanf(fin, "%s\n", fname); printf("Opened %s in orbital process loop\n", fname); fflush(stdout); input = fopen(fname, "r"); if ( input == (FILE *) NULL ) { printf("Failed to open %s from orbit process loop\n", fname); fflush(stdout); break; } /* Begin stripping out header */ fseek(input, (long)(0) , SEEK_SET); /* Begin stripping out data */ for (i = 0; i < nproc ; i++) while (!feof(input) ) { /* Add data only if data is within window -- Now assumed that the*/ /* date-time filtering has been done prior to input */ j = 0; while (!feof(input) ) { #ifdef VERBOSE2 printf("Process record number %d\n",j); fflush(stdout); #endif j++; fread(&a_bufr_line, sizeof(ssmisupt), 1, input); #ifdef VERBOSE2 printf("%d Satno = %d\n",j, a_bufr_line.satid); #endif if (a_bufr_line.satid == satno ) { nerrs = process_bufr(&a_bufr_line); lines_used += 1; err_spots += nerrs; #ifdef VERBOSE if (nerrs != 0) { printf("%2d errors in record %6d of orbit %2d\n",nerrs, j, i); fflush(stdout); } #endif ice_add_bufr(north_tmp, south_tmp, &a_bufr_line, arctic, antarctic); } } } /* end processing the orbits */ fclose(input); } /* End while - reading files */ printf("seaice_seaissmisu Done reading all input files\n"); fflush(stdout); /* Compute averaged tb and concentrations */ ice_avg_data(north_tmp, south_tmp, north, south, north_pts, south_pts, arctic, antarctic); #ifdef HIRES /* Now that we have a clean grid, work on the high resolution version */ hires(north, NX_NORTH, NY_NORTH, (int) (0.5+polei_NORTH), (int) (0.5+polej_NORTH), &fld_n[0][0]); hires(south, NX_SOUTH, NY_SOUTH, (int) (0.5+polei_SOUTH), (int) (0.5+polej_SOUTH), &fld_s[0][0]); ssmis_getfld(north, north_pts, &fld_n[0][0], &rfld_n[0][0], SSMIS_HIRES_CONC); ssmis_getfld(south, south_pts, &fld_s[0][0], &rfld_s[0][0], SSMIS_HIRES_CONC); #else ssmis_getfld(north, north_pts, &fld_n[0][0], &rfld_n[0][0], SSMIS_BAR_CONC); ssmis_getfld(south, south_pts, &fld_s[0][0], &rfld_s[0][0], SSMIS_BAR_CONC); #endif printf("past getting concentration fields in to fld grids \n"); fflush(stdout); /* Fill in pole for output purposes */ ssmis_pole_fill(&fld_n[0][0], 1); ssmis_pole_fill(&fld_s[0][0], 2); fwrite(fld_n, sizeof(unsigned char), north_pts, nraw); fwrite(fld_s, sizeof(unsigned char), south_pts, sraw); /* Mask out the land now */ ice_mask(north, south, north_pts, south_pts, &fld_n[0][0], &fld_s[0][0]); /* Write out full ssmis and concentration data files */ fwrite(north, sizeof(ssmis), north_pts, outn); fwrite(south, sizeof(ssmis), south_pts, outs); fclose(outn); fclose(outs); // From here on is doing some system monitoring/checking: /* Expt: Write out the error counts by type */ if (lines_used > 0) { printf("Total lores spots: %8d Error spots %8d fraction %6.4f\n", lines_used, err_spots, (float) err_spots / (float) (lines_used) ); printf("Error stats\n"); printf(" lat %d lon %d\n",err_lat, err_lon); printf(" range 19v %5d 19h %5d 22v %5d 37v %5d 37h %5d 92v %5d 92h %5d 150h %5d\n", err_19v_range, err_19h_range, err_22v_range, err_37v_range, err_37h_range, err_92v_range, err_92h_range, err_150h_range); printf(" polar 19GHz %5d 37GHz %5d 92GHz %6d \n", err_19_polar, err_37_polar, err_92_polar); printf(" bad_low %4d crop_low %4d bad high %4d\n",bad_low, crop_low, bad_high); printf(" filt37 %d filt22 %d efilt_37_n %6d efilt_37_s %6d\n", filt37, filt22, efilt_37_n, efilt_37_s); } else { printf("Did not find data for satellite %3d\n",satno ); } /* Expt: Write out the data count maps */ #ifdef VERBOSE outn = fopen("count.n", "w"); outs = fopen("count.s", "w"); getfld(north, north_pts, &fld_n[0][0], &rfld_n[0][0], COUNT); getfld(south, south_pts, &fld_s[0][0], &rfld_s[0][0], COUNT); fwrite(fld_n, sizeof(unsigned char), north_pts, outn); fwrite(fld_s, sizeof(unsigned char), south_pts, outs); fclose(outn); fclose(outs); #endif return 0; } #include "algorithm.C" #include "hires.C" #include "getfld.C" #include "icetools.C" #include "process_bufr.C" #include "filt.C"