(file) Return to sharp.c CVS log (file) (dir) Up to [Development] / JSOC / proj / sharp / apps

Diff for /JSOC/proj/sharp/apps/sharp.c between version 1.2 and 1.16

version 1.2, 2012/08/27 19:56:00 version 1.16, 2013/07/04 02:14:14
Line 19 
Line 19 
  *      Version:  *      Version:
  *              v0.0    Jul 02 2012  *              v0.0    Jul 02 2012
  *              v0.1    Jul 23 2012  *              v0.1    Jul 23 2012
    *              v0.2    Sep 04 2012
    *              v0.3    Dec 18 2012
    *              v0.4    Jan 02 2013
    *      v0.5    Jan 23 2012
  *  *
  *      Notes:  *      Notes:
  *              v0.0  *              v0.0
Line 31 
Line 35 
  *              Fixed char I/O thanks to Art  *              Fixed char I/O thanks to Art
  *              SW indices fixed  *              SW indices fixed
  *              Added doppler and continuum  *              Added doppler and continuum
    *              Added other keywords: HEADER (populated by cvs build version), DATE_B
    *              v0.3
    *              Fixed memory leakage of 0.15G per rec; denoted with "Dec 18"
    *              v0.4
    *              Took out convert_inplace(). Was causing all the images to be int
    *      v0.5
    *      Corrected ephemeris keywords, added argument mInfo for setKeys()
  *  *
  *      Example:  *      Example:
  *      sharp "mharp=hmi.Mharp_720s[1404][2012.02.20_10:00]" \  *      sharp "mharp=hmi.Mharp_720s[1404][2012.02.20_10:00]" \
Line 62 
Line 73 
 #include "errorprop.c" #include "errorprop.c"
 #include "sw_functions.c" #include "sw_functions.c"
  
   //#include <mkl.h> // Comment out mkl.h, which can only run on solar3
 #include <mkl_blas.h> #include <mkl_blas.h>
 #include <mkl_service.h> #include <mkl_service.h>
 #include <mkl_lapack.h> #include <mkl_lapack.h>
Line 80 
Line 92 
 // Nyqvist rate at disk center is 0.03 degree. Oversample above 0.015 degree // Nyqvist rate at disk center is 0.03 degree. Oversample above 0.015 degree
 #define NYQVIST         (0.015) #define NYQVIST         (0.015)
  
   // Some other things
 #ifndef MIN #ifndef MIN
 #define MIN(a,b) (((a)<(b)) ? (a) : (b)) #define MIN(a,b) (((a)<(b)) ? (a) : (b))
 #endif #endif
Line 107 
Line 120 
 #define INTERP                  0 #define INTERP                  0
 #define dpath    "/home/jsoc/cvs/Development/JSOC" #define dpath    "/home/jsoc/cvs/Development/JSOC"
  
   
 /* ========================================================================================================== */ /* ========================================================================================================== */
  
 // Space weather keywords // Space weather keywords
 struct swIndex { struct swIndex {
         float mean_vf;         float mean_vf;
       float count_mask;
         float absFlux;         float absFlux;
         float mean_hf;         float mean_hf;
         float mean_gamma;         float mean_gamma;
Line 131  struct swIndex {
Line 146  struct swIndex {
         float meanshear_angle;         float meanshear_angle;
         float area_w_shear_gt_45h;         float area_w_shear_gt_45h;
         float meanshear_angleh;         float meanshear_angleh;
       float mean_derivative_btotal_err;
       float mean_vf_err;
       float mean_gamma_err;
       float mean_derivative_bh_err;
       float mean_derivative_bz_err;
       float mean_jz_err;
       float us_i_err;
       float mean_alpha_err;
       float mean_ih_err;
       float total_us_ih_err;
       float total_abs_ih_err;
       float totaljz_err;
       float meanpot_err;
       float totpot_err;
       float meanshear_angle_err;
 }; };
  
 // Mapping method // Mapping method
Line 147  enum projection {
Line 177  enum projection {
         lambert         lambert
 }; };
  
 // Ephemeris  // WSC code
   char *wcsCode[] = {"CAR", "CAS", "MER", "CEA", "GLS", "TAN", "ARC", "STG",
           "SIN", "ZEA"};
   
   // Ephemeris information
 struct ephemeris { struct ephemeris {
         double disk_lonc, disk_latc;         double disk_lonc, disk_latc;
         double disk_xc, disk_yc;         double disk_xc, disk_yc;
Line 165  struct mapInfo {
Line 199  struct mapInfo {
         float *xi_out, *zeta_out;       // coordinate on full disk image to sample at         float *xi_out, *zeta_out;       // coordinate on full disk image to sample at
 }; };
  
   
 /* ========================================================================================================== */ /* ========================================================================================================== */
  
 /* Get all input data series */ /* Get all input data series */
Line 181  int getInputRS_aux(DRMS_RecordSet_t **in
Line 214  int getInputRS_aux(DRMS_RecordSet_t **in
 /* Find record from record set with given T_rec */ /* Find record from record set with given T_rec */
 int getInputRec_aux(DRMS_Record_t **inRec_ptr, DRMS_RecordSet_t *inRS, TIME trec); int getInputRec_aux(DRMS_Record_t **inRec_ptr, DRMS_RecordSet_t *inRS, TIME trec);
  
 // ===================  
   
 /* Create CEA record */ /* Create CEA record */
 int createCeaRecord(DRMS_Record_t *mharpRec, DRMS_Record_t *bharpRec, int createCeaRecord(DRMS_Record_t *mharpRec, DRMS_Record_t *bharpRec,
                                         DRMS_Record_t *dopRec, DRMS_Record_t *contRec,                                         DRMS_Record_t *dopRec, DRMS_Record_t *contRec,
Line 245  void computeSWIndex(struct swIndex *swKe
Line 276  void computeSWIndex(struct swIndex *swKe
 void setSWIndex(DRMS_Record_t *outRec, struct swIndex *swKeys_ptr); void setSWIndex(DRMS_Record_t *outRec, struct swIndex *swKeys_ptr);
  
 /* Set all keywords, no error checking for now */ /* Set all keywords, no error checking for now */
 void setKeys(DRMS_Record_t *outRec, DRMS_Record_t *inRec);  void setKeys(DRMS_Record_t *outRec, DRMS_Record_t *inRec, struct mapInfo *mInfo);
  
 // =================== // ===================
  
Line 287  char *CutSegs[] = {"magnetogram", "bitma
Line 318  char *CutSegs[] = {"magnetogram", "bitma
         "disambig", "conf_disambig"};         "disambig", "conf_disambig"};
 char *CEASegs[] = {"magnetogram", "bitmap", "Dopplergram", "continuum", "disambig", char *CEASegs[] = {"magnetogram", "bitmap", "Dopplergram", "continuum", "disambig",
                                         BR_SEG_CEA, BT_SEG_CEA, BP_SEG_CEA, BR_ERR_SEG_CEA, BT_ERR_SEG_CEA, BP_ERR_SEG_CEA};                                         BR_SEG_CEA, BT_SEG_CEA, BP_SEG_CEA, BR_ERR_SEG_CEA, BT_ERR_SEG_CEA, BP_ERR_SEG_CEA};
 /* ========================================================================================================== */  
   
  
   /* ========================================================================================================== */
  
 char *module_name = "sharp"; char *module_name = "sharp";
 char *version_id = "2012 Jul 02";  /* Version number */  char *version_id = "2013 Jun 26";  /* Version number */
   int seed;
  
 ModuleArgs_t module_args[] = ModuleArgs_t module_args[] =
 { {
Line 302  ModuleArgs_t module_args[] =
Line 333  ModuleArgs_t module_args[] =
         {ARG_STRING, "cont", kNotSpecified, "Input Continuum series."},         {ARG_STRING, "cont", kNotSpecified, "Input Continuum series."},
         {ARG_STRING, "sharp_cea", kNotSpecified, "Output Sharp CEA series."},         {ARG_STRING, "sharp_cea", kNotSpecified, "Output Sharp CEA series."},
         {ARG_STRING, "sharp_cut", kNotSpecified, "Output Sharp cutout series."},         {ARG_STRING, "sharp_cut", kNotSpecified, "Output Sharp cutout series."},
       {ARG_INT,    "seed", "987654", "Seed for the random number generator."},
         {ARG_END}         {ARG_END}
 }; };
  
 int DoIt(void) int DoIt(void)
 { {
       int errbufstat=setvbuf(stderr, NULL, _IONBF, BUFSIZ);
       int outbufstat=setvbuf(stdout, NULL, _IONBF, BUFSIZ);
  
         int status = DRMS_SUCCESS;         int status = DRMS_SUCCESS;
         int nrecs, irec;         int nrecs, irec;
Line 327  int DoIt(void)
Line 361  int DoIt(void)
         sharpCeaQuery = (char *) params_get_str(&cmdparams, "sharp_cea");         sharpCeaQuery = (char *) params_get_str(&cmdparams, "sharp_cea");
         sharpCutQuery = (char *) params_get_str(&cmdparams, "sharp_cut");         sharpCutQuery = (char *) params_get_str(&cmdparams, "sharp_cut");
  
       seed = params_get_int(&cmdparams, "seed");
   
         /* Get input data, check everything */         /* Get input data, check everything */
  
         if (getInputRS(&mharpRS, &bharpRS, mharpQuery, bharpQuery))         if (getInputRS(&mharpRS, &bharpRS, mharpQuery, bharpQuery))
Line 403  int DoIt(void)
Line 439  int DoIt(void)
  
         } // irec         } // irec
  
   
         drms_close_records(mharpRS, DRMS_FREE_RECORD);         drms_close_records(mharpRS, DRMS_FREE_RECORD);
         drms_close_records(bharpRS, DRMS_FREE_RECORD);         drms_close_records(bharpRS, DRMS_FREE_RECORD);
           drms_close_records(dopRS, DRMS_FREE_RECORD);                            // Dec 18 2012
           drms_close_records(contRS, DRMS_FREE_RECORD);                           // Dec 18 2012
  
         return 0;         return 0;
  
Line 564  int createCeaRecord(DRMS_Record_t *mharp
Line 603  int createCeaRecord(DRMS_Record_t *mharp
  
         // Get ephemeris         // Get ephemeris
  
         if (getEphemeris(mharpRec, &(mInfo.ephem))) return 1;          if (getEphemeris(mharpRec, &(mInfo.ephem))) {
                   SHOW("CEA: get ephemeris error\n");
                   return 1;
           }
  
         // Find position         // Find position
  
         if (findPosition(mharpRec, &mInfo)) return 1;          if (findPosition(mharpRec, &mInfo)) {
                   SHOW("CEA: find position error\n");
                   return 1;
           }
  
         // Create xi_out, zeta_out array in mInfo:         // Create xi_out, zeta_out array in mInfo:
         // Coordinates to sample in original full disk image         // Coordinates to sample in original full disk image
Line 584  int createCeaRecord(DRMS_Record_t *mharp
Line 629  int createCeaRecord(DRMS_Record_t *mharp
  
         // Mapping single segment: Mharp, etc.         // Mapping single segment: Mharp, etc.
  
         if (mapScaler(sharpRec, mharpRec, mharpRec, &mInfo, "magnetogram")) return 1;          if (mapScaler(sharpRec, mharpRec, mharpRec, &mInfo, "magnetogram")) {
         if (mapScaler(sharpRec, mharpRec, mharpRec, &mInfo, "bitmap")) return 1;                  SHOW("CEA: mapping magnetogram error\n");
                   return 1;
           }
           if (mapScaler(sharpRec, mharpRec, mharpRec, &mInfo, "bitmap")) {
                   SHOW("CEA: mapping bitmap error\n");
                   return 1;
           }
         printf("Magnetogram mapping done.\n");         printf("Magnetogram mapping done.\n");
           if (mapScaler(sharpRec, dopRec, mharpRec, &mInfo, "Dopplergram")) {
         if (mapScaler(sharpRec, dopRec, mharpRec, &mInfo, "Dopplergram")) return 1;                  SHOW("CEA: mapping dopplergram error\n");
                   return 1;
           }
         printf("Dopplergram mapping done.\n");         printf("Dopplergram mapping done.\n");
  
         if (mapScaler(sharpRec, contRec, mharpRec, &mInfo, "continuum")) return 1;          if (mapScaler(sharpRec, contRec, mharpRec, &mInfo, "continuum")) {
                   SHOW("CEA: mapping continuum error\n");
                   return 1;
           }
         printf("Intensitygram mapping done.\n");         printf("Intensitygram mapping done.\n");
  
         if (mapScaler(sharpRec, bharpRec, mharpRec, &mInfo, "conf_disambig")) return 1;          if (mapScaler(sharpRec, bharpRec, mharpRec, &mInfo, "conf_disambig")) {
                   SHOW("CEA: mapping conf_disambig error\n");
                   return 1;
           }
         printf("Conf disambig mapping done.\n");         printf("Conf disambig mapping done.\n");
  
         // Mapping vector B         // Mapping vector B
  
         if (mapVectorB(sharpRec, bharpRec, &mInfo)) return 1;          if (mapVectorB(sharpRec, bharpRec, &mInfo)) {
                   SHOW("CEA: mapping vector B error\n");
                   return 1;
           }
         printf("Vector B mapping done.\n");         printf("Vector B mapping done.\n");
  
         // Mapping vector B errors         // Mapping vector B errors
  
         if (mapVectorBErr(sharpRec, bharpRec, &mInfo)) return 1;          if (mapVectorBErr(sharpRec, bharpRec, &mInfo)) {
                   SHOW("CEA: mapping vector B uncertainty error\n");
                   return 1;
           }
         printf("Vector B error done.\n");         printf("Vector B error done.\n");
  
         // Keywords & Links         // Keywords & Links
Line 617  int createCeaRecord(DRMS_Record_t *mharp
Line 682  int createCeaRecord(DRMS_Record_t *mharp
         DRMS_Link_t *bHarpLink = hcon_lookup_lower(&sharpRec->links, "BHARP");         DRMS_Link_t *bHarpLink = hcon_lookup_lower(&sharpRec->links, "BHARP");
         if (bHarpLink) drms_link_set("BHARP", sharpRec, bharpRec);         if (bHarpLink) drms_link_set("BHARP", sharpRec, bharpRec);
  
         setKeys(sharpRec, bharpRec);            // Set all other keywords           setKeys(sharpRec, bharpRec, &mInfo);            // Set all other keywords
           drms_copykey(sharpRec, mharpRec, "QUALITY");            // copied from los records
  
         // Space weather         // Space weather
  
Line 700  int mapScaler(DRMS_Record_t *sharpRec, D
Line 766  int mapScaler(DRMS_Record_t *sharpRec, D
         outSeg = drms_segment_lookup(sharpRec, segName);         outSeg = drms_segment_lookup(sharpRec, segName);
         if (!outSeg) return 1;         if (!outSeg) return 1;
  
         DRMS_Type_t arrayType = outSeg->info->type;      //  DRMS_Type_t arrayType = outSeg->info->type;
         DRMS_Array_t *outArray = drms_array_create(DRMS_TYPE_FLOAT, 2, dims, map, &status);         DRMS_Array_t *outArray = drms_array_create(DRMS_TYPE_FLOAT, 2, dims, map, &status);
         if (status) {if (inArray) drms_free_array(inArray); free(map); return 1;}         if (status) {if (inArray) drms_free_array(inArray); free(map); return 1;}
  
         // convert to needed data type         // convert to needed data type
  
         drms_array_convert_inplace(outSeg->info->type, 0, 1, outArray);      //  drms_array_convert_inplace(outSeg->info->type, 0, 1, outArray);         // Jan 02 2013
  
         outSeg->axis[0] = outArray->axis[0]; outSeg->axis[1] = outArray->axis[1];         outSeg->axis[0] = outArray->axis[0]; outSeg->axis[1] = outArray->axis[1];
         outArray->parent_segment = outSeg;      //  outArray->parent_segment = outSeg;
         outArray->israw = 0;            // always compressed         outArray->israw = 0;            // always compressed
         outArray->bzero = outSeg->bzero;         outArray->bzero = outSeg->bzero;
         outArray->bscale = outSeg->bscale;         outArray->bscale = outSeg->bscale;
Line 718  int mapScaler(DRMS_Record_t *sharpRec, D
Line 784  int mapScaler(DRMS_Record_t *sharpRec, D
         if (status) return 0;         if (status) return 0;
  
         if (inArray) drms_free_array(inArray);         if (inArray) drms_free_array(inArray);
           if ((xsz != FOURK) || (ysz != FOURK)) free(inData);                     // Dec 18 2012
         if (outArray) drms_free_array(outArray);         if (outArray) drms_free_array(outArray);
         return 0;         return 0;
  
Line 784  int mapVectorB(DRMS_Record_t *sharpRec,
Line 851  int mapVectorB(DRMS_Record_t *sharpRec,
                 outSeg = drms_segment_lookup(sharpRec, segName[iSeg]);                 outSeg = drms_segment_lookup(sharpRec, segName[iSeg]);
                 outArray = drms_array_create(DRMS_TYPE_FLOAT, 2, dims, data_prt[iSeg], &status);                 outArray = drms_array_create(DRMS_TYPE_FLOAT, 2, dims, data_prt[iSeg], &status);
                 outSeg->axis[0] = outArray->axis[0]; outSeg->axis[1] = outArray->axis[1];                 outSeg->axis[0] = outArray->axis[0]; outSeg->axis[1] = outArray->axis[1];
                 outArray->parent_segment = outSeg;          //              outArray->parent_segment = outSeg;
                 outArray->israw = 0;                 outArray->israw = 0;
                 outArray->bzero = outSeg->bzero;                 outArray->bzero = outSeg->bzero;
                 outArray->bscale = outSeg->bscale;                 outArray->bscale = outSeg->bscale;
Line 836  int mapVectorBErr(DRMS_Record_t *sharpRe
Line 903  int mapVectorBErr(DRMS_Record_t *sharpRe
                 outSeg = drms_segment_lookup(sharpRec, segName[iSeg]);                 outSeg = drms_segment_lookup(sharpRec, segName[iSeg]);
                 outArray = drms_array_create(DRMS_TYPE_FLOAT, 2, dims, data_prt[iSeg], &status);                 outArray = drms_array_create(DRMS_TYPE_FLOAT, 2, dims, data_prt[iSeg], &status);
                 outSeg->axis[0] = outArray->axis[0]; outSeg->axis[1] = outArray->axis[1];                 outSeg->axis[0] = outArray->axis[0]; outSeg->axis[1] = outArray->axis[1];
                 outArray->parent_segment = outSeg;          //              outArray->parent_segment = outSeg;
                 outArray->israw = 0;                 outArray->israw = 0;
                 outArray->bzero = outSeg->bzero;                 outArray->bzero = outSeg->bzero;
                 outArray->bscale = outSeg->bscale;                 outArray->bscale = outSeg->bscale;
Line 881  int findPosition(DRMS_Record_t *inRec, s
Line 948  int findPosition(DRMS_Record_t *inRec, s
         // We compute minlon & minlat then by         // We compute minlon & minlat then by
         // LONDTMIN(t) = LONDTMIN(t0) + (t - t0) * OMEGA_DT         // LONDTMIN(t) = LONDTMIN(t0) + (t - t0) * OMEGA_DT
  
         float psize = drms_getkey_float(inRec, "SIZE", &status);      //  float psize = drms_getkey_float(inRec, "SIZE", &status);
         if (psize != psize) {      //  if (psize != psize) {
                 TIME t0 = drms_getkey_time(inRec, "T_FRST", &status); if (status) return 1;  
       if (minlon != minlon || maxlon != maxlon) {         // check lons instead of SIZE
                   TIME t0 = drms_getkey_time(inRec, "T_FRST1", &status); if (status) return 1;                    // changed from T_FRST to T_FRST1, T_FRST may not exist
                 double omega = drms_getkey_double(inRec, "OMEGA_DT", &status); if (status) return 1;                 double omega = drms_getkey_double(inRec, "OMEGA_DT", &status); if (status) return 1;
                 char firstRecQuery[100], t0_str[100];                 char firstRecQuery[100], t0_str[100];
                 sprint_time(t0_str, t0, "TAI", 0);                 sprint_time(t0_str, t0, "TAI", 0);
Line 935  int getEphemeris(DRMS_Record_t *inRec, s
Line 1004  int getEphemeris(DRMS_Record_t *inRec, s
         float crpix1 = drms_getkey_float(inRec, "IMCRPIX1", &status);         float crpix1 = drms_getkey_float(inRec, "IMCRPIX1", &status);
         float crpix2 = drms_getkey_float(inRec, "IMCRPIX2", &status);         float crpix2 = drms_getkey_float(inRec, "IMCRPIX2", &status);
         float cdelt = drms_getkey_float(inRec, "CDELT1", &status);  // in arcsec, assumimg dx=dy         float cdelt = drms_getkey_float(inRec, "CDELT1", &status);  // in arcsec, assumimg dx=dy
         printf("cdelt=%f\n",cdelt);  
         ephem->disk_xc = PIX_X(0.0,0.0) - 1.0;          // Center of disk in pixel, starting at 0         ephem->disk_xc = PIX_X(0.0,0.0) - 1.0;          // Center of disk in pixel, starting at 0
         ephem->disk_yc = PIX_Y(0.0,0.0) - 1.0;         ephem->disk_yc = PIX_Y(0.0,0.0) - 1.0;
  
Line 1077  int performSampling(float *outData, floa
Line 1145  int performSampling(float *outData, floa
         // Rebinning, smoothing         // Rebinning, smoothing
  
         frebin(outData0, outData, ncol0, nrow0, mInfo->nbin, 1);                // Gaussian         frebin(outData0, outData, ncol0, nrow0, mInfo->nbin, 1);                // Gaussian
           free(outData0);         // Dec 18 2012
  
         //         //
  
Line 1454  int createCutRecord(DRMS_Record_t *mharp
Line 1523  int createCutRecord(DRMS_Record_t *mharp
                         break;                         break;
                 }                 }
         }         }
         if (iHarpSeg != nMharpSegs) return 1;           // if failed          if (iHarpSeg != nMharpSegs) {
                   SHOW("Cutout: segment number unmatch\n");
                   return 1;               // if failed
           }
         printf("Magnetogram cutout done.\n");         printf("Magnetogram cutout done.\n");
  
         // Cutout Doppler         // Cutout Doppler
Line 1495  int createCutRecord(DRMS_Record_t *mharp
Line 1567  int createCutRecord(DRMS_Record_t *mharp
         if (bHarpLink) drms_link_set("BHARP", sharpRec, bharpRec);         if (bHarpLink) drms_link_set("BHARP", sharpRec, bharpRec);
  
         setSWIndex(sharpRec, swKeys_ptr);       // Set space weather indices         setSWIndex(sharpRec, swKeys_ptr);       // Set space weather indices
         setKeys(sharpRec, bharpRec);            // Set all other keywords          setKeys(sharpRec, bharpRec, NULL);              // Set all other keywords, NULL specifies cutout
  
         // Stats         // Stats
  
Line 1556  int writeCutout(DRMS_Record_t *outRec, D
Line 1628  int writeCutout(DRMS_Record_t *outRec, D
  
 #if DISAMB_AZI #if DISAMB_AZI
         if (!strcmp(SegName, "azimuth")) {         if (!strcmp(SegName, "azimuth")) {
                 DRMS_Segment_t *disambSeg = drms_segment_lookup(inRec, "disambig");                  DRMS_Segment_t *disambSeg = NULL;
                   disambSeg = drms_segment_lookup(inRec, "disambig");
                 if (!disambSeg) {drms_free_array(cutoutArray); return 1;}                 if (!disambSeg) {drms_free_array(cutoutArray); return 1;}
                 DRMS_Array_t *disambArray;                 DRMS_Array_t *disambArray;
                 if (disambSeg->axis[0] == nx && disambSeg->axis[1] == ny) {                 if (disambSeg->axis[0] == nx && disambSeg->axis[1] == ny) {
                         disambArray = drms_segment_read(disambSeg, DRMS_TYPE_CHAR, &status);                         disambArray = drms_segment_read(disambSeg, DRMS_TYPE_CHAR, &status);
                         if (status) {drms_free_array(cutoutArray); return 1;}                         if (status) {drms_free_array(cutoutArray); return 1;}
                 } else {                 } else {
               drms_free_array(cutoutArray);
                         return 1;                         return 1;
                 }                 }
                 double *azimuth = (double *) cutoutArray->data;                 double *azimuth = (double *) cutoutArray->data;
Line 1578  int writeCutout(DRMS_Record_t *outRec, D
Line 1652  int writeCutout(DRMS_Record_t *outRec, D
  
         outSeg = drms_segment_lookup(outRec, SegName);         outSeg = drms_segment_lookup(outRec, SegName);
         if (!outSeg) return 1;         if (!outSeg) return 1;
         drms_array_convert_inplace(outSeg->info->type, 0, 1, cutoutArray);      //  drms_array_convert_inplace(outSeg->info->type, 0, 1, cutoutArray);      // Jan 02 2013
         outSeg->axis[0] = cutoutArray->axis[0];         outSeg->axis[0] = cutoutArray->axis[0];
         outSeg->axis[1] = cutoutArray->axis[1];         outSeg->axis[1] = cutoutArray->axis[1];
         cutoutArray->parent_segment = outSeg;      //  cutoutArray->parent_segment = outSeg;
         cutoutArray->israw = 0;         // always compressed         cutoutArray->israw = 0;         // always compressed
     cutoutArray->bzero = outSeg->bzero;     cutoutArray->bzero = outSeg->bzero;
     cutoutArray->bscale = outSeg->bscale;               // Same as inArray's     cutoutArray->bscale = outSeg->bscale;               // Same as inArray's
Line 1608  void computeSWIndex(struct swIndex *swKe
Line 1682  void computeSWIndex(struct swIndex *swKe
         int nx = mInfo->ncol, ny = mInfo->nrow;         int nx = mInfo->ncol, ny = mInfo->nrow;
         int nxny = nx * ny;         int nxny = nx * ny;
         int dims[2] = {nx, ny};         int dims[2] = {nx, ny};
   
         // Get bx, by, bz, mask         // Get bx, by, bz, mask
  
         // Use HARP (Turmon) bitmap as a threshold on spaceweather quantities         // Use HARP (Turmon) bitmap as a threshold on spaceweather quantities
         //DRMS_Segment_t *maskSeg = drms_segment_lookup(inRec, "bitmap");          DRMS_Segment_t *bitmaskSeg = drms_segment_lookup(inRec, "bitmap");
         //DRMS_Array_t *maskArray = drms_segment_read(maskSeg, DRMS_TYPE_INT, &status);          DRMS_Array_t *bitmaskArray = drms_segment_read(bitmaskSeg, DRMS_TYPE_INT, &status);
         //int *mask = (int *) maskArray->data;          // get the previously made mask array          int *bitmask = (int *) bitmaskArray->data;              // get the previously made mask array
  
         //Use conf_disambig map as a threshold on spaceweather quantities         //Use conf_disambig map as a threshold on spaceweather quantities
         DRMS_Segment_t *maskSeg = drms_segment_lookup(inRec, "conf_disambig");         DRMS_Segment_t *maskSeg = drms_segment_lookup(inRec, "conf_disambig");
Line 1633  void computeSWIndex(struct swIndex *swKe
Line 1708  void computeSWIndex(struct swIndex *swKe
         DRMS_Array_t *bzArray = drms_segment_read(bzSeg, DRMS_TYPE_FLOAT, &status);         DRMS_Array_t *bzArray = drms_segment_read(bzSeg, DRMS_TYPE_FLOAT, &status);
         float *bz = (float *) bzArray->data;            // bz         float *bz = (float *) bzArray->data;            // bz
  
         // Get emphemeris          DRMS_Segment_t *bz_errSeg = drms_segment_lookup(inRec, BR_ERR_SEG_CEA);
           DRMS_Array_t *bz_errArray = drms_segment_read(bz_errSeg, DRMS_TYPE_FLOAT, &status);
           float *bz_err = (float *) bz_errArray->data;            // bz_err
   
           DRMS_Segment_t *by_errSeg = drms_segment_lookup(inRec, BT_ERR_SEG_CEA);
           DRMS_Array_t *by_errArray = drms_segment_read(by_errSeg, DRMS_TYPE_FLOAT, &status);
           float *by_err = (float *) by_errArray->data;            // by_err
           //for (int i = 0; i < nxny; i++) by_err[i] *= -1;
   
           DRMS_Segment_t *bx_errSeg = drms_segment_lookup(inRec, BP_ERR_SEG_CEA);
           DRMS_Array_t *bx_errArray = drms_segment_read(bx_errSeg, DRMS_TYPE_FLOAT, &status);
           float *bx_err = (float *) bx_errArray->data;            // bx_err
  
         //float cdelt1_orig = drms_getkey_float(inRec, "CDELT1",   &status);          // Get emphemeris
         float cdelt1      = drms_getkey_float(inRec, "CDELT1",   &status);          float  cdelt1_orig = drms_getkey_float(inRec, "CDELT1",   &status);
         float dsun_obs    = drms_getkey_float(inRec, "DSUN_OBS",   &status);         float dsun_obs    = drms_getkey_float(inRec, "DSUN_OBS",   &status);
         double rsun_ref   = drms_getkey_double(inRec, "RSUN_REF", &status);         double rsun_ref   = drms_getkey_double(inRec, "RSUN_REF", &status);
         double rsun_obs   = drms_getkey_double(inRec, "RSUN_OBS", &status);         double rsun_obs   = drms_getkey_double(inRec, "RSUN_OBS", &status);
Line 1645  void computeSWIndex(struct swIndex *swKe
Line 1731  void computeSWIndex(struct swIndex *swKe
         float crpix1      = drms_getkey_float(inRec, "CRPIX1", &status);         float crpix1      = drms_getkey_float(inRec, "CRPIX1", &status);
         float crpix2      = drms_getkey_float(inRec, "CRPIX2", &status);         float crpix2      = drms_getkey_float(inRec, "CRPIX2", &status);
  
         //float cdelt1=( (rsun_ref*cdelt1_orig*PI/180.) / (dsun_obs) )*(180./PI)*(3600.); //convert cdelt1 from degrees to arcsec (approximately)          // convert cdelt1_orig from degrees to arcsec
           float cdelt1       = (atan((rsun_ref*cdelt1_orig*RADSINDEG)/(dsun_obs)))*(1/RADSINDEG)*(3600.);
         printf("cdelt1=%f\n",cdelt1);  
         printf("rsun_ref=%f\n",rsun_ref);  
         printf("rsun_obs=%f\n",rsun_obs);  
         printf("dsun_obs=%f\n",dsun_obs);  
  
         // Temp arrays         // Temp arrays
   
         float *bh = (float *) (malloc(nxny * sizeof(float)));         float *bh = (float *) (malloc(nxny * sizeof(float)));
         float *bt = (float *) (malloc(nxny * sizeof(float)));         float *bt = (float *) (malloc(nxny * sizeof(float)));
         float *jz = (float *) (malloc(nxny * sizeof(float)));         float *jz = (float *) (malloc(nxny * sizeof(float)));
           float *jz_smooth = (float *) (malloc(nxny * sizeof(float)));
         float *bpx = (float *) (malloc(nxny * sizeof(float)));         float *bpx = (float *) (malloc(nxny * sizeof(float)));
         float *bpy = (float *) (malloc(nxny * sizeof(float)));         float *bpy = (float *) (malloc(nxny * sizeof(float)));
         float *bpz = (float *) (malloc(nxny * sizeof(float)));         float *bpz = (float *) (malloc(nxny * sizeof(float)));
Line 1668  void computeSWIndex(struct swIndex *swKe
Line 1750  void computeSWIndex(struct swIndex *swKe
         float *dery_bh = (float *) (malloc(nxny * sizeof(float)));         float *dery_bh = (float *) (malloc(nxny * sizeof(float)));
         float *derx_bz = (float *) (malloc(nxny * sizeof(float)));         float *derx_bz = (float *) (malloc(nxny * sizeof(float)));
         float *dery_bz = (float *) (malloc(nxny * sizeof(float)));         float *dery_bz = (float *) (malloc(nxny * sizeof(float)));
           float *bt_err  = (float *) (malloc(nxny * sizeof(float)));
           float *bh_err  = (float *) (malloc(nxny * sizeof(float)));
           float *jz_err  = (float *) (malloc(nxny * sizeof(float)));
           float *jz_err_squared = (float *) (malloc(nxny * sizeof(float)));
           float *jz_err_squared_smooth = (float *) (malloc(nxny * sizeof(float)));
           float *jz_rms_err = (float *) (malloc(nxny * sizeof(float)));
           //spaceweather quantities computed
  
         // Compute  
  
         if (computeAbsFlux(bz, dims, &(swKeys_ptr->absFlux), &(swKeys_ptr->mean_vf),          if (computeAbsFlux(bz_err, bz , dims, &(swKeys_ptr->absFlux), &(swKeys_ptr->mean_vf),  &(swKeys_ptr->mean_vf_err),
                                            mask, cdelt1, rsun_ref, rsun_obs)){                         &(swKeys_ptr->count_mask), mask, bitmask, cdelt1, rsun_ref, rsun_obs))
       {
                 swKeys_ptr->absFlux = DRMS_MISSING_FLOAT;               // If fail, fill in NaN                 swKeys_ptr->absFlux = DRMS_MISSING_FLOAT;               // If fail, fill in NaN
                 swKeys_ptr->mean_vf = DRMS_MISSING_FLOAT;                 swKeys_ptr->mean_vf = DRMS_MISSING_FLOAT;
           swKeys_ptr->mean_vf_err = DRMS_MISSING_FLOAT;
           swKeys_ptr->count_mask  = DRMS_MISSING_INT;
         }         }
  
         for (int i = 0; i < nxny; i++) bpz[i] = bz[i];         for (int i = 0; i < nxny; i++) bpz[i] = bz[i];
         greenpot(bpx, bpy, bpz, nx, ny);         greenpot(bpx, bpy, bpz, nx, ny);
  
         computeBh(bx, by, bz, bh, dims, &(swKeys_ptr->mean_hf), mask);          computeBh(bx_err, by_err, bh_err, bx, by, bz, bh, dims, &(swKeys_ptr->mean_hf), mask, bitmask);
  
         if (computeGamma(bx, by, bz, bh, dims, &(swKeys_ptr->mean_gamma), mask))          if (computeGamma(bz_err, bh_err, bx, by, bz, bh, dims, &(swKeys_ptr->mean_gamma), &(swKeys_ptr->mean_gamma_err),mask, bitmask))
           {
                 swKeys_ptr->mean_gamma = DRMS_MISSING_FLOAT;                 swKeys_ptr->mean_gamma = DRMS_MISSING_FLOAT;
           swKeys_ptr->mean_gamma_err =  DRMS_MISSING_FLOAT;
       }
  
         computeB_total(bx, by, bz, bt, dims, mask);          computeB_total(bx_err, by_err, bz_err, bt_err, bx, by, bz, bt, dims, mask, bitmask);
  
         if (computeBtotalderivative(bt, dims, &(swKeys_ptr->mean_derivative_btotal), mask, derx_bt, dery_bt))          if (computeBtotalderivative(bt, dims, &(swKeys_ptr->mean_derivative_btotal), mask, bitmask, derx_bt, dery_bt, bt_err, &(swKeys_ptr->mean_derivative_btotal_err)))
       {
                 swKeys_ptr->mean_derivative_btotal = DRMS_MISSING_FLOAT;                 swKeys_ptr->mean_derivative_btotal = DRMS_MISSING_FLOAT;
                   swKeys_ptr->mean_derivative_btotal_err = DRMS_MISSING_FLOAT;
       }
  
         if (computeBhderivative(bh, dims, &(swKeys_ptr->mean_derivative_bh), mask, derx_bh, dery_bh))          if (computeBhderivative(bh, bh_err, dims, &(swKeys_ptr->mean_derivative_bh), &(swKeys_ptr->mean_derivative_bh_err), mask, bitmask, derx_bh, dery_bh))
       {
                 swKeys_ptr->mean_derivative_bh = DRMS_MISSING_FLOAT;                 swKeys_ptr->mean_derivative_bh = DRMS_MISSING_FLOAT;
           swKeys_ptr->mean_derivative_bh_err = DRMS_MISSING_FLOAT;
           }
  
         if (computeBzderivative(bz, dims, &(swKeys_ptr->mean_derivative_bz), mask, derx_bz, dery_bz))          if (computeBzderivative(bz, bz_err, dims, &(swKeys_ptr->mean_derivative_bz), &(swKeys_ptr->mean_derivative_bz_err), mask, bitmask, derx_bz, dery_bz))
       {
                 swKeys_ptr->mean_derivative_bz = DRMS_MISSING_FLOAT; // If fail, fill in NaN                 swKeys_ptr->mean_derivative_bz = DRMS_MISSING_FLOAT; // If fail, fill in NaN
           swKeys_ptr->mean_derivative_bz_err = DRMS_MISSING_FLOAT;
       }
  
           computeJz(bx_err, by_err, bx, by, dims, jz, jz_err, jz_err_squared, mask, bitmask, cdelt1, rsun_ref, rsun_obs,
                 derx, dery);
  
  
         if(computeJz(bx, by, dims, jz, &(swKeys_ptr->mean_jz), &(swKeys_ptr->us_i), mask,      if(computeJzsmooth(bx, by, dims, jz, jz_smooth, jz_err, jz_rms_err, jz_err_squared_smooth, &(swKeys_ptr->mean_jz),
                      cdelt1, rsun_ref, rsun_obs, derx, dery)) {                         &(swKeys_ptr->mean_jz_err), &(swKeys_ptr->us_i), &(swKeys_ptr->us_i_err), mask, bitmask, cdelt1,
                          rsun_ref, rsun_obs, derx, dery))
       {
                 swKeys_ptr->mean_jz = DRMS_MISSING_FLOAT;                 swKeys_ptr->mean_jz = DRMS_MISSING_FLOAT;
                 swKeys_ptr->us_i = DRMS_MISSING_FLOAT;                 swKeys_ptr->us_i = DRMS_MISSING_FLOAT;
           swKeys_ptr->mean_jz_err        = DRMS_MISSING_FLOAT;
           swKeys_ptr->us_i_err           = DRMS_MISSING_FLOAT;
         }         }
  
                         printf("swKeys_ptr->mean_jz=%f\n",swKeys_ptr->mean_jz);          if (computeAlpha(jz_err, bz_err, bz, dims, jz, jz_smooth, &(swKeys_ptr->mean_alpha), &(swKeys_ptr->mean_alpha_err), mask, bitmask, cdelt1, rsun_ref, rsun_obs))
       {
         if (computeAlpha(bz, dims, jz, &(swKeys_ptr->mean_alpha), mask, cdelt1, rsun_ref, rsun_obs))  
                 swKeys_ptr->mean_alpha = DRMS_MISSING_FLOAT;                 swKeys_ptr->mean_alpha = DRMS_MISSING_FLOAT;
           swKeys_ptr->mean_alpha_err     = DRMS_MISSING_FLOAT;
       }
  
         if (computeHelicity(bz, dims, jz, &(swKeys_ptr->mean_ih),          if (computeHelicity(jz_err, jz_rms_err, bz_err, bz, dims, jz, &(swKeys_ptr->mean_ih), &(swKeys_ptr->mean_ih_err), &(swKeys_ptr->total_us_ih), &(swKeys_ptr->total_abs_ih),
                                                 &(swKeys_ptr->total_us_ih), &(swKeys_ptr->total_abs_ih),                          &(swKeys_ptr->total_us_ih_err), &(swKeys_ptr->total_abs_ih_err), mask, bitmask, cdelt1, rsun_ref, rsun_obs))
                                                 mask, cdelt1, rsun_ref, rsun_obs)) {      {
                 swKeys_ptr->mean_ih = DRMS_MISSING_FLOAT;                 swKeys_ptr->mean_ih = DRMS_MISSING_FLOAT;
                 swKeys_ptr->total_us_ih = DRMS_MISSING_FLOAT;                 swKeys_ptr->total_us_ih = DRMS_MISSING_FLOAT;
                 swKeys_ptr->total_abs_ih = DRMS_MISSING_FLOAT;                 swKeys_ptr->total_abs_ih = DRMS_MISSING_FLOAT;
           swKeys_ptr->mean_ih_err        = DRMS_MISSING_FLOAT;
           swKeys_ptr->total_us_ih_err    = DRMS_MISSING_FLOAT;
           swKeys_ptr->total_abs_ih_err   = DRMS_MISSING_FLOAT;
         }         }
  
         if (computeSumAbsPerPolarity(bz, jz, dims, &(swKeys_ptr->totaljz),          if (computeSumAbsPerPolarity(jz_err, bz_err, bz, jz, dims, &(swKeys_ptr->totaljz), &(swKeys_ptr->totaljz_err),
                                                                  mask, cdelt1, rsun_ref, rsun_obs))                                                                   mask, bitmask, cdelt1, rsun_ref, rsun_obs))
       {
                 swKeys_ptr->totaljz = DRMS_MISSING_FLOAT;                 swKeys_ptr->totaljz = DRMS_MISSING_FLOAT;
           swKeys_ptr->totaljz_err        = DRMS_MISSING_FLOAT;
           }
  
           if (computeFreeEnergy(bx_err, by_err, bx, by, bpx, bpy, dims,
         if (computeFreeEnergy(bx, by, bpx, bpy, dims,                                                    &(swKeys_ptr->meanpot), &(swKeys_ptr->meanpot_err), &(swKeys_ptr->totpot), &(swKeys_ptr->totpot_err),
                                                   &(swKeys_ptr->meanpot), &(swKeys_ptr->totpot),                                                    mask, bitmask, cdelt1, rsun_ref, rsun_obs))
                                                   mask, cdelt1, rsun_ref, rsun_obs)) {      {
                 swKeys_ptr->meanpot = DRMS_MISSING_FLOAT; // If fail, fill in NaN                 swKeys_ptr->meanpot = DRMS_MISSING_FLOAT; // If fail, fill in NaN
                 swKeys_ptr->totpot = DRMS_MISSING_FLOAT;                 swKeys_ptr->totpot = DRMS_MISSING_FLOAT;
           swKeys_ptr->meanpot_err        = DRMS_MISSING_FLOAT;
           swKeys_ptr->totpot_err         = DRMS_MISSING_FLOAT;
         }         }
  
         if (computeShearAngle(bx, by, bz, bpx, bpy, bpz, dims,          if (computeShearAngle(bx_err, by_err, bh_err, bx, by, bz, bpx, bpy, bpz, dims,
                                                   &(swKeys_ptr->meanshear_angle), &(swKeys_ptr->area_w_shear_gt_45),                                                    &(swKeys_ptr->meanshear_angle), &(swKeys_ptr->meanshear_angle_err), &(swKeys_ptr->area_w_shear_gt_45),
                                                   &(swKeys_ptr->meanshear_angleh), &(swKeys_ptr->area_w_shear_gt_45h),                                                    mask, bitmask)) {
                                                   mask)) {  
                 swKeys_ptr->meanshear_angle = DRMS_MISSING_FLOAT; // If fail, fill in NaN                 swKeys_ptr->meanshear_angle = DRMS_MISSING_FLOAT; // If fail, fill in NaN
                 swKeys_ptr->area_w_shear_gt_45 = DRMS_MISSING_FLOAT;                 swKeys_ptr->area_w_shear_gt_45 = DRMS_MISSING_FLOAT;
                 swKeys_ptr->meanshear_angleh = DRMS_MISSING_FLOAT; // If fail, fill in NaN          swKeys_ptr->meanshear_angle_err= DRMS_MISSING_FLOAT;
                 swKeys_ptr->area_w_shear_gt_45h = DRMS_MISSING_FLOAT;  
         }         }
  
         // Clean up          // Clean up the arrays
  
           drms_free_array(bitmaskArray);          // Dec 18 2012 Xudong
         drms_free_array(maskArray);         drms_free_array(maskArray);
         drms_free_array(bxArray);         drms_free_array(bxArray);
         drms_free_array(byArray);         drms_free_array(byArray);
         drms_free_array(bzArray);         drms_free_array(bzArray);
  
         free(bh); free(bt); free(jz);          free(bh); free(bt); free(jz); free(jz_smooth);
         free(bpx); free(bpy); free(bpz);         free(bpx); free(bpy); free(bpz);
         free(derx); free(dery);         free(derx); free(dery);
         free(derx_bt); free(dery_bt);         free(derx_bt); free(dery_bt);
         free(derx_bz); free(dery_bz);         free(derx_bz); free(dery_bz);
         free(derx_bh); free(dery_bh);         free(derx_bh); free(dery_bh);
           free(bt_err); free(bh_err);  free(jz_err);
       free(jz_err_squared); free(jz_rms_err);
       free(jz_err_squared_smooth);
 } }
  
   
 /* /*
  * Set space weather indices, no error checking for now  * Set space weather indices, no error checking for now
  *  *
Line 1779  void setSWIndex(DRMS_Record_t *outRec, s
Line 1897  void setSWIndex(DRMS_Record_t *outRec, s
         drms_setkey_float(outRec, "TOTPOT", swKeys_ptr->totpot);         drms_setkey_float(outRec, "TOTPOT", swKeys_ptr->totpot);
         drms_setkey_float(outRec, "MEANSHR", swKeys_ptr->meanshear_angle);         drms_setkey_float(outRec, "MEANSHR", swKeys_ptr->meanshear_angle);
         drms_setkey_float(outRec, "SHRGT45", swKeys_ptr->area_w_shear_gt_45);         drms_setkey_float(outRec, "SHRGT45", swKeys_ptr->area_w_shear_gt_45);
       drms_setkey_float(outRec, "CMASK",   swKeys_ptr->count_mask);
       drms_setkey_float(outRec, "ERRBT",   swKeys_ptr->mean_derivative_btotal_err);
       drms_setkey_float(outRec, "ERRVF",   swKeys_ptr->mean_vf_err);
       drms_setkey_float(outRec, "ERRGAM",  swKeys_ptr->mean_gamma_err);
       drms_setkey_float(outRec, "ERRBH",   swKeys_ptr->mean_derivative_bh_err);
       drms_setkey_float(outRec, "ERRBZ",   swKeys_ptr->mean_derivative_bz_err);
       drms_setkey_float(outRec, "ERRJZ",   swKeys_ptr->mean_jz_err);
       drms_setkey_float(outRec, "ERRUSI",  swKeys_ptr->us_i_err);
       drms_setkey_float(outRec, "ERRALP",  swKeys_ptr->mean_alpha_err);
       drms_setkey_float(outRec, "ERRMIH",  swKeys_ptr->mean_ih_err);
       drms_setkey_float(outRec, "ERRTUI",  swKeys_ptr->total_us_ih_err);
       drms_setkey_float(outRec, "ERRTAI",  swKeys_ptr->total_abs_ih_err);
       drms_setkey_float(outRec, "ERRJHT",  swKeys_ptr->totaljz_err);
       drms_setkey_float(outRec, "ERRMPOT", swKeys_ptr->meanpot_err);
       drms_setkey_float(outRec, "ERRTPOT", swKeys_ptr->totpot_err);
       drms_setkey_float(outRec, "ERRMSHA", swKeys_ptr->meanshear_angle_err);
 }; };
  
   
 /* /*
  * Set all keywords, no error checking for now  * Set all keywords, no error checking for now
  *  *
  */  */
  
 void setKeys(DRMS_Record_t *outRec, DRMS_Record_t *inRec)  void setKeys(DRMS_Record_t *outRec, DRMS_Record_t *inRec, struct mapInfo *mInfo)
 { {
    copy_me_keys(inRec, outRec);    copy_me_keys(inRec, outRec);
    copy_patch_keys(inRec, outRec);    copy_patch_keys(inRec, outRec);
    copy_geo_keys(inRec, outRec);    copy_geo_keys(inRec, outRec);
    copy_ambig_keys(inRec, outRec);    copy_ambig_keys(inRec, outRec);
   
       int status = 0;
   
           // Change a few geometry keywords for CEA records
           if (mInfo != NULL) {
   
           drms_setkey_float(outRec, "CRPIX1", mInfo->ncol/2. + 0.5);
                   drms_setkey_float(outRec, "CRPIX2", mInfo->nrow/2. + 0.5);
   
                   drms_setkey_float(outRec, "CRVAL1", mInfo->xc);
                   drms_setkey_float(outRec, "CRVAL2", mInfo->yc);
                   drms_setkey_float(outRec, "CDELT1", mInfo->xscale);
                   drms_setkey_float(outRec, "CDELT2", mInfo->yscale);
                   drms_setkey_string(outRec, "CUNIT1", "degree");
                   drms_setkey_string(outRec, "CUNIT2", "degree");
   
                   char key[64];
                   snprintf (key, 64, "CRLN-%s", wcsCode[(int) mInfo->proj]);
                   drms_setkey_string(outRec, "CTYPE1", key);
                   snprintf (key, 64, "CRLT-%s", wcsCode[(int) mInfo->proj]);
                   drms_setkey_string(outRec, "CTYPE2", key);
                   drms_setkey_float(outRec, "CROTA2", 0.0);
   
           } else {
   
           float disk_xc = drms_getkey_float(inRec, "IMCRPIX1", &status);
           float disk_yc = drms_getkey_float(inRec, "IMCRPIX2", &status);
           float x_ll = drms_getkey_float(inRec, "CRPIX1", &status);
           float y_ll = drms_getkey_float(inRec, "CRPIX2", &status);
           // Defined as disk center's pixel address wrt lower-left of cutout
           drms_setkey_float(outRec, "CRPIX1", disk_xc - x_ll + 1.);
                   drms_setkey_float(outRec, "CRPIX2", disk_yc - y_ll + 1.);
                   // Always 0.
                   drms_setkey_float(outRec, "CRVAL1", 0);
                   drms_setkey_float(outRec, "CRVAL2", 0);
   
           }
   
           char timebuf[1024];
           float UNIX_epoch = -220924792.000; /* 1970.01.01_00:00:00_UTC */
           double val;
   
           val = drms_getkey_double(inRec, "DATE",&status);
           drms_setkey_double(outRec, "DATE_B", val);
           sprint_time(timebuf, (double)time(NULL) + UNIX_epoch, "ISO", 0);
           drms_setkey_string(outRec, "DATE", timebuf);
   
           // set cvs commit version into keyword HEADER
           char *cvsinfo = strdup("$Id$");
           char *cvsinfo2 = sw_functions_version();
           char cvsinfoall[2048];
           strcat(cvsinfoall,cvsinfo);
           strcat(cvsinfoall,"\n");
           strcat(cvsinfoall,cvsinfo2);
           status = drms_setkey_string(outRec, "CODEVER7", cvsinfoall);
   
 }; };
  
 // //
Line 1832  void frebin (float *image_in, float *ima
Line 2021  void frebin (float *image_in, float *ima
         fresize(&fresizes, image_in, image_out, nx, ny, nlead, nxout, nyout, nxout, xoff, yoff, DRMS_MISSING_FLOAT);         fresize(&fresizes, image_in, image_out, nx, ny, nlead, nxout, nyout, nxout, xoff, yoff, DRMS_MISSING_FLOAT);
  
 } }
   


Legend:
Removed from v.1.2  
changed lines
  Added in v.1.16

Karen Tian
Powered by
ViewCVS 0.9.4