Changeset 26467 for branches/eam_branches/20091201/Ohana/src/opihi/cmd.data/hermitian2d.c

Timestamp:

Dec 22, 2009, 9:52:30 AM (17 years ago)

Author:

eugene

Message:

adding hermitian polynomials and functions

File:

• branches/eam_branches/20091201/Ohana/src/opihi/cmd.data/hermitian2d.c (modified) (4 diffs)

branches/eam_branches/20091201/Ohana/src/opihi/cmd.data/hermitian2d.c

@@ -3 +3 @@
 int hermitian2d (int argc, char **argv) {
   
-  int i, j, Nx, Ny, N;
+  int i, j, N, Nx, Ny, Xorder, Yorder, Poly;
   float *in;
-  double Sig_x, Sig_y, Theta;
-  double root1, root2, R, A1, A2, A3;
-  double Sx, Sy, Sxy;
-  double x, y, r, f, Xo, Yo;
+  double Xo, Yo, sigma;
+  double x, y, xf, yf, value, norm, weight;
+  double nf, mf;
   Buffer *buf;
 
-  Poly = FALSE;
+  // optionally return the Hermitian Polynomial or the Hermitian Function
+  Poly = FALSE; 
   if ((N = get_argument (argc, argv, "-poly"))) {
     Poly = TRUE;
@@ -26 +26 @@
   }
 
-  if (argc > 5) {
-    gprint (GP_ERR, "USAGE: hermitian2d (buffer) (sigma) Xorder Yorder\n");
-    return (FALSE);
-  }
+  if (argc != 5) goto usage;
 
   /* select input / output buffers */
@@ -43 +40 @@
   Xorder = atof (argv[3]);
   Yorder = atof (argv[4]);
+  if (Xorder < 0) goto usage;
+  if (Yorder < 0) goto usage;
+  if (Xorder > 10) goto usage;
+  if (Yorder > 10) goto usage;
+
+  nf = exp(lgamma(Xorder + 1));
+  mf = exp(lgamma(Yorder + 1));
+  norm = 1.0 / sqrt(nf*mf*2*M_PI) / sigma;
 
   in = (float *) buf[0].matrix.buffer;
@@ -58 +63 @@
       value = xf*yf;
       if (!Poly) {
-        value *= xf*yf;
+        weight = norm*exp(-0.25*(x*x + y*y));
+        value *= weight;
       }
 
-      *in += xf*yf;
+      *in += value;
     }
   }
 
   return (TRUE);
+
+usage:
+  gprint (GP_ERR, "USAGE: hermitian2d (buffer) (sigma) Xorder Yorder\n");
+  gprint (GP_ERR, "  Note : only orders up to 10 are supported\n");
+  return (FALSE);
 }