JSOC/proj/sharp/apps/update_sharp_keys.c - diff

Return to update_sharp_keys.c CVS log

Up to [Development] / JSOC / proj / sharp / apps

Diff for /JSOC/proj/sharp/apps/update_sharp_keys.c between version 1.12 and 1.19

version 1.12, 2015/02/27 19:50:51

version 1.19, 2021/05/26 19:25:17

Line 1

* MODULE NAME: update_sharp_keys.c

* DESCRIPTION: This module recalculates SHARP keywords. This is accomplished by

* DESCRIPTION: This module recalculates SHARP keywords.

* cloning a record, recalculating keywords of choice (i.e., user input),

* This is accomplished by cloning a record, recalculating keywords of choice

* and pointing to the same segments as the old record. Associated error keys

* (i.e., user input), and pointing to the same segments as the old record.

* are computed by default.

* Associated error keys are computed by default.

* This module accounts for versioning by appending to the keyword CODEVER7 and HISTORY:

* CODEVER7 will contain multiple lines: the production build of sharp.c, the

Line 24

* -- DRMS SHARP CEA series

* -- HARPNUM

* -- comma separated list of keywords to recalculate

* -- DEBUG flag (use like this debug=debug)

* -- DEBUG flag

* AUTHOR : Monica Bobra

* Version : v0.0 Jun 14 2013

* v0.1 Feb 11 2014 Added different input and output series

* EXAMPLE :

* update_sharp_keys sharpseriesin=hmi.sharp_720s sharpceaseriesin=hmi.sharp_cea_720s //

* HARPNUM=1 sharpseriesout=hmi.sharp_720s sharpceaseriesout=hmi.sharp_cea_720s keylist=USFLUX,TOTPOT debug=debug

Line 66

Line 63

#define ARRLENGTH(ARR) (sizeof(ARR) / sizeof(ARR[0]))

#define DIE(msg) {fflush(stdout); fprintf(stderr,"%s, status=%d\n", msg, status); return(status);}

#define SHOW(msg) {printf("%s", msg); fflush(stdout);}

#define DEBUG(true)

char *module_name = "update_sharp_keys"; /* Module name */

char *version_id = "2014 Feb 12"; /* Version number */

char *version_id = "2021 May 25"; /* Version number */

ModuleArgs_t module_args[] =

{

Line 77 ModuleArgs_t module_args[] =

Line 75 ModuleArgs_t module_args[] =

{ARG_STRING, "sharpseriesout", NULL, "Output Sharp dataseries"},

{ARG_STRING, "sharpceaseriesout", NULL, "Output Sharp CEA dataseries"},

{ARG_INT, "HARPNUM", NULL, "HARP number"},

{ARG_STRING, "keylist", NULL, "comma separated list of keywords to update"},

{ARG_STRING, "keylist", NULL, "comma separated list of keywords to update (for LOS keywords enter MEANGBZLOS or USFLUXLOS)"},

{ARG_STRING, "debug", NULL, "debug mode functionality. use like this: debug=debug"},

{ARG_END}

};

Line 93 int DoIt(void)

Line 91 int DoIt(void)

/* Keywords */

float mean_vf;

float mean_vf_los;

float absFlux;

float absFlux_los;

float count_mask;

float count_mask_los;

float mean_hf;

float mean_gamma;

float mean_derivative_btotal;

float mean_derivative_bh;

float mean_derivative_bz;

float mean_derivative_los;

float mean_jz;

float us_i;

float mean_alpha;

Line 191 int DoIt(void)

Line 193 int DoIt(void)

DIE("Problem creating sharp cea records.");

}

char *keylist = (char *) params_get_str(&cmdparams, "keylist");

char *debug = (char *) params_get_str(&cmdparams, "debug");

// Flags to indicate which keyword will be recalculated

int meanvfflag = (strstr(keylist,"USFLUX") != NULL); // generalize so that lowercase is acceptable

int meanvflosflag = (strstr(keylist,"USFLUXL") != NULL);

int meanglosflag = (strstr(keylist,"MEANGBL") != NULL);

int meanvfflag = (strstr(keylist,"USFLUX") != NULL);

int totpotflag = (strstr(keylist,"TOTPOT") != NULL);

int meangamflag = (strstr(keylist,"MEANGAM") != NULL);

int meangbtflag = (strstr(keylist,"MEANGBT") != NULL);

Line 222 int DoIt(void)

Line 225 int DoIt(void)

int epszflag = (strstr(keylist,"EPSZ") != NULL);

int debugflag = (strstr(debug,"debug") != NULL);

for (irec=0;irec<nrecs;irec++)

{

DRMS_Record_t *sharpinrec = sharpinrecset->records[irec];

DRMS_Record_t *sharpceainrec = sharpceainrecset->records[irec];

DRMS_Segment_t *inseg = drms_segment_lookup(sharpceainrec, "Br");

Line 257 int DoIt(void)

float *jz_rms_err = (float *) (malloc(nxny * sizeof(float)));

float *jz_err_squared_smooth = (float *) (malloc(nxny * sizeof(float)));

float *jz_smooth = (float *) (malloc(nxny * sizeof(float)));

float *err_term1 = (float *) (malloc(nxny * sizeof(float)));

float *err_term1 = (float *) (calloc(nxny, sizeof(float)));

float *err_term2 = (float *) (malloc(nxny * sizeof(float)));

float *err_term2 = (float *) (calloc(nxny, sizeof(float)));

float *err_termA = (float *) (calloc(nxny, sizeof(float)));

float *fx = (float *) (malloc(nxny * sizeof(float)));

float *err_termB = (float *) (calloc(nxny, sizeof(float)));

float *fy = (float *) (malloc(nxny * sizeof(float)));

float *err_termAt = (float *) (calloc(nxny, sizeof(float)));

float *fz = (float *) (malloc(nxny * sizeof(float)));

float *err_termBt = (float *) (calloc(nxny, sizeof(float)));

float *derx_los = (float *) (malloc(nxny * sizeof(float)));

float *err_termAh = (float *) (calloc(nxny, sizeof(float)));

float *dery_los = (float *) (malloc(nxny * sizeof(float)));

float *err_termBh = (float *) (calloc(nxny, sizeof(float)));

// ephemeris variables

float cdelt1_orig, cdelt1, dsun_obs, imcrpix1, imcrpix2, crpix1, crpix2;

Line 290 int DoIt(void)

Line 289 int DoIt(void)

float UNIX_epoch = -220924792.000;

sprint_time(timebuf, (double)time(NULL) + UNIX_epoch, "ISO", 0);

sprintf(historyofthemodule,"The following keywords were re-computed on %s: %s.",timebuf,keylist);

printf("historyofthemodule=%s\n",historyofthemodule);

//printf("historyofthemodule=%s\n",historyofthemodule);

history0 = drms_getkey_string(sharpinrec, "HISTORY", &status);

strcat(historyofthemodule,"\n");

strcat(historyofthemodule,history0);

Line 301 int DoIt(void)

Line 300 int DoIt(void)

dsun_obs = drms_getkey_float(testrec, "DSUN_OBS", &status);

rsun_ref = drms_getkey_double(testrec, "RSUN_REF", &status);

cdelt1_orig = drms_getkey_float(testrec, "CDELT1", &status);

cdelt1 = (atan((rsun_ref*cdelt1_orig*RADSINDEG)/(dsun_obs)))*(1/RADSINDEG)*(3600.);

int scale = round(2.0/cdelt1);

int nx1 = nx/scale;

int ny1 = ny/scale;

Line 313 int DoIt(void)

Line 311 int DoIt(void)

DIE("X-dimension of output array in fsample() is too large.");

if (ny1 > floor((ny-1)/scale + 1) )

DIE("Y-dimension of output array in fsample() is too large.");

// malloc some arrays for the R parameter calculation

float *rim = (float *)malloc(nx1*ny1*sizeof(float));

float *p1p0 = (float *)malloc(nx1*ny1*sizeof(float));

Line 325 int DoIt(void)

Line 322 int DoIt(void)

float *p1pad = (float *)malloc(nxp*nyp*sizeof(float));

float *pmapn = (float *)malloc(nx1*ny1*sizeof(float));

// define some arrays for the lorentz force calculation

float *fx = (float *) (malloc(nxny * sizeof(float)));

float *fy = (float *) (malloc(nxny * sizeof(float)));

float *fz = (float *) (malloc(nxny * sizeof(float)));

// Get emphemeris

sharpinrec = sharpinrecset->records[irec];

sharpceainrec = sharpceainrecset->records[irec];

Line 341 int DoIt(void)

Line 333 int DoIt(void)

imcrpix2 = drms_getkey_float(sharpinrec, "IMCRPIX2", &status);

crpix1 = drms_getkey_float(sharpinrec, "CRPIX1", &status);

crpix2 = drms_getkey_float(sharpinrec, "CRPIX2", &status);

cdelt1 = (atan((rsun_ref*cdelt1_orig*RADSINDEG)/(dsun_obs)))*(1/RADSINDEG)*(3600.);

sharpoutrec = sharpoutrecset->records[irec];

sharpceaoutrec = sharpceaoutrecset->records[irec];

Line 367 int DoIt(void)

Line 360 int DoIt(void)

drms_setkey_time(sharpoutrec,"DATE",tnow);

drms_setkey_time(sharpceaoutrec,"DATE",tnow);

// Get bx, by, bz, mask

// Get all the segments if computing vector magnetic field keywords

// Use HARP (Turmon) bitmap as a threshold on spaceweather quantities

DRMS_Segment_t *bitmaskSeg = drms_segment_lookup(sharpceainrec, "bitmap");

DRMS_Array_t *bitmaskArray = drms_segment_read(bitmaskSeg, DRMS_TYPE_INT, &status);

int *bitmask = (int *) bitmaskArray->data; // get the previously made mask array

//Use conf_disambig map as a threshold on spaceweather quantities

DRMS_Segment_t *maskSeg = drms_segment_lookup(sharpceainrec, "conf_disambig");

DRMS_Array_t *maskArray = drms_segment_read(maskSeg, DRMS_TYPE_INT, &status);

if (status != DRMS_SUCCESS)

DIE("No CONF_DISAMBIG segment.");

int *mask = (int *) maskArray->data; // get the previously made mask array

//Use magnetogram map to compute R

DRMS_Segment_t *losSeg = drms_segment_lookup(sharpceainrec, "magnetogram");

DRMS_Array_t *losArray = drms_segment_read(losSeg, DRMS_TYPE_FLOAT, &status);

if (status != DRMS_SUCCESS)

DIE("No magnetogram segment.");

float *los = (float *) losArray->data; // los

//printf("Line 387");

DRMS_Segment_t *bxSeg = drms_segment_lookup(sharpceainrec, "Bp");

DRMS_Array_t *bxArray = drms_segment_read(bxSeg, DRMS_TYPE_FLOAT, &status);

if (status != DRMS_SUCCESS)

DIE("No Bp segment.");

float *bx = (float *) bxArray->data; // bx

//printf("Line 393");

DRMS_Segment_t *bySeg = drms_segment_lookup(sharpceainrec, "Bt");

DRMS_Array_t *byArray = drms_segment_read(bySeg, DRMS_TYPE_FLOAT, &status);

if (status != DRMS_SUCCESS)

DIE("No Bt segment.");

float *by = (float *) byArray->data; // by

for (int i = 0; i < nxny; i++) by[i] *= -1;

//printf("Line 399");

DRMS_Segment_t *bzSeg = drms_segment_lookup(sharpceainrec, "Br");

//printf("Line 408\n");

if (bzSeg == NULL)

DIE("somekind of pointer problem.");

DRMS_Array_t *bzArray = drms_segment_read(bzSeg, DRMS_TYPE_FLOAT, &status);

//printf("status=%d\n",status);

if (status != DRMS_SUCCESS)

DIE("No Br segment.");

float *bz = (float *) bzArray->data; // bz

//printf("Line 414");

DRMS_Segment_t *bz_errSeg = drms_segment_lookup(sharpceainrec, "Br_err");

DRMS_Array_t *bz_errArray = drms_segment_read(bz_errSeg, DRMS_TYPE_FLOAT, &status);

if (status != DRMS_SUCCESS)

DIE("No Br_err segment.");

float *bz_err = (float *) bz_errArray->data; // bz_err

//printf("Line 408");

DRMS_Segment_t *by_errSeg = drms_segment_lookup(sharpceainrec, "Bt_err");

DRMS_Array_t *by_errArray = drms_segment_read(by_errSeg, DRMS_TYPE_FLOAT, &status);

if (status != DRMS_SUCCESS)

DIE("No Bt_err segment.");

float *by_err = (float *) by_errArray->data; // by_err

//printf("Line 412");

DRMS_Segment_t *bx_errSeg = drms_segment_lookup(sharpceainrec, "Bp_err");

DRMS_Array_t *bx_errArray = drms_segment_read(bx_errSeg, DRMS_TYPE_FLOAT, &status);

if (status != DRMS_SUCCESS)

DIE("No Bp_err segment.");

float *bx_err = (float *) bx_errArray->data; // bx_err

/***** USFLUX, Example Function 1 *************************************/

Line 443 int DoIt(void)

Line 459 int DoIt(void)

drms_setkey_float(sharpoutrec, "R_VALUE", Rparam);

drms_setkey_float(sharpceaoutrec, "R_VALUE", Rparam);

}

/***** MEANGBZ, Example Function 7 ************************************/

if (meangbzflag)

{

// Compute Bz derivative

if (computeBzderivative(bz, bz_err, dims, &mean_derivative_bz, &mean_derivative_bz_err, mask, bitmask, derx_bz, dery_bz, err_term1, err_term2))

{

mean_derivative_bz = DRMS_MISSING_FLOAT;

mean_derivative_bz_err = DRMS_MISSING_FLOAT;

}

/***** MEANPOT and TOTPOT, Example Function 13 (Requires Keiji's code) **/

// Set sharp keys

drms_setkey_float(sharpoutrec, "MEANGBZ", mean_derivative_bz);

drms_setkey_float(sharpoutrec, "ERRBZ", mean_derivative_bz_err);

// Set sharp cea keys

drms_setkey_float(sharpceaoutrec, "MEANGBZ", mean_derivative_bz);

drms_setkey_float(sharpceaoutrec, "ERRBZ", mean_derivative_bz_err);

}

/***** MEANPOT and TOTPOT, Example Function 13 (Requires Keiji's code) **/

if (totpotflag || meanpotflag)

{

// First compute potential field

Line 479 int DoIt(void)

Line 512 int DoIt(void)

}

/***** MEANGAM, Example Function 3 ************************************/

if (meangamflag)

{

// First compute horizontal field

Line 502 int DoIt(void)

Line 534 int DoIt(void)

}

/***** MEANGBT, Example Function 5 (Requires Function 4) *************/

if (meangbtflag)

{

// First, compute Bt

Line 510 int DoIt(void)

Line 541 int DoIt(void)

// Then take the derivative of Bt

if (computeBtotalderivative(bt, dims, &mean_derivative_btotal, mask, bitmask, derx_bt, dery_bt, bt_err,

&mean_derivative_btotal_err, err_termAt, err_termBt))

&mean_derivative_btotal_err, err_term1, err_term2))

{

mean_derivative_btotal = DRMS_MISSING_FLOAT;

mean_derivative_btotal_err = DRMS_MISSING_FLOAT;

Line 526 int DoIt(void)

Line 557 int DoIt(void)

}

/***** MEANGBH, Example Function 6 (Requires Function 2) *************/

if (meangbhflag)

{

// First, compute Bh

computeBh(bx_err, by_err, bh_err, bx, by, bz, bh, dims, &mean_hf, mask, bitmask);

// Then take the derivative of Bh

if (computeBhderivative(bh, bh_err, dims, &mean_derivative_bh, &mean_derivative_bh_err, mask, bitmask, derx_bh, dery_bh, err_termAh, err_termBh))

if (computeBhderivative(bh, bh_err, dims, &mean_derivative_bh, &mean_derivative_bh_err, mask, bitmask, derx_bh, dery_bh, err_term1, err_term2))

{

mean_derivative_bh = DRMS_MISSING_FLOAT;

mean_derivative_bh_err = DRMS_MISSING_FLOAT;

Line 548 int DoIt(void)

Line 578 int DoIt(void)

drms_setkey_float(sharpceaoutrec, "ERRBH", mean_derivative_bh_err);

}

/***** MEANGBZ, Example Function 7 ************************************/

if (meangbzflag)

{

// Compute Bz derivative

if (computeBzderivative(bz, bz_err, dims, &mean_derivative_bz, &mean_derivative_bz_err, mask, bitmask, derx_bz, dery_bz, err_termA, err_termB))

{

mean_derivative_bz = DRMS_MISSING_FLOAT;

mean_derivative_bz_err = DRMS_MISSING_FLOAT;

}

// Set sharp keys

drms_setkey_float(sharpoutrec, "MEANGBZ", mean_derivative_bz);

drms_setkey_float(sharpoutrec, "ERRBZ", mean_derivative_bz_err);

// Set sharp cea keys

drms_setkey_float(sharpceaoutrec, "MEANGBZ", mean_derivative_bz);

drms_setkey_float(sharpceaoutrec, "ERRBZ", mean_derivative_bz_err);

}

/***** MEANJZD and TOTUSJZ, Example Function 9 (Requires Function 8) ***/

if (meanjzdflag || totusjzflag)

{

// First, compute Jz on all pixels

Line 581 int DoIt(void)

Line 590 int DoIt(void)

&mean_jz_err, &us_i, &us_i_err, mask, bitmask, cdelt1,

rsun_ref, rsun_obs, derx, dery))

{

mean_jz = DRMS_MISSING_FLOAT;

us_i = DRMS_MISSING_FLOAT;

Line 603 int DoIt(void)

Line 611 int DoIt(void)

}

/***** MEANALP, Example Function 10 (Requires Function 8)*********/

if (meanalpflag)

{

// First, compute Jz on all pixels

Line 629 int DoIt(void)

Line 636 int DoIt(void)

}

/***** MEANJZH, TOTUSJH, ABSNJZH, Example Function 11 (Requires Function 8) ***/

if (meanjzhflag || totusjhflag || absnjzhflag)

{

// First, compute Jz on all pixels

Line 666 int DoIt(void)

Line 672 int DoIt(void)

}

/***** SAVNCPP, Example Function 12 (Requires Function 8) *******************/

if (savncppflag)

{

// First, compute Jz on all pixels

Line 691 int DoIt(void)

Line 696 int DoIt(void)

}

/***** MEANSHR and SHRGT45, Example Function 14 (Requires Keiji's code) **********/

if (meanshrflag || shrgt45flag)

{

// First compute potential field

Line 720 int DoIt(void)

Line 724 int DoIt(void)

}

/***** TOTFX, TOTFY, TOTFX, TOTBSQ, EPSX, EPSY, EPSZ , Example Function 16 (Lorentz forces) **********/

if (totfxflag || totfyflag || totfzflag || totbsqflag || epsxflag || epsyflag || epszflag)

{

// Compute Lorentz forces

if (computeLorentz(bx, by, bz, fx, fy, fz, dims, &totfx, &totfy, &totfz, &totbsq,

&epsx, &epsy, &epsz, mask, bitmask, cdelt1, rsun_ref, rsun_obs))

Line 757 int DoIt(void)

Line 759 int DoIt(void)

}

/******************************* END FLAGS **********************************/

/***** USFLUX, Example Function 17 *************************************/

if (meanvflosflag)

{

//printf("772");

if (computeAbsFlux_los(los, dims, &absFlux_los, &mean_vf_los,

&count_mask_los, bitmask, cdelt1, rsun_ref, rsun_obs))

{

absFlux_los = DRMS_MISSING_FLOAT; // If fail, fill in NaN

mean_vf_los = DRMS_MISSING_FLOAT;

count_mask_los = DRMS_MISSING_INT;

}

//printf("line757\n");

//printf("mean_vf_los=%f\n",mean_vf_los);

//printf("cdelt1=%f\n",cdelt1);

//printf("rsun_ref=%f\n",rsun_ref);

//printf("rsun_obs=%f\n",rsun_obs);

//printf("Key #%d done.\n", irec);

drms_setkey_float(sharpoutrec, "USFLUXL", mean_vf_los);

drms_setkey_float(sharpoutrec, "CMASKL", count_mask_los);

drms_setkey_float(sharpceaoutrec, "USFLUXL", mean_vf_los);

drms_setkey_float(sharpceaoutrec, "CMASKL", count_mask_los);

}

/***** MEANGBZ, Example Function 18 ************************************/

if (meanglosflag)

{

//printf("line 794\n");

// Compute Bz derivative

if (computeLOSderivative(los, dims, &mean_derivative_los, bitmask, derx_los, dery_los))

{

//printf("line 799\n");

mean_derivative_los = DRMS_MISSING_FLOAT;

}

// Set sharp keys

//printf("mean_derivative_los=%f\n",mean_derivative_los);

drms_setkey_float(sharpoutrec, "MEANGBL", mean_derivative_los);

// Set sharp cea keys

drms_setkey_float(sharpceaoutrec, "MEANGBL", mean_derivative_los);

}

/******************************* END FLAGS **********************************/

//printf("line809\n");

drms_free_array(bitmaskArray);

drms_free_array(maskArray);

drms_free_array(bxArray);

Line 769 int DoIt(void)

Line 813 int DoIt(void)

drms_free_array(by_errArray);

drms_free_array(bz_errArray);

drms_free_array(losArray);

// free the arrays that are related to the lorentz calculation

free(fx); free(fy); free(fz);

//used for select calculations

free(bh); free(bt); free(jz);

free(bpx); free(bpy); free(bpz);

free(derx); free(dery);

free(derx_bt); free(dery_bt);

free(derx_bz); free(dery_bz);

free(derx_bh); free(dery_bh);

free(bt_err); free(bh_err); free(jz_err);

free(jz_err_squared); free(jz_rms_err);

free(cvsinfoall);

free(rim);

free(p1p0);

free(p1n0);

Line 792 int DoIt(void)

Line 823 int DoIt(void)

free(pmap);

free(p1pad);

free(pmapn);

free(fx); free(fy); free(fz);

// free the arrays that are related to the numerical derivatives

free(bh); free(bt); free(jz);

free(bpx); free(bpy); free(bpz);

free(derx); free(dery);

free(derx_bt); free(dery_bt);

free(derx_los); free(dery_los);

free(derx_bz); free(dery_bz);

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);

free(jz_smooth);

free(err_term2);

free(err_term1);

free(err_termB);

free(err_termA);

free(err_termBt);

free(err_termAt);

free(err_termBh);

free(err_termAh);

} //endfor

// Close all the records

drms_close_records(sharpinrecset, DRMS_FREE_RECORD);

drms_close_records(sharpceainrecset, DRMS_FREE_RECORD);

Line 812 int DoIt(void)

Line 848 int DoIt(void)

drms_close_records(sharpceaoutrecset, DRMS_INSERT_RECORD);

printf("Success=%d\n",sameflag);

return 0;

} // DoIt

Legend:

Removed from v.1.12
changed lines
	Added in v.1.19

Karen Tian