version 1.20, 2014/02/13 04:41:06
|
version 1.21, 2014/02/18 23:35:31
|
|
|
/* | /* |
* sharp.c | * sharp.c |
* | * |
* This module creates the pipeline space weather harps |
* This module creates the pipeline Space Weather Active Region Patches (SHARPs). |
* It is a hard-coded strip-down version of bmap.c |
* It is a hard-coded strip-down version of bmap.c. |
* It takes the Mharp and Bharp series and crete the following quantities |
* It takes the Mharp and Bharp series and create the following quantities: |
|
* |
* Series 1: Sharp_CEA | * Series 1: Sharp_CEA |
* CEA remapped magnetogram, bitmap, continuum, doppler (same size in map coordinate, need manual spec?) | * CEA remapped magnetogram, bitmap, continuum, doppler (same size in map coordinate, need manual spec?) |
* CEA remapped vector field (Br, Bt, Bp) (same as above) | * CEA remapped vector field (Br, Bt, Bp) (same as above) |
* Space weather indices based on vector cutouts (step 2) | * Space weather indices based on vector cutouts (step 2) |
|
* |
* Series 2: Sharp_cutout: | * Series 2: Sharp_cutout: |
* cutouts of magnetogram, bitmap, continuum, doppler (HARP defined, various sizes in CCD pixels) | * cutouts of magnetogram, bitmap, continuum, doppler (HARP defined, various sizes in CCD pixels) |
* cutouts of all vector data segments (same as above) | * cutouts of all vector data segments (same as above) |
* Series 3: Other remaps |
|
* | * |
* Author: | * Author: |
* Xudong Sun; Monica Bobra | * Xudong Sun; Monica Bobra |
|
|
* Fixed disambig cutout for patch: 0 for even, 7 for odd | * Fixed disambig cutout for patch: 0 for even, 7 for odd |
* | * |
* | * |
* Example: |
* Example Calls: |
B (full disk disambiguation) |
* [I] B (full disk disambiguation) |
sharp "mharp=hmi.Mharp_720s[1832][2012.07.12_15:24]" "b=hmi_test.B_720s[2012.07.12_15:24]" "dop=hmi.V_720s[2012.07.12_15:24]" "cont=hmi.Ic_720s[2012.07.12_15:24]" "sharp_cea=su_xudong.sharp_cea_720s" "sharp_cut=su_xudong.sharp_720s" |
* > sharp "mharp=hmi.Mharp_720s[1832][2012.07.12_15:24]" "b=hmi_test.B_720s[2012.07.12_15:24]" "dop=hmi.V_720s[2012.07.12_15:24]" "cont=hmi.Ic_720s[2012.07.12_15:24]" "sharp_cea=su_xudong.sharp_cea_720s" "sharp_cut=su_xudong.sharp_720s" |
BHARP (patch disambiguation) |
* [II] BHARP (patch disambiguation) |
sharp "mharp=hmi.Mharp_720s[1832][2012.07.12_15:24]" "bharp=hmi.Bharp_720s[1832][2012.07.12_15:24]" "dop=hmi.V_720s[2012.07.12_15:24]" "cont=hmi.Ic_720s[2012.07.12_15:24]" "sharp_cea=su_xudong.sharp_cea_720s" "sharp_cut=su_xudong.sharp_720s" |
* > sharp "mharp=hmi.Mharp_720s[1832][2012.07.12_15:24]" "bharp=hmi.Bharp_720s[1832][2012.07.12_15:24]" "dop=hmi.V_720s[2012.07.12_15:24]" "cont=hmi.Ic_720s[2012.07.12_15:24]" "sharp_cea=su_xudong.sharp_cea_720s" "sharp_cut=su_xudong.sharp_720s" |
* |
|
* For comparison: |
|
* bmap "in=hmi_test.Bharp_720s_fd10[1404][2012.02.20_10:00]" \ |
|
"out=hmi_test.B_720s_CEA" -s -a "map=cyleqa" |
|
* | * |
* | * |
*/ | */ |
Line 178 struct swIndex { |
|
Line 175 struct swIndex { |
|
float meanpot_err; | float meanpot_err; |
float totpot_err; | float totpot_err; |
float meanshear_angle_err; | float meanshear_angle_err; |
|
float Rparam; |
}; | }; |
| |
// Mapping method | // Mapping method |
Line 1886 void computeSWIndex(struct swIndex *swKe |
|
Line 1884 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 |
| |
|
//Use magnetogram map to compute R |
|
DRMS_Segment_t *losSeg = drms_segment_lookup(inRec, "magnetogram"); |
|
DRMS_Array_t *losArray = drms_segment_read(losSeg, DRMS_TYPE_FLOAT, &status); |
|
float *los = (float *) losArray->data; // los |
|
|
DRMS_Segment_t *bz_errSeg = drms_segment_lookup(inRec, BR_ERR_SEG_CEA); | 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); | DRMS_Array_t *bz_errArray = drms_segment_read(bz_errSeg, DRMS_TYPE_FLOAT, &status); |
float *bz_err = (float *) bz_errArray->data; // bz_err | float *bz_err = (float *) bz_errArray->data; // bz_err |
Line 1911 void computeSWIndex(struct swIndex *swKe |
|
Line 1914 void computeSWIndex(struct swIndex *swKe |
|
| |
// convert cdelt1_orig from degrees to arcsec | // convert cdelt1_orig from degrees to arcsec |
float cdelt1 = (atan((rsun_ref*cdelt1_orig*RADSINDEG)/(dsun_obs)))*(1/RADSINDEG)*(3600.); | float cdelt1 = (atan((rsun_ref*cdelt1_orig*RADSINDEG)/(dsun_obs)))*(1/RADSINDEG)*(3600.); |
|
int nx1 = nx*cdelt1/2; |
|
int ny1 = ny*cdelt1/2; |
| |
// Temp arrays | // Temp arrays |
float *bh = (float *) (malloc(nxny * sizeof(float))); | float *bh = (float *) (malloc(nxny * sizeof(float))); |
Line 1934 void computeSWIndex(struct swIndex *swKe |
|
Line 1939 void computeSWIndex(struct swIndex *swKe |
|
float *jz_err_squared = (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_err_squared_smooth = (float *) (malloc(nxny * sizeof(float))); |
float *jz_rms_err = (float *) (malloc(nxny * sizeof(float))); | float *jz_rms_err = (float *) (malloc(nxny * sizeof(float))); |
//spaceweather quantities computed |
|
| |
|
// malloc some arrays for the R parameter calculation |
|
float *rim = (float *)malloc(nx1*ny1*sizeof(float)); |
|
float *p1p0 = (float *)malloc(nx1*ny1*sizeof(float)); |
|
float *p1n0 = (float *)malloc(nx1*ny1*sizeof(float)); |
|
float *p1p = (float *)malloc(nx1*ny1*sizeof(float)); |
|
float *p1n = (float *)malloc(nx1*ny1*sizeof(float)); |
|
float *p1 = (float *)malloc(nx1*ny1*sizeof(float)); |
|
float *pmap = (float *)malloc(nx1*ny1*sizeof(float)); |
| |
|
//spaceweather quantities computed |
if (computeAbsFlux(bz_err, bz , dims, &(swKeys_ptr->absFlux), &(swKeys_ptr->mean_vf), &(swKeys_ptr->mean_vf_err), | if (computeAbsFlux(bz_err, bz , dims, &(swKeys_ptr->absFlux), &(swKeys_ptr->mean_vf), &(swKeys_ptr->mean_vf_err), |
&(swKeys_ptr->count_mask), mask, bitmask, cdelt1, rsun_ref, rsun_obs)) | &(swKeys_ptr->count_mask), mask, bitmask, cdelt1, rsun_ref, rsun_obs)) |
{ | { |
Line 2033 void computeSWIndex(struct swIndex *swKe |
|
Line 2046 void computeSWIndex(struct swIndex *swKe |
|
swKeys_ptr->meanshear_angle_err= DRMS_MISSING_FLOAT; | swKeys_ptr->meanshear_angle_err= DRMS_MISSING_FLOAT; |
} | } |
| |
|
if (computeR(bz_err, los , dims, &(swKeys_ptr->Rparam), cdelt1, rim, p1p0, p1n0, |
|
p1p, p1n, p1, pmap, nx1, ny1)) |
|
{ |
|
swKeys_ptr->Rparam = DRMS_MISSING_FLOAT; // If fail, fill in NaN |
|
} |
|
|
|
|
// Clean up the arrays | // Clean up the arrays |
| |
drms_free_array(bitmaskArray); // Dec 18 2012 Xudong | drms_free_array(bitmaskArray); // Dec 18 2012 Xudong |
Line 2050 void computeSWIndex(struct swIndex *swKe |
|
Line 2070 void computeSWIndex(struct swIndex *swKe |
|
free(bt_err); free(bh_err); free(jz_err); | free(bt_err); free(bh_err); free(jz_err); |
free(jz_err_squared); free(jz_rms_err); | free(jz_err_squared); free(jz_rms_err); |
free(jz_err_squared_smooth); | free(jz_err_squared_smooth); |
|
|
|
free(rim); |
|
free(p1p0); |
|
free(p1n0); |
|
free(p1p); |
|
free(p1n); |
|
free(p1); |
|
free(pmap); |
|
|
} | } |
| |
/* | /* |
Line 2091 void setSWIndex(DRMS_Record_t *outRec, s |
|
Line 2120 void setSWIndex(DRMS_Record_t *outRec, s |
|
drms_setkey_float(outRec, "ERRMPOT", swKeys_ptr->meanpot_err); | drms_setkey_float(outRec, "ERRMPOT", swKeys_ptr->meanpot_err); |
drms_setkey_float(outRec, "ERRTPOT", swKeys_ptr->totpot_err); | drms_setkey_float(outRec, "ERRTPOT", swKeys_ptr->totpot_err); |
drms_setkey_float(outRec, "ERRMSHA", swKeys_ptr->meanshear_angle_err); | drms_setkey_float(outRec, "ERRMSHA", swKeys_ptr->meanshear_angle_err); |
|
drms_setkey_float(outRec, "R_VALUE", swKeys_ptr->Rparam); |
}; | }; |
| |
/* | /* |
Line 2176 void setKeys(DRMS_Record_t *outRec, DRMS |
|
Line 2206 void setKeys(DRMS_Record_t *outRec, DRMS |
|
| |
} | } |
| |
char timebuf[1024]; |
TIME val, trec, tnow, UNIX_epoch = -220924792.000; /* 1970.01.01_00:00:00_UTC */ |
float UNIX_epoch = -220924792.000; /* 1970.01.01_00:00:00_UTC */ |
tnow = (double)time(NULL); |
double val; |
tnow += UNIX_epoch; |
|
|
val = drms_getkey_double(bharpRec, "DATE", &status); |
val = drms_getkey_time(bharpRec, "DATE", &status); |
drms_setkey_double(outRec, "DATE_B", val); |
drms_setkey_time(outRec, "DATE_B", val); |
sprint_time(timebuf, (double)time(NULL) + UNIX_epoch, "ISO", 0); |
drms_setkey_time(outRec, "DATE", tnow); |
drms_setkey_string(outRec, "DATE", timebuf); |
|
| |
// set cvs commit version into keyword HEADER | // set cvs commit version into keyword HEADER |
char *cvsinfo = strdup("$Id$"); | char *cvsinfo = strdup("$Id$"); |