r520: no message
[ctsim.git] / libctsim / filter.cpp
index a88fad5a51c3cb1a663908a81b22dc95f8f7d07b..f4a10f3186dbbeb6da08a3fbdc810ee11a4b8e27 100644 (file)
@@ -7,9 +7,9 @@
 **     Date Started:           Aug 1984
 **
 **  This is part of the CTSim program
-**  Copyright (C) 1983-2000 Kevin Rosenberg
+**  Copyright (c) 1983-2000 Kevin Rosenberg
 **
-**  $Id: filter.cpp,v 1.25 2000/08/19 22:59:06 kevin Exp $
+**  $Id: filter.cpp,v 1.36 2001/02/11 04:56:37 kevin Exp $
 **
 **  This program is free software; you can redistribute it and/or modify
 **  it under the terms of the GNU General Public License (version 2) as
 int SignalFilter::N_INTEGRAL=500;  //static member
 
 const int SignalFilter::FILTER_INVALID = -1 ;
-const int SignalFilter::FILTER_ABS_BANDLIMIT = 0;      // filter times |x = |
-const int SignalFilter::FILTER_ABS_SINC = 1;
-const int SignalFilter::FILTER_ABS_G_HAMMING = 2;
-const int SignalFilter::FILTER_ABS_COSINE = 3;
+const int SignalFilter::FILTER_ABS_BANDLIMIT = 0;      // filter times |x|
+const int SignalFilter::FILTER_ABS_G_HAMMING = 1;
+const int SignalFilter::FILTER_ABS_COSINE = 2;
+const int SignalFilter::FILTER_ABS_SINC = 3;
 const int SignalFilter::FILTER_SHEPP = 4;
 const int SignalFilter::FILTER_BANDLIMIT = 5;
 const int SignalFilter::FILTER_SINC = 6;
@@ -41,11 +41,13 @@ const int SignalFilter::FILTER_G_HAMMING = 7;
 const int SignalFilter::FILTER_COSINE = 8;
 const int SignalFilter::FILTER_TRIANGLE = 9;
 
+const int SignalFilter::s_iReconstructFilterCount = 5;
+
 const char* SignalFilter::s_aszFilterName[] = {
   {"abs_bandlimit"},
-  {"abs_sinc"},
   {"abs_hamming"},
   {"abs_cosine"},
+  {"abs_sinc"},
   {"shepp"},
   {"bandlimit"},
   {"sinc"},
@@ -56,9 +58,9 @@ const char* SignalFilter::s_aszFilterName[] = {
 
 const char* SignalFilter::s_aszFilterTitle[] = {
   {"Abs(w) * Bandlimit"},
-  {"Abs(w) * Sinc"},
   {"Abs(w) * Hamming"},
   {"Abs(w) * Cosine"},
+  {"Abs(w) * Sinc"},
   {"Shepp"},
   {"Bandlimit"},
   {"Sinc"},
@@ -364,12 +366,13 @@ SignalFilter::spatialResponseCalc (int filterID, double bw, double x, double par
   double zinc = (zmax - zmin) / (n - 1);
 
   double z = zmin;
-  double q [n];
+  double* q = new double [n];
   for (int i = 0; i < n; i++, z += zinc)
     q[i] = frequencyResponse (filterID, bw, z, param) * cos (TWOPI * z * x);
   
   double y = 2 * integrateSimpson (zmin, zmax, q, n);
-  
+  delete q;
+
   return (y);
 }
 
@@ -401,31 +404,31 @@ SignalFilter::frequencyResponse (int filterID, double bw, double u, double param
 
   switch (filterID) {
   case FILTER_BANDLIMIT:
-    if (au >= bw / 2)
+    if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
       q = 0.;
     else
       q = 1;
     break;
   case FILTER_ABS_BANDLIMIT:
-    if (au >= bw / 2)
+    if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
       q = 0.;
     else
       q = au;
     break;
   case FILTER_TRIANGLE:
-    if (au >= bw)
+    if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
       q = 0;
     else
       q = 1 - au / bw;
     break;
   case FILTER_COSINE:
-    if (au >= bw / 2)
+    if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
       q = 0;
     else
       q = cos(PI * u / bw);
     break;
   case FILTER_ABS_COSINE:
-    if (au >= bw / 2)
+    if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
       q = 0;
     else
       q = au * cos(PI * u / bw);
@@ -437,13 +440,13 @@ SignalFilter::frequencyResponse (int filterID, double bw, double u, double param
     q = au * bw * sinc (PI * bw * u, 1.);
     break;
   case FILTER_G_HAMMING:
-    if (au >= bw / 2)
+    if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
       q = 0;
     else
       q = param + (1 - param) * cos (TWOPI * u / bw);
     break;
   case FILTER_ABS_G_HAMMING:
-    if (au >= bw / 2)
+    if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
       q = 0;
     else
       q = au * (param + (1 - param) * cos(TWOPI * u / bw));
@@ -558,30 +561,21 @@ SignalFilter::spatialResponseAnalytic (int filterID, double bw, double x, double
 }
 
 
-/* NAME
- *   sinc                      Return sin(x)/x function
- *
- * SYNOPSIS
- *   v = sinc (x, mult)
- *   double v                  sinc value
- *   double x, mult
- *
- * DESCRIPTION
- *   v = sin(x * mult) / x;
- */
 
+// Functions that are inline in filter.h
+
+
+//  sinc                       Return sin(x)/x function 
+//   v = sinc (x, mult)
+// Calculates sin(x * mult) / x;
+
+//  integral_abscos    Returns integral of u*cos(u)
+//
+//   q = integral_abscos (u, w)
+//   double q                  Integral value
+//   double u                  Integration variable
+//   double w                  Upper integration boundary
+// Returns the value of integral of u*cos(u)*dV for V = 0 to w
 
-/* NAME
- *   integral_abscos                   Returns integral of u*cos(u)
- *
- * SYNOPSIS
- *   q = integral_abscos (u, w)
- *   double q                  Integral value
- *   double u                  Integration variable
- *   double w                  Upper integration boundary
- *
- * DESCRIPTION
- *   Returns the value of integral of u*cos(u)*dV for V = 0 to w
- */