IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Jan 17, 2015, 6:11:34 AM (12 years ago)
Author:
eugene
Message:

allow dvomath to handle 3D (and higher) buffers

File:
1 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20150112/Ohana/src/opihi/lib.shell/stack_math.c

    r36528 r37857  
    127127int MMM_trinary (StackVar *OUT, StackVar *V1, StackVar *V2, StackVar *V3, char *op) {
    128128
    129   int i, Nx, Ny;
     129  int i, Npix;
    130130  float *out, *M1, *M2, *M3;
    131131  char line[512]; // this is only used to report an error
    132132 
    133   Nx = V1[0].buffer[0].matrix.Naxis[0];
    134   Ny = V1[0].buffer[0].matrix.Naxis[1];
    135 
     133  Npix = 1;
     134  for (i = 0; i < V1[0].buffer[0].matrix.Naxes; i++) {
     135    if (V1[0].buffer[0].matrix.Naxis[i] == 0) break;
     136    Npix *= V1[0].buffer[0].matrix.Naxis[i];
     137  }
     138 
    136139  if (V1[0].type == ST_MATRIX_TMP) {  /** use V1 as temp buffer **/
    137140    OUT[0].buffer = V1[0].buffer;
    … …  
    154157
    155158# define MMM_FUNC(OP)                                   \
    156   for (i = 0; i < Nx*Ny; i++, out++, M1++, M2++, M3++) {        \
     159  for (i = 0; i < Npix; i++, out++, M1++, M2++, M3++) { \
    157160    *out = OP;                                          \
    158161  }                                                     \
    … …  
    531534}
    532535
    533 // the vector is applied to each column
     536// the vector is applied to each column (currently only valid for 2D matrix)
    534537int MV_binary (StackVar *OUT, StackVar *V1, StackVar *V2, char *op) {
    535538
    … …  
    626629}
    627630
    628 // the vector is applied to each row
     631// the vector is applied to each row (currently only valid for 2D matrix)
    629632int VM_binary (StackVar *OUT, StackVar *V1, StackVar *V2, char *op) {
    630633
    … …  
    722725int MM_binary (StackVar *OUT, StackVar *V1, StackVar *V2, char *op) {
    723726
    724   int i, Nx, Ny;
     727  int i, Npix;
    725728  float *out, *M1, *M2;
    726729  char line[512]; // this is only used to report an error
    727730 
    728   Nx = V1[0].buffer[0].matrix.Naxis[0];
    729   Ny = V1[0].buffer[0].matrix.Naxis[1];
    730 
     731  Npix = 1;
     732  for (i = 0; i < V1[0].buffer[0].matrix.Naxes; i++) {
     733    if (V1[0].buffer[0].matrix.Naxis[i] == 0) break;
     734    Npix *= V1[0].buffer[0].matrix.Naxis[i];
     735  }
     736 
    731737  if (V1[0].type == ST_MATRIX_TMP) {  /** use V1 as temp buffer **/
    732738    OUT[0].buffer = V1[0].buffer;
    … …  
    748754
    749755# define MM_FUNC(OP)                                    \
    750   for (i = 0; i < Nx*Ny; i++, out++, M1++, M2++) {      \
     756  for (i = 0; i < Npix; i++, out++, M1++, M2++) {       \
    751757    *out = OP;                                          \
    752758  }                                                     \
    … …  
    803809int MS_binary (StackVar *OUT, StackVar *V1, StackVar *V2, char *op) {
    804810
    805   int i, Nx, Ny;
     811  int i, Npix;
    806812  char line[512]; // this is only used to report an error
    807813 
    808   Nx = V1[0].buffer[0].matrix.Naxis[0];
    809   Ny = V1[0].buffer[0].matrix.Naxis[1];
     814  Npix = 1;
     815  for (i = 0; i < V1[0].buffer[0].matrix.Naxes; i++) {
     816    if (V1[0].buffer[0].matrix.Naxis[i] == 0) break;
     817    Npix *= V1[0].buffer[0].matrix.Naxis[i];
     818  }
    810819
    811820  /* if possible, use V1 as temp buffer, otherwise create new one */
    … …  
    825834    if (V2->type == ST_SCALAR_FLT)  {                           \
    826835      opihi_flt M2 = V2[0].FltValue;                    \
    827       for (i = 0; i < Nx*Ny; i++, out++, M1++) {        \
     836      for (i = 0; i < Npix; i++, out++, M1++) { \
    828837        *out = OP;                                      \
    829838      }                                                 \
    … …  
    832841    if (V2->type == ST_SCALAR_INT)  {                           \
    833842      opihi_int M2 = V2[0].IntValue;                    \
    834       for (i = 0; i < Nx*Ny; i++, out++, M1++) {        \
     843      for (i = 0; i < Npix; i++, out++, M1++) { \
    835844        *out = OP;                                      \
    836845      }                                                 \
    … …  
    881890int SM_binary (StackVar *OUT, StackVar *V1, StackVar *V2, char *op) {
    882891
    883   int i, Nx, Ny;
     892  int i, Npix;
    884893  char line[512]; // this is only used to report an error
    885894 
    886   Nx = V2[0].buffer[0].matrix.Naxis[0];
    887   Ny = V2[0].buffer[0].matrix.Naxis[1];
    888 
     895  Npix = 1;
     896  for (i = 0; i < V2[0].buffer[0].matrix.Naxes; i++) {
     897    if (V2[0].buffer[0].matrix.Naxis[i] == 0) break;
     898    Npix *= V2[0].buffer[0].matrix.Naxis[i];
     899  }
     900 
    889901  if (V2[0].type == ST_MATRIX_TMP) {  /* V2[0] is NOT temporary, we can't use it for storage */
    890902    OUT[0].buffer = V2[0].buffer;
    … …  
    902914    if (V1->type == ST_SCALAR_FLT)  {                           \
    903915      opihi_flt M1 = V1[0].FltValue;                    \
    904       for (i = 0; i < Nx*Ny; i++, out++, M2++) {        \
     916      for (i = 0; i < Npix; i++, out++, M2++) { \
    905917        *out = OP;                                      \
    906918      }                                                 \
    … …  
    909921    if (V1->type == ST_SCALAR_INT)  {                           \
    910922      opihi_int M1 = V1[0].IntValue;                    \
    911       for (i = 0; i < Nx*Ny; i++, out++, M2++) {        \
     923      for (i = 0; i < Npix; i++, out++, M2++) { \
    912924        *out = OP;                                      \
    913925      }                                                 \
    … …  
    12591271int M_unary (StackVar *OUT, StackVar *V1, char *op) {
    12601272
    1261   int i, j, Nx, Ny;
     1273  int i, j, Npix;
    12621274  float *out, *M1;
    12631275 
    1264   Nx = V1[0].buffer[0].matrix.Naxis[0];
    1265   Ny = V1[0].buffer[0].matrix.Naxis[1];
    1266 
     1276  Npix = 1;
     1277  for (i = 0; i < V1[0].buffer[0].matrix.Naxes; i++) {
     1278    if (V1[0].buffer[0].matrix.Naxis[i] == 0) break;
     1279    Npix *= V1[0].buffer[0].matrix.Naxis[i];
     1280  }
     1281 
    12671282  if (V1[0].type == ST_MATRIX_TMP) {
    12681283    OUT[0].buffer = V1[0].buffer;
    … …  
    12771292
    12781293  if (!strcmp (op, "="))     { }
    1279   if (!strcmp (op, "abs"))   { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = fabs(*M1);         }}
    1280   if (!strcmp (op, "int"))   { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = (opihi_flt)(long long)(*M1); }}
    1281 
    1282   if (!strcmp (op, "floor")) { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = floor (*M1); }}
    1283   if (!strcmp (op, "ceil"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = ceil (*M1); }}
    1284   // if (!strcmp (op, "rint"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = nearbyint (*M1); }}
    1285 
    1286   if (!strcmp (op, "exp"))   { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = exp(*M1);          }}
    1287   if (!strcmp (op, "ten"))   { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = pow(10.0,*M1);     }}
    1288   if (!strcmp (op, "log"))   { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = log10(*M1);        }}
    1289   if (!strcmp (op, "ln"))    { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = log(*M1);          }}
    1290   if (!strcmp (op, "sqrt"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = sqrt(*M1);         }}
    1291   if (!strcmp (op, "erf"))   { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = erf(*M1);          }}
    1292 
    1293   if (!strcmp (op, "sinh"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = sinh(*M1);         }}
    1294   if (!strcmp (op, "cosh"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = cosh(*M1);         }}
    1295   if (!strcmp (op, "asinh")) { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = asinh(*M1);        }}
    1296   if (!strcmp (op, "acosh")) { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = acosh(*M1);        }}
    1297   if (!strcmp (op, "lgamma")) { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = lgamma(*M1);      }}
    1298 
    1299   if (!strcmp (op, "sin"))   { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = sin(*M1);          }}
    1300   if (!strcmp (op, "cos"))   { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = cos(*M1);          }}
    1301   if (!strcmp (op, "tan"))   { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = tan(*M1);          }}
    1302   if (!strcmp (op, "dsin"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = sin(*M1*RAD_DEG);  }}
    1303   if (!strcmp (op, "dcos"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = cos(*M1*RAD_DEG);  }}
    1304   if (!strcmp (op, "dtan"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = tan(*M1*RAD_DEG);  }}
    1305   if (!strcmp (op, "asin"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = asin(*M1);         }}
    1306   if (!strcmp (op, "acos"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = acos(*M1);         }}
    1307   if (!strcmp (op, "atan"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = atan(*M1);         }}
    1308   if (!strcmp (op, "dasin")) { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = asin(*M1)*DEG_RAD; }}
    1309   if (!strcmp (op, "dacos")) { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = acos(*M1)*DEG_RAD; }}
    1310   if (!strcmp (op, "datan")) { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = atan(*M1)*DEG_RAD; }}
    1311   if (!strcmp (op, "not"))   { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = !(*M1);            }}
    1312   if (!strcmp (op, "--"))    { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = -(*M1);            }}
    1313   if (!strcmp (op, "rnd"))   { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = drand48();         }}
    1314   if (!strcmp (op, "ramp"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = i;                 }}
    1315   if (!strcmp (op, "zero"))  { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = 0;                 }}
    1316   if (!strcmp (op, "isinf")) { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = !finite(*M1);      }}
    1317   if (!strcmp (op, "isnan")) { for (i = 0; i < Nx*Ny; i++, out++, M1++) { *out = isnan(*M1);        }}
    1318 
    1319   /* xrm and yrm only make sense in for matrices. see special meaning for vectors */
     1294  if (!strcmp (op, "abs"))   { for (i = 0; i < Npix; i++, out++, M1++) { *out = fabs(*M1);         }}
     1295  if (!strcmp (op, "int"))   { for (i = 0; i < Npix; i++, out++, M1++) { *out = (opihi_flt)(long long)(*M1); }}
     1296
     1297  if (!strcmp (op, "floor")) { for (i = 0; i < Npix; i++, out++, M1++) { *out = floor (*M1); }}
     1298  if (!strcmp (op, "ceil"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = ceil (*M1); }}
     1299  // if (!strcmp (op, "rint"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = nearbyint (*M1); }}
     1300
     1301  if (!strcmp (op, "exp"))   { for (i = 0; i < Npix; i++, out++, M1++) { *out = exp(*M1);          }}
     1302  if (!strcmp (op, "ten"))   { for (i = 0; i < Npix; i++, out++, M1++) { *out = pow(10.0,*M1);     }}
     1303  if (!strcmp (op, "log"))   { for (i = 0; i < Npix; i++, out++, M1++) { *out = log10(*M1);        }}
     1304  if (!strcmp (op, "ln"))    { for (i = 0; i < Npix; i++, out++, M1++) { *out = log(*M1);          }}
     1305  if (!strcmp (op, "sqrt"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = sqrt(*M1);         }}
     1306  if (!strcmp (op, "erf"))   { for (i = 0; i < Npix; i++, out++, M1++) { *out = erf(*M1);          }}
     1307
     1308  if (!strcmp (op, "sinh"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = sinh(*M1);         }}
     1309  if (!strcmp (op, "cosh"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = cosh(*M1);         }}
     1310  if (!strcmp (op, "asinh")) { for (i = 0; i < Npix; i++, out++, M1++) { *out = asinh(*M1);        }}
     1311  if (!strcmp (op, "acosh")) { for (i = 0; i < Npix; i++, out++, M1++) { *out = acosh(*M1);        }}
     1312  if (!strcmp (op, "lgamma")) { for (i = 0; i < Npix; i++, out++, M1++) { *out = lgamma(*M1);      }}
     1313
     1314  if (!strcmp (op, "sin"))   { for (i = 0; i < Npix; i++, out++, M1++) { *out = sin(*M1);          }}
     1315  if (!strcmp (op, "cos"))   { for (i = 0; i < Npix; i++, out++, M1++) { *out = cos(*M1);          }}
     1316  if (!strcmp (op, "tan"))   { for (i = 0; i < Npix; i++, out++, M1++) { *out = tan(*M1);          }}
     1317  if (!strcmp (op, "dsin"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = sin(*M1*RAD_DEG);  }}
     1318  if (!strcmp (op, "dcos"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = cos(*M1*RAD_DEG);  }}
     1319  if (!strcmp (op, "dtan"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = tan(*M1*RAD_DEG);  }}
     1320  if (!strcmp (op, "asin"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = asin(*M1);         }}
     1321  if (!strcmp (op, "acos"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = acos(*M1);         }}
     1322  if (!strcmp (op, "atan"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = atan(*M1);         }}
     1323  if (!strcmp (op, "dasin")) { for (i = 0; i < Npix; i++, out++, M1++) { *out = asin(*M1)*DEG_RAD; }}
     1324  if (!strcmp (op, "dacos")) { for (i = 0; i < Npix; i++, out++, M1++) { *out = acos(*M1)*DEG_RAD; }}
     1325  if (!strcmp (op, "datan")) { for (i = 0; i < Npix; i++, out++, M1++) { *out = atan(*M1)*DEG_RAD; }}
     1326  if (!strcmp (op, "not"))   { for (i = 0; i < Npix; i++, out++, M1++) { *out = !(*M1);            }}
     1327  if (!strcmp (op, "--"))    { for (i = 0; i < Npix; i++, out++, M1++) { *out = -(*M1);            }}
     1328  if (!strcmp (op, "rnd"))   { for (i = 0; i < Npix; i++, out++, M1++) { *out = drand48();         }}
     1329  if (!strcmp (op, "ramp"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = i;                 }}
     1330  if (!strcmp (op, "zero"))  { for (i = 0; i < Npix; i++, out++, M1++) { *out = 0;                 }}
     1331  if (!strcmp (op, "isinf")) { for (i = 0; i < Npix; i++, out++, M1++) { *out = !finite(*M1);      }}
     1332  if (!strcmp (op, "isnan")) { for (i = 0; i < Npix; i++, out++, M1++) { *out = isnan(*M1);        }}
     1333
     1334  /* xrm and yrm only make sense for 2D matrices. see special meaning for vectors */
    13201335  if (!strcmp (op, "xramp")) {
     1336    int Nx = V1[0].buffer[0].matrix.Naxis[0];
     1337    int Ny = V1[0].buffer[0].matrix.Naxis[1];
    13211338    for (j = 0; j < Ny; j++) {
    13221339      for (i = 0; i < Nx; i++, out++, M1++) {
    … …  
    13261343  }
    13271344  if (!strcmp (op, "yramp")) {
     1345    int Nx = V1[0].buffer[0].matrix.Naxis[0];
     1346    int Ny = V1[0].buffer[0].matrix.Naxis[1];
    13281347    for (j = 0; j < Ny; j++) {
    13291348      for (i = 0; i < Nx; i++, out++, M1++) {
Note: See TracChangeset for help on using the changeset viewer.