Return to sw_functions.c CVS log Up to [Development] / JSOC / proj / sharp / apps

### Diff for /JSOC/proj/sharp/apps/sw_functions.c between version 1.31 and 1.32

version 1.31, 2014/06/05 21:27:04 version 1.32, 2014/09/05 21:59:48
 Line 832  int computeHelicity(float *jz_err, float
 Line 832  int computeHelicity(float *jz_err, float
/* Example function 12:  Sum of Absolute Value per polarity  */ /* Example function 12:  Sum of Absolute Value per polarity  */

//  The Sum of the Absolute Value per polarity is defined as the following: //  The Sum of the Absolute Value per polarity is defined as the following:
//  fabs(sum(jz gt 0)) + fabs(sum(jz lt 0)) and the units are in Amperes.  //  fabs(sum(jz gt 0)) + fabs(sum(jz lt 0)) and the units are in Amperes per arcsecond.
//  The units of jz are in G/pix. In this case, we would have the following: //  The units of jz are in G/pix. In this case, we would have the following:
//  Jz = (Gauss/pix)(1/CDELT1)(0.00010)(1/MUNAUGHT)(RSUN_REF/RSUN_OBS)(RSUN_REF/RSUN_OBS)(RSUN_OBS/RSUN_REF), //  Jz = (Gauss/pix)(1/CDELT1)(0.00010)(1/MUNAUGHT)(RSUN_REF/RSUN_OBS)(RSUN_REF/RSUN_OBS)(RSUN_OBS/RSUN_REF),
//     = (Gauss/pix)(1/CDELT1)(0.00010)(1/MUNAUGHT)(RSUN_REF/RSUN_OBS) //     = (Gauss/pix)(1/CDELT1)(0.00010)(1/MUNAUGHT)(RSUN_REF/RSUN_OBS)
 Line 869  int computeSumAbsPerPolarity(float *jz_e
 Line 869  int computeSumAbsPerPolarity(float *jz_e
}         }
}     }

*totaljzptr    = fabs(sum1) + fabs(sum2);  /* Units are A */      *totaljzptr    = fabs(sum1) + fabs(sum2);  /* Units are Amperes per arcsecond */
*totaljz_err_ptr = sqrt(err)*(1/cdelt1)*fabs((0.00010)*(1/MUNAUGHT)*(rsun_ref/rsun_obs));     *totaljz_err_ptr = sqrt(err)*(1/cdelt1)*fabs((0.00010)*(1/MUNAUGHT)*(rsun_ref/rsun_obs));
//printf("SAVNCPP=%g\n",*totaljzptr);     //printf("SAVNCPP=%g\n",*totaljzptr);
//printf("SAVNCPP_err=%g\n",*totaljz_err_ptr);     //printf("SAVNCPP_err=%g\n",*totaljz_err_ptr);
 Line 1160  int computeR(float *bz_err, float *los,
 Line 1160  int computeR(float *bz_err, float *los,
for (j = 0; j < ny1; j++)         for (j = 0; j < ny1; j++)
{         {
index = j * nx1 + i;             index = j * nx1 + i;
if isnan(pmapn[index]) continue;
if isnan(rim[index]) continue;
sum += pmapn[index]*abs(rim[index]);             sum += pmapn[index]*abs(rim[index]);
}         }
}     }
 Line 1169  int computeR(float *bz_err, float *los,
 Line 1171  int computeR(float *bz_err, float *los,
else     else
*Rparam = log10(sum);         *Rparam = log10(sum);

//printf("R_VALUE=%f\n",*Rparam);

free_fresize(&fresboxcar);     free_fresize(&fresboxcar);
free_fresize(&fresgauss);     free_fresize(&fresgauss);

 Line 1201  int computeLorentz(float *bx, float *by
 Line 1205  int computeLorentz(float *bx, float *by
double k_h = -1.0 * area / (4. * PI) / 1.0e20;     double k_h = -1.0 * area / (4. * PI) / 1.0e20;
double k_z = area / (8. * PI) / 1.0e20;     double k_z = area / (8. * PI) / 1.0e20;

/* Multiplier */
float vectorMulti[] = {1.,-1.,1.};

if (nx <= 0 || ny <= 0) return 1;     if (nx <= 0 || ny <= 0) return 1;

for (int i = 0; i < nxny; i++)     for (int i = 0; i < nxny; i++)
 Line 1212  int computeLorentz(float *bx, float *by
 Line 1213  int computeLorentz(float *bx, float *by
if isnan(bx[i]) continue;        if isnan(bx[i]) continue;
if isnan(by[i]) continue;        if isnan(by[i]) continue;
if isnan(bz[i]) continue;        if isnan(bz[i]) continue;
fx[i]  = (bx[i] * vectorMulti[0]) * (bz[i] * vectorMulti[2]) * k_h;         fx[i]  = bx[i] * bz[i] * k_h;
fy[i]  = (by[i] * vectorMulti[1]) * (bz[i] * vectorMulti[2]) * k_h;         fy[i]  = by[i] * bz[i] * k_h;
fz[i]  = (bx[i] * bx[i] + by[i] * by[i] - bz[i] * bz[i]) * k_z;        fz[i]  = (bx[i] * bx[i] + by[i] * by[i] - bz[i] * bz[i]) * k_z;
bsq    = bx[i] * bx[i] + by[i] * by[i] + bz[i] * bz[i];        bsq    = bx[i] * bx[i] + by[i] * by[i] + bz[i] * bz[i];
totfx  += fx[i]; totfy += fy[i]; totfz += fz[i];        totfx  += fx[i]; totfy += fy[i]; totfz += fz[i];
 Line 1228  int computeLorentz(float *bx, float *by
 Line 1229  int computeLorentz(float *bx, float *by
*epsy_ptr   = (totfy / k_h) / totbsq;     *epsy_ptr   = (totfy / k_h) / totbsq;
*epsz_ptr   = (totfz / k_z) / totbsq;     *epsz_ptr   = (totfz / k_z) / totbsq;

//printf("TOTBSQ=%f\n",*totbsq_ptr);

return 0;     return 0;

} }

Legend:
 Removed from v.1.31 changed lines Added in v.1.32