r186: *** empty log message ***

diff --git a/ChangeLog b/ChangeLog
index bd59494c285d850eacb3349767e53211b667b42e..e77b5944d954b29da019af3d7d6b4bd9f0d8432a 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,10 @@
+2.0.0-b11 - 8/31/00
+   Added reconstruction for Equilinear and Equiangular geometries
+   Changed theta to be -PI/2 to make compliant with Kak-Slaney formulas
+   ctsim: Added projection graph to animation of projection collection 
+   ctsim: Added single stepping to projection collection animation
+   ctsim: improved File/Properties display for projection files
+       

 2.0.0-b10 - 8/25/00
    ctsim: Added animation of projection collection processs
    ctsim: Added Auto Scaling for image windows
 2.0.0-b10 - 8/25/00
    ctsim: Added animation of projection collection processs
    ctsim: Added Auto Scaling for image windows

diff --git a/acconfig.h b/acconfig.h
index 9c06cd55c1fec6d76f7360aab0fe96c9c1ae7243..c7dffb528072c7e94951d69b1658ad3249e021e2 100644 (file)
--- a/acconfig.h
+++ b/acconfig.h
@@ -74,3 +74,5 @@
 
 /* Define for HDF5 library */
 #undef HAVE_HDF5
 
 /* Define for HDF5 library */
 #undef HAVE_HDF5

+
+#undef HAVE_USLEEP

diff --git a/config.h.in b/config.h.in
index 31cc9fd0d134d92aca5e950c5471e6d0984287e3..97f445e22721b277981b20360b0f1c06854a7d42 100644 (file)
--- a/config.h.in
+++ b/config.h.in
@@ -90,6 +90,9 @@
 /* Define if you have the strtol function.  */
 #undef HAVE_STRTOL
 
 /* Define if you have the strtol function.  */
 #undef HAVE_STRTOL
 

+/* Define if you have the usleep function.  */
+#undef HAVE_USLEEP
+

 /* Define if you have the <X11/extensions/Xdbe.h> header file.  */
 #undef HAVE_X11_EXTENSIONS_XDBE_H
 
 /* Define if you have the <X11/extensions/Xdbe.h> header file.  */
 #undef HAVE_X11_EXTENSIONS_XDBE_H
 

diff --git a/configure b/configure
index b46197cbe668229da4d8b7897a998072b53eac57..299bf4348756e26f01cf987ea1397049503d315a 100755 (executable)
--- a/configure
+++ b/configure
@@ -712,7 +712,7 @@ fi
 
 PACKAGE=ctsim
 
 
 PACKAGE=ctsim
 

-VERSION=2.0.0b10
+VERSION=2.0.0b11

 
 if test "`cd $srcdir && pwd`" != "`pwd`" && test -f $srcdir/config.status; then
   { echo "configure: error: source directory already configured; run "make distclean" there first" 1>&2; exit 1; }
 
 if test "`cd $srcdir && pwd`" != "`pwd`" && test -f $srcdir/config.status; then
   { echo "configure: error: source directory already configured; run "make distclean" there first" 1>&2; exit 1; }


@@ -2614,7 +2614,7 @@ fi
 
 fi
 
 
 fi
 

-for ac_func in strtod strtol snprintf htonl
+for ac_func in strtod strtol snprintf htonl usleep

 do
 echo $ac_n "checking for $ac_func""... $ac_c" 1>&6
 echo "configure:2621: checking for $ac_func" >&5
 do
 echo $ac_n "checking for $ac_func""... $ac_c" 1>&6
 echo "configure:2621: checking for $ac_func" >&5

diff --git a/configure.in b/configure.in
index 09d8ec94d7cfb31bb1830543c66c6aa9a5b80881..0b5e2705f8e4cfa53a4a1e4eed87603375630ea1 100644 (file)
--- a/configure.in
+++ b/configure.in
@@ -4,7 +4,7 @@ dnl Must reset CDPATH so that bash's cd does not print to stdout
 dnl CDPATH=
 
 AC_INIT(src/ctsim.cpp)
 dnl CDPATH=
 
 AC_INIT(src/ctsim.cpp)

-AM_INIT_AUTOMAKE(ctsim,2.0.0b10)
+AM_INIT_AUTOMAKE(ctsim,2.0.0b11)

 AM_CONFIG_HEADER(config.h)
 
 dnl Checks for programs.
 AM_CONFIG_HEADER(config.h)
 
 dnl Checks for programs.


@@ -99,7 +99,7 @@ AC_STRUCT_TM
 
 dnl Checks for library functions.
 AC_FUNC_VPRINTF
 
 dnl Checks for library functions.
 AC_FUNC_VPRINTF

-AC_CHECK_FUNCS(strtod strtol snprintf htonl)
+AC_CHECK_FUNCS(strtod strtol snprintf htonl usleep)

 AC_CHECK_FUNC(basename)
 AC_CHECK_FUNC(setjmp)
 if test "${OSTYPE}" = "cygwin" ; then
 AC_CHECK_FUNC(basename)
 AC_CHECK_FUNC(setjmp)
 if test "${OSTYPE}" = "cygwin" ; then

diff --git a/html/simulate.html.in b/html/simulate.html.in
index 5039901700e4bc47beb0faece40cf40d2646d1b8..bd0745a34768c4510e0d9993c56888f9731b8e10 100644 (file)
--- a/html/simulate.html.in
+++ b/html/simulate.html.in
@@ -28,7 +28,12 @@ MPI Supercomputing:<br>
 <INPUT type="radio" name="MPI" value="no" checked>No (Single CPU)<br>
 </td><td valign="top">
 <h3>Simulate X-Ray acquistion</h3>
 <INPUT type="radio" name="MPI" value="no" checked>No (Single CPU)<br>
 </td><td valign="top">
 <h3>Simulate X-Ray acquistion</h3>

-Number of detectors: <input type="text" name="PJ_NDet" size="4" value="367"><p>
+Geometry:<br>
+<INPUT type="radio" name="geometry" value="parallel" checked>Parallel<br>
+<INPUT type="radio" name="geometry" value="equiangular">Equiangular<br>
+<INPUT type="radio" name="geometry" value="collinear">Collinear<br>
+<p>
+Number of Detectors: <input type="text" name="PJ_NDet" size="4" value="367"><p>

 Number of Rotations: <input type="text" name="PJ_NRot" size="4" value="320"><p>
 Number of Rays<br>(samples) per detector: <input type="text" name="PJ_NRay" size="2" value="3"><p>
 Rotation Angle<br>as a multiple of PI: <input type="text" name="PJ_RotAngle" size="3" value="1.0"><br>
 Number of Rotations: <input type="text" name="PJ_NRot" size="4" value="320"><p>
 Number of Rays<br>(samples) per detector: <input type="text" name="PJ_NRay" size="2" value="3"><p>
 Rotation Angle<br>as a multiple of PI: <input type="text" name="PJ_RotAngle" size="3" value="1.0"><br>

diff --git a/include/backprojectors.h b/include/backprojectors.h
index 1781fe18b5f2c60a519edd6153eaade5c7c3f755..b36dfb83cefbb7ef776ef5d4e2da3871e378a9ff 100644 (file)
--- a/include/backprojectors.h
+++ b/include/backprojectors.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: backprojectors.h,v 1.13 2000/08/25 15:59:13 kevin Exp $
+**  $Id: backprojectors.h,v 1.14 2000/08/31 08:38:58 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
 **
 **  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


@@ -132,6 +132,7 @@ class Backproject
     double xMin, xMax, yMin, yMax;     // Retangular coords of phantom
     double xInc, yInc; // size of cells
     int m_interpFactor;
     double xMin, xMax, yMin, yMax;     // Retangular coords of phantom
     double xInc, yInc; // size of cells
     int m_interpFactor;

+    double m_dFocalLength;

 
  private:
     Backproject (const Backproject& rhs);
 
  private:
     Backproject (const Backproject& rhs);


@@ -158,7 +159,7 @@ class BackprojectTable : public Backproject
 
   void BackprojectView (const double* const t, double view_angle);
 
 
   void BackprojectView (const double* const t, double view_angle);
 

- private:
+ protected:

   Array2d<kfloat64> arrayR;
   Array2d<kfloat64> arrayPhi;
   kfloat64** r;
   Array2d<kfloat64> arrayR;
   Array2d<kfloat64> arrayPhi;
   kfloat64** r;


@@ -211,5 +212,25 @@ class BackprojectIntDiff3 : public BackprojectDiff
   void BackprojectView (const double* const t, double view_angle);
 };
 
   void BackprojectView (const double* const t, double view_angle);
 };
 

+class BackprojectEquilinear : public BackprojectTable
+{
+ public:
+  BackprojectEquilinear (const Projections& proj, ImageFile& im, int interpID, const int interpFactor)
+      : BackprojectTable::BackprojectTable (proj, im, interpID, interpFactor)
+      {}
+
+  void BackprojectView (const double* const t, double view_angle);
+};
+
+class BackprojectEquiangular : public BackprojectTable
+{
+ public:
+  BackprojectEquiangular (const Projections& proj, ImageFile& im, int interpID, const int interpFactor)
+      : BackprojectTable::BackprojectTable (proj, im, interpID, interpFactor)
+      {}
+
+  void BackprojectView (const double* const t, double view_angle);
+};
+

 
 #endif
 
 #endif

diff --git a/include/ezplot.h b/include/ezplot.h
index 22c2360873ad0ead0dc1d38032fa5b072fe9c269..716fd3e2b361c8946448a7ecb5ea12d5edcbac6c 100644 (file)
--- a/include/ezplot.h
+++ b/include/ezplot.h
@@ -7,7 +7,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: ezplot.h,v 1.12 2000/07/29 19:50:08 kevin Exp $
+**  $Id: ezplot.h,v 1.13 2000/08/31 08:38:58 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
 **
 **  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


@@ -221,6 +221,7 @@ class EZPlot {
     void addCurve (const double* x, const float* y, int num);
     void addCurve (const double* x, const double* y, int num);
     void addCurve (const double* y, int n);
     void addCurve (const double* x, const float* y, int num);
     void addCurve (const double* x, const double* y, int num);
     void addCurve (const double* y, int n);

+    void addCurve (const float* y, int n);

 
     void plot ();
 };     
 
     void plot ();
 };     

diff --git a/include/procsignal.h b/include/procsignal.h
index e8e9679e3d247d814d601ac1518002f2ea6e922d..ad2858161a4e1158522c9a22c78eb24076974cd3 100644 (file)
--- a/include/procsignal.h
+++ b/include/procsignal.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: procsignal.h,v 1.3 2000/08/27 20:32:54 kevin Exp $
+**  $Id: procsignal.h,v 1.4 2000/08/31 08:38:58 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
 **
 **  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


@@ -58,14 +58,11 @@ class ProcessSignal {
     static const int FILTER_GENERATION_DIRECT;
     static const int FILTER_GENERATION_INVERSE_FOURIER;
 
     static const int FILTER_GENERATION_DIRECT;
     static const int FILTER_GENERATION_INVERSE_FOURIER;
 

-    ProcessSignal (const char* szFilterName, const char* szFilterMethodName,double bw, double signalIncrement, int n, double param, const char* szDomainName, const char* szFilterGenerationName, const int zeropad = 0, const int preinterpolationFactor = 1, const int iTraceLevel = Trace::TRACE_NONE, int iGeometry = Scanner::GEOMETRY_PARALLEL);
-
-    ProcessSignal (const int idFilter, int idFilterMethod, double dBandwidth, double dSignalIncrement, int nSignalPoints, double dFilterParam, const int idDomain, int idFilterGeneration, const int iZeropad = 0, const int iPreinterpolationFactor = 1, const int iTraceLevel = Trace::TRACE_NONE);
+    ProcessSignal (const char* szFilterName, const char* szFilterMethodName,double bw, double signalIncrement, int n, double param, const char* szDomainName, const char* szFilterGenerationName, const int zeropad = 0, const int preinterpolationFactor = 1, const int iTraceLevel = Trace::TRACE_NONE, int iGeometry = Scanner::GEOMETRY_PARALLEL, double dFocalLength = 1.);

 
     ~ProcessSignal();
 
 
     ~ProcessSignal();
 

-    void filterSignal (const double input[], double output[], int view =0) const;
-    void filterSignal (const float input[], double output[], int view = 0) const;
+    void filterSignal (const float input[], double output[]) const;

 
     bool fail(void) const      {return m_fail;}
     const string& failMessage(void) const {return m_failMessage;}
 
     bool fail(void) const      {return m_fail;}
     const string& failMessage(void) const {return m_failMessage;}


@@ -133,6 +130,7 @@ class ProcessSignal {
     int m_nOutputPoints;
     int m_iPreinterpolationFactor;
     int m_idGeometry;
     int m_nOutputPoints;
     int m_iPreinterpolationFactor;
     int m_idGeometry;

+    double m_dFocalLength;

 
     bool m_fail;
     string m_failMessage;
 
     bool m_fail;
     string m_failMessage;


@@ -151,7 +149,7 @@ class ProcessSignal {
     fftw_plan m_complexPlanForward, m_complexPlanBackward;
 #endif
 
     fftw_plan m_complexPlanForward, m_complexPlanBackward;
 #endif
 

-    void init (const int idFilter, int idFilterMethod, double dBandwidth, double dSignalIncrement, int nSignalPoints, double dFilterParam, const int idDomain, int idFilterGeneration, const int iZeropad, const int iPreinterpolationFactor, const int iTraceLevel, const int iGeometry);
+    void init (const int idFilter, int idFilterMethod, double dBandwidth, double dSignalIncrement, int nSignalPoints, double dFilterParam, const int idDomain, int idFilterGeneration, const int iZeropad, const int iPreinterpolationFactor, const int iTraceLevel, const int iGeometry, double dFocalLength);

 
     // transforms that use precalculated trig tables, therefore don't 
     // require number of data points (n) as an argument
 
     // transforms that use precalculated trig tables, therefore don't 
     // require number of data points (n) as an argument


@@ -159,6 +157,7 @@ class ProcessSignal {
     void finiteFourierTransform (const complex<double> input[], complex<double> output[], const int direction) const;
     void finiteFourierTransform (const complex<double> input[], double output[], const int direction) const;
 
     void finiteFourierTransform (const complex<double> input[], complex<double> output[], const int direction) const;
     void finiteFourierTransform (const complex<double> input[], double output[], const int direction) const;
 

+    double convolve (const double func[], const double filter[], const double dx, const int n, const int np) const;

     double convolve (const double f[], const double dx, const int n, const int np) const;
     double convolve (const float f[], const double dx, const int n, const int np) const;
 
     double convolve (const double f[], const double dx, const int n, const int np) const;
     double convolve (const float f[], const double dx, const int n, const int np) const;
 

diff --git a/include/projections.h b/include/projections.h
index 38828a1532428590a7f6d28e854921453f39b78b..8a916166170327ee92adcde5dffe93d6b7993c21 100644 (file)
--- a/include/projections.h
+++ b/include/projections.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: projections.h,v 1.12 2000/08/27 20:32:55 kevin Exp $
+**  $Id: projections.h,v 1.13 2000/08/31 08:38:58 kevin Exp $

 **
 **
 **  This program is free software; you can redistribute it and/or modify
 **
 **
 **  This program is free software; you can redistribute it and/or modify


@@ -32,7 +32,7 @@
 class Scanner;
 class DetectorArray;
 class Array2dFileLabel;
 class Scanner;
 class DetectorArray;
 class Array2dFileLabel;

-
+class ostringstream;

 
 // Projections
 class Projections
 
 // Projections
 class Projections


@@ -40,16 +40,18 @@ class Projections
  public:
   Projections (const int nView, const int nDet);
   Projections (const Scanner& scanner);
  public:
   Projections (const int nView, const int nDet);
   Projections (const Scanner& scanner);

-  Projections (void);
-  ~Projections (void);
+  Projections ();
+  ~Projections ();

 
   void initFromScanner (const Scanner& scanner);
 
 
   void initFromScanner (const Scanner& scanner);
 

-  void printProjectionData (void);
-  void printScanInfo (void) const;
+  void printProjectionData ();
+  void printScanInfo (ostringstream& os) const;

 
 

+  bool read (const string& fname);

   bool read (const char* fname);
   bool write (const char* fname);
   bool read (const char* fname);
   bool write (const char* fname);

+  bool write (const string& fname);

   bool detarrayRead (fnetorderstream& fs, DetectorArray& darray, const int view_num);
   bool detarrayWrite (fnetorderstream& fs, const DetectorArray& darray, const int view_num);
 
   bool detarrayRead (fnetorderstream& fs, DetectorArray& darray, const int view_num);
   bool detarrayWrite (fnetorderstream& fs, const DetectorArray& darray, const int view_num);
 


@@ -61,20 +63,24 @@ class Projections
   void setPhmLen (double phmLen) { m_phmLen = phmLen;}
   void setCalcTime (double calcTime) {m_calcTime = calcTime;}
   void setRemark (const char* remark) {m_remark = remark;}
   void setPhmLen (double phmLen) { m_phmLen = phmLen;}
   void setCalcTime (double calcTime) {m_calcTime = calcTime;}
   void setRemark (const char* remark) {m_remark = remark;}

-  void setRemark (const string remark) {m_remark = remark;}
-
-  const double detStart(void) const {return m_detStart;}
-  const double rotStart(void) const {return m_rotStart;}
-  const double calcTime(void) const {return m_calcTime;}
-  const double phmLen(void) const {return m_phmLen;}
-  const char*  remark(void) const {return m_remark.c_str();}
-  const double detInc(void) const {return m_detInc;}
-  const double rotInc(void) const {return m_rotInc;}
-  const int nDet(void) const {return m_nDet;}
-  const int nView(void) const {return m_nView;}
-  const string& getFilename(void) const {return m_filename;}
-  Array2dFileLabel& getLabel(void) {return m_label;}
-  const Array2dFileLabel& getLabel(void) const {return m_label;}
+  void setRemark (const string& remark) {m_remark = remark;}
+
+  double detStart() const {return m_detStart;}
+  double rotStart() const {return m_rotStart;}
+  double calcTime() const {return m_calcTime;}
+  double phmLen() const {return m_phmLen;}
+  const char*  remark() const {return m_remark.c_str();}
+  double detInc() const {return m_detInc;}
+  double rotInc() const {return m_rotInc;}
+  int nDet() const {return m_nDet;}
+  int nView() const {return m_nView;}
+  int geometry() const {return m_geometry;}
+  double focalLength() const {return m_focalLength;}
+  double fieldOfView() const {return m_fieldOfView;}
+
+  const string& getFilename() const {return m_filename;}
+  Array2dFileLabel& getLabel() {return m_label;}
+  const Array2dFileLabel& getLabel() const {return m_label;}

 
   DetectorArray& getDetectorArray (const int iview)
       { return (*m_projData[iview]); }
 
   DetectorArray& getDetectorArray (const int iview)
       { return (*m_projData[iview]); }


@@ -108,12 +114,12 @@ class Projections
 
   const static kuint16 m_signature = ('P'*256 + 'J');
 
 
   const static kuint16 m_signature = ('P'*256 + 'J');
 

-  bool headerRead (void);
-  bool headerWrite (void);
+  bool headerRead ();
+  bool headerWrite ();

   bool headerRead (fnetorderstream& fs);
   bool headerWrite (fnetorderstream& fs);
   bool headerRead (fnetorderstream& fs);
   bool headerWrite (fnetorderstream& fs);

-  void newProjData (void);
-  void deleteProjData (void);
+  void newProjData ();
+  void deleteProjData ();

 
   void init (const int nView, const int nDet);
 
 
   void init (const int nView, const int nDet);
 

diff --git a/include/scanner.h b/include/scanner.h
index fa6ea1d8633010b7e26a61e473d724cf524f1329..7f8f95e7e2aafcb29a4e8323167b0c66bf91bfb6 100644 (file)
--- a/include/scanner.h
+++ b/include/scanner.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: scanner.h,v 1.13 2000/08/27 20:32:55 kevin Exp $
+**  $Id: scanner.h,v 1.14 2000/08/31 08:38:58 kevin Exp $

 **
 **
 **  This program is free software; you can redistribute it and/or modify
 **
 **
 **  This program is free software; you can redistribute it and/or modify


@@ -93,6 +93,7 @@ class Scanner
   double detLen() const {return m_detLen;}
   double focalLength() const {return m_dFocalLength;}
   double fieldOfView() const {return m_dFieldOfView;}
   double detLen() const {return m_detLen;}
   double focalLength() const {return m_dFocalLength;}
   double fieldOfView() const {return m_dFieldOfView;}

+  int geometry() const {return m_idGeometry;}

 
   static int getGeometryCount() {return s_iGeometryCount;}
   static const char** getGeometryNameArray() {return s_aszGeometryName;}
 
   static int getGeometryCount() {return s_iGeometryCount;}
   static const char** getGeometryNameArray() {return s_aszGeometryName;}


@@ -119,12 +120,15 @@ class Scanner
   double m_phmLen;              // Maximum Length of phantom or area of interest
   double m_dXCenter;            // Center of Phantom
   double m_dYCenter;            
   double m_phmLen;              // Maximum Length of phantom or area of interest
   double m_dXCenter;            // Center of Phantom
   double m_dYCenter;            

+  double m_dAngularDetIncrement;
+  double m_dAngularDetLen;

 
   int m_trace;
   struct {
     double xd1,yd1,xd2,yd2;    /* Coordinates of detector endpoints */
     double xs1,ys1,xs2,ys2;    /* Coordinates of source endpoints */
     double angle;              /* Starting angle */
 
   int m_trace;
   struct {
     double xd1,yd1,xd2,yd2;    /* Coordinates of detector endpoints */
     double xs1,ys1,xs2,ys2;    /* Coordinates of source endpoints */
     double angle;              /* Starting angle */

+    double dAngularDet;

   } m_initPos;
 
   GRFMTX_2D m_rotmtxIncrement;
   } m_initPos;
 
   GRFMTX_2D m_rotmtxIncrement;

diff --git a/include/sgp.h b/include/sgp.h
index 0a5fe0344ea7a9c9bd5359d750a7106f703b1a08..b026cdc6f1634cdcb67f2b08fae8b4633d28406a 100644 (file)
--- a/include/sgp.h
+++ b/include/sgp.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: sgp.h,v 1.15 2000/08/27 20:32:55 kevin Exp $
+**  $Id: sgp.h,v 1.16 2000/08/31 08:38:58 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
 **
 **  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


@@ -72,7 +72,7 @@ public:
 #endif
 
   SGPDriver (const char* szWinTitle = "", int xsize = 640, int ysize = 480);
 #endif
 
   SGPDriver (const char* szWinTitle = "", int xsize = 640, int ysize = 480);

-
+  

   ~SGPDriver ();
 
   int getPhysicalXSize () const
   ~SGPDriver ();
 
   int getPhysicalXSize () const


@@ -96,6 +96,9 @@ public:
 #ifdef HAVE_WXWINDOWS
   wxDC* idWX () const
     { return m_pDC; }
 #ifdef HAVE_WXWINDOWS
   wxDC* idWX () const
     { return m_pDC; }

+
+  void setDC (wxDC* dc)
+      { m_pDC = dc; }

 #endif
 };
 
 #endif
 };
 


@@ -105,7 +108,7 @@ class SGP {
 private:
   int m_iPhysicalXSize;   // Physical Window size 
   int m_iPhysicalYSize;
 private:
   int m_iPhysicalXSize;   // Physical Window size 
   int m_iPhysicalYSize;

-  const SGPDriver m_driver;
+  SGPDriver m_driver;

 
   double xw_min;    // Window extents 
   double yw_min;
 
   double xw_min;    // Window extents 
   double yw_min;


@@ -199,6 +202,7 @@ public:
   void setMarker (int idMarker, int color);
   void setRasterOp (int ro);
 
   void setMarker (int idMarker, int color);
   void setRasterOp (int ro);
 

+  void getViewport (double &xmin, double &ymin, double& xmax, double& ymax);

   void getTextExtent (const char *szText, double* x, double* y);
   double getCharHeight ();
 
   void getTextExtent (const char *szText, double* x, double* y);
   double getCharHeight ();
 


@@ -219,6 +223,10 @@ public:
   void stylusNDC (double x, double y, bool beam);
   void pointNDC (double x, double y);
   void markerNDC (double x, double y);
   void stylusNDC (double x, double y, bool beam);
   void pointNDC (double x, double y);
   void markerNDC (double x, double y);

+
+#if HAVE_WXWINDOWS
+  void setDC (wxDC* pDC);
+#endif

 };
 
 
 };
 
 

diff --git a/libctgraphics/ezplot.cpp b/libctgraphics/ezplot.cpp
index 7038132c128db2a6ed0f4ef68457676594aadbda..c29575737f47870ca14c934b2c7e218b8236405e 100644 (file)
--- a/libctgraphics/ezplot.cpp
+++ b/libctgraphics/ezplot.cpp
@@ -6,7 +6,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: ezplot.cpp,v 1.11 2000/08/22 07:02:48 kevin Exp $
+**  $Id: ezplot.cpp,v 1.12 2000/08/31 08:38:58 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
 **
 **  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


@@ -71,6 +71,18 @@ EZPlot::addCurve (const double *y, int n)
 }
 
 
 }
 
 

+void 
+EZPlot::addCurve (const float *y, int n)
+{
+  double yDouble [n];
+
+  for (int i = 0; i < n; i++)
+    yDouble[i] = y[i];
+
+  addCurve (yDouble, n);
+}
+
+

 void
 EZPlot::addCurve (const float x[], const double y[], int num)
 {
 void
 EZPlot::addCurve (const float x[], const double y[], int num)
 {


@@ -300,7 +312,8 @@ EZPlot::plot ()
   xtl_wid = x_fldwid * charwidth;              /* calc size of tick labels */
   ytl_wid = y_fldwid * charwidth;
   tl_height = charheight;
   xtl_wid = x_fldwid * charwidth;              /* calc size of tick labels */
   ytl_wid = y_fldwid * charwidth;
   tl_height = charheight;

-  
+
+  //  rSGP.getViewport (xp_min, yp_min, xp_max, yp_max);

   /* calculate the extent of the plot frame */
   xp_min = o_xporigin;
   yp_min = o_yporigin;
   /* calculate the extent of the plot frame */
   xp_min = o_xporigin;
   yp_min = o_yporigin;

diff --git a/libctgraphics/sgp.cpp b/libctgraphics/sgp.cpp
index 52a2430ad942d98d0f63d2b934a671db10b5f6d0..c1004d699bc2a720b9c2e4e1ab80c16d8a53a0ca 100644 (file)
--- a/libctgraphics/sgp.cpp
+++ b/libctgraphics/sgp.cpp
@@ -7,7 +7,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: sgp.cpp,v 1.11 2000/08/27 20:32:55 kevin Exp $
+**  $Id: sgp.cpp,v 1.12 2000/08/31 08:38:58 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
 **
 **  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


@@ -53,9 +53,8 @@ int SGP::s_iRGBColorCount = sizeof(s_aRGBColor) / sizeof(class RGBColor);
 
 #ifdef HAVE_WXWINDOWS
 SGPDriver::SGPDriver (wxDC* pDC, int xsize = 640, int ysize = 480)
 
 #ifdef HAVE_WXWINDOWS
 SGPDriver::SGPDriver (wxDC* pDC, int xsize = 640, int ysize = 480)

-  : m_iPhysicalXSize(xsize), m_iPhysicalYSize(ysize), m_idDriver(0)
+  : m_iPhysicalXSize(xsize), m_iPhysicalYSize(ysize), m_idDriver(0), m_pDC(pDC)

 {
 {

-  m_pDC = pDC;

   m_idDriver |= SGPDRIVER_WXWINDOWS;
 }
 #endif
   m_idDriver |= SGPDRIVER_WXWINDOWS;
 }
 #endif


@@ -103,7 +102,6 @@ SGP::SGP (const SGPDriver& driver)
   setTextAngle (0.);
   setTextSize (1. / 25.);
   setColor (C_BLACK);
   setTextAngle (0.);
   setTextSize (1. / 25.);
   setColor (C_BLACK);

-

 }
 
 
 }
 
 


@@ -199,6 +197,15 @@ SGP::setViewport (double xmin, double ymin, double xmax, double ymax)
   viewNDC[3] = ymax;
 }
 
   viewNDC[3] = ymax;
 }
 

+void
+SGP::getViewport (double& xmin, double& ymin, double& xmax, double& ymax)
+{
+    xmin = xv_min;
+    ymin = yv_min;
+    xmax = xv_max;
+    ymax = yv_max;
+}
+

 
 // NAME
 //     frameViewport           draw box around viewport
 
 // NAME
 //     frameViewport           draw box around viewport


@@ -566,47 +573,46 @@ SGP::drawRect (double xmin, double ymin, double xmax, double ymax)
 void 
 SGP::drawCircle (const double r)
 {
 void 
 SGP::drawCircle (const double r)
 {

-       drawArc (r, 0.0, 7.0);
+       drawArc (r, 0.0, TWOPI);

 }
 
 }
 

-// =============================================================
-// draw arc around current center.  pass angles and radius     
+//==============================================================
+// draw arc around current center.  angles in radius   

 //==============================================================
 
 void 
 SGP::drawArc (const double r, double start, double stop)
 {
 //==============================================================
 
 void 
 SGP::drawArc (const double r, double start, double stop)
 {

-  if ((stop-start) > 2 * PI)
-    stop = start + 2 * PI;
-  if ((start-stop) > 2 * PI)
-    stop = start + 2 * PI;
-  while (start >= stop)
-    stop += 2*PI;
+  if (start > stop) {
+    double temp = start;
+    start = stop;
+    stop = temp;
+  }

 
   double x = r * cos ((double) start);
   double y = r * sin ((double) start);
   moveRel (x, y);          // move from center to start of arc 
 
 
   double x = r * cos ((double) start);
   double y = r * sin ((double) start);
   moveRel (x, y);          // move from center to start of arc 
 

-  double theta = 5 * PI / 180;
-  double c = cos(theta);
-  double s = sin(theta);
+  const double thetaIncrement = (5 * (TWOPI / 360)); 
+  double cosTheta = cos(thetaIncrement);
+  double sinTheta = sin(thetaIncrement);

 
   double angle, xp, yp;
 
   double angle, xp, yp;

-  for (angle = start; angle < stop - theta; angle += theta) {
-    xp = c * x - s * y;
-    yp = s * x + c * y;
-    lineRel (xp - x, yp - y);
+  for (angle = start; angle < stop - thetaIncrement; angle += thetaIncrement) {
+    xp = cosTheta * x - sinTheta * y; // translate point by thetaIncrement
+    yp = sinTheta * x + cosTheta * y;
+    lineAbs (xp, yp);

     x = xp; y = yp;
   }
 
     x = xp; y = yp;
   }
 

-  c = cos (stop - angle);
-  s = sin (stop - angle);
+  double c = cos (stop - angle);
+  double s = sin (stop - angle);

   xp = c * x - s * y;
   yp = s * x + c * y;
   xp = c * x - s * y;
   yp = s * x + c * y;

-  lineRel (xp - x, yp - y);
+  lineAbs (xp, yp);

 
 

-  x = r * cos ((double) stop);
-  y = r * sin ((double) stop);
+  x = r * cos (stop);
+  y = r * sin (stop);

   moveRel (-x, -y);            // move back to center of circle 
 }
 
   moveRel (-x, -y);            // move back to center of circle 
 }
 


@@ -783,3 +789,13 @@ const unsigned char SGP::MARKER_BITMAP[MARK_COUNT][5] =
     {'\076', '\042', '\042', '\042', '\076'},    // big open square 
     {'\010', '\024', '\042', '\024', '\010'},    // big open diamond 
 };
     {'\076', '\042', '\042', '\042', '\076'},    // big open square 
     {'\010', '\024', '\042', '\024', '\010'},    // big open diamond 
 };

+
+
+#if HAVE_WXWINDOWS
+void
+SGP::setDC (wxDC* pDC)
+{
+  if (m_driver.isWX())
+    m_driver.setDC(pDC);
+}
+#endif

diff --git a/libctsim/backprojectors.cpp b/libctsim/backprojectors.cpp
index 71aa8d4824722d2940d874b0b87666007121954f..0bcd4077cd527281d11769b4084c3d29e07cce1d 100644 (file)
--- a/libctsim/backprojectors.cpp
+++ b/libctsim/backprojectors.cpp
@@ -8,7 +8,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: backprojectors.cpp,v 1.12 2000/08/25 15:59:13 kevin Exp $
+**  $Id: backprojectors.cpp,v 1.13 2000/08/31 08:38:58 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
 **
 **  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


@@ -148,22 +148,27 @@ Backprojector::initBackprojector (const Projections& proj, ImageFile& im, const
     return false;
   }
 
     return false;
   }
 

-  if (m_idBackproject == BPROJ_TRIG)
-    m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectTrig (proj, im, m_idInterpolation, interpFactor));
-  else if (m_idBackproject == BPROJ_TABLE)
-    m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectTable (proj, im, m_idInterpolation, interpFactor));
-  else if (m_idBackproject == BPROJ_DIFF)
-    m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectDiff (proj, im, m_idInterpolation, interpFactor));
-  else if (m_idBackproject == BPROJ_DIFF2)
-    m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectDiff2 (proj, im, m_idInterpolation, interpFactor));
-  else if (m_idBackproject == BPROJ_IDIFF2)
-    m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectIntDiff2 (proj, im, m_idInterpolation, interpFactor));
-  else if (m_idBackproject == BPROJ_IDIFF3)
-    m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectIntDiff3 (proj, im, m_idInterpolation, interpFactor));
-  else {
-    m_fail = true;
-    m_failMessage = "Unable to select a backprojection method [Backprojector::initBackprojector]";
-    return false;
+  if (proj.geometry() == Scanner::GEOMETRY_EQUILINEAR)
+      m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectEquilinear(proj, im, m_idInterpolation, interpFactor));
+  else if (proj.geometry() == Scanner::GEOMETRY_EQUIANGULAR) 
+      m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectEquiangular(proj, im, m_idInterpolation, interpFactor));
+  else if (proj.geometry() == Scanner::GEOMETRY_PARALLEL) {
+      if (m_idBackproject == BPROJ_TRIG)
+         m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectTrig (proj, im, m_idInterpolation, interpFactor));
+      else if (m_idBackproject == BPROJ_TABLE)
+         m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectTable (proj, im, m_idInterpolation, interpFactor));
+      else if (m_idBackproject == BPROJ_DIFF)
+         m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectDiff (proj, im, m_idInterpolation, interpFactor));
+      else if (m_idBackproject == BPROJ_DIFF2)
+         m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectDiff2 (proj, im, m_idInterpolation, interpFactor));
+      else if (m_idBackproject == BPROJ_IDIFF2)
+         m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectIntDiff2 (proj, im, m_idInterpolation, interpFactor));
+      else if (m_idBackproject == BPROJ_IDIFF3)
+         m_pBackprojectImplem = static_cast<Backproject*>(new BackprojectIntDiff3 (proj, im, m_idInterpolation, interpFactor));
+  } else {
+      m_fail = true;
+      m_failMessage = "Unable to select a backprojection method [Backprojector::initBackprojector]";
+      return false;

   }
 
   return true;
   }
 
   return true;


@@ -256,7 +261,7 @@ Backproject::Backproject (const Projections& proj, ImageFile& im, const int inte
 {
   detInc = proj.detInc();
   nDet = proj.nDet();
 {
   detInc = proj.detInc();
   nDet = proj.nDet();

-  iDetCenter = (nDet - 1) / 2; // index refering to L=0 projection 
+  iDetCenter = nDet / 2;       // index refering to L=0 projection 

   rotInc = proj.rotInc();
 
   v = im.getArray();
   rotInc = proj.rotInc();
 
   v = im.getArray();


@@ -271,6 +276,8 @@ Backproject::Backproject (const Projections& proj, ImageFile& im, const int inte
 
   xInc = (xMax - xMin) / nx;   // size of cells
   yInc = (yMax - yMin) / ny;
 
   xInc = (xMax - xMin) / nx;   // size of cells
   yInc = (yMax - yMin) / ny;

+
+  m_dFocalLength = proj.focalLength();

 }
 
 Backproject::~Backproject ()
 }
 
 Backproject::~Backproject ()


@@ -312,12 +319,12 @@ void Backproject::errorIndexOutsideDetector (int ix, int iy, double theta, doubl
 void
 BackprojectTrig::BackprojectView (const double* const filteredProj, const double view_angle)
 {
 void
 BackprojectTrig::BackprojectView (const double* const filteredProj, const double view_angle)
 {

-  double theta = HALFPI + view_angle;   // Add PI/2 to get perpendicular angle to detector     
-  int ix, iy;
-  double x, y;                 // Rectang coords of center of pixel 
+  double theta = view_angle;

 
 

-  for (x = xMin + xInc / 2, ix = 0; ix < nx; x += xInc, ix++)
-    for (y = yMin + yInc / 2, iy = 0; iy < ny; y += yInc, iy++) {
+  double x = xMin + xInc / 2;  // Rectang coords of center of pixel 
+  for (int ix = 0; ix < nx; x += xInc, ix++) {
+    double y = yMin + yInc / 2;
+    for (int iy = 0; iy < ny; y += yInc, iy++) {

       double r = sqrt (x * x + y * y);   // distance of cell from center
       double phi = atan2 (y, x);         // angle of cell from center
       double L = r * cos (theta - phi);  // position on detector
       double r = sqrt (x * x + y * y);   // distance of cell from center
       double phi = atan2 (y, x);         // angle of cell from center
       double L = r * cos (theta - phi);  // position on detector


@@ -340,6 +347,7 @@ BackprojectTrig::BackprojectView (const double* const filteredProj, const double
            v[ix][iy] += rotInc * ((1-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1]);
       }
     }
            v[ix][iy] += rotInc * ((1-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1]);
       }
     }

+  }

 }  
 
 
 }  
 
 


@@ -374,7 +382,7 @@ BackprojectTable::~BackprojectTable ()
 void
 BackprojectTable::BackprojectView (const double* const filteredProj, const double view_angle)
 {
 void
 BackprojectTable::BackprojectView (const double* const filteredProj, const double view_angle)
 {

-  double theta = HALFPI + view_angle;  // add half PI to view angle to get perpendicular theta angle
+  double theta = view_angle;

 
   for (int ix = 0; ix < nx; ix++) {
     ImageFileColumn pImCol = v[ix];
 
   for (int ix = 0; ix < nx; ix++) {
     ImageFileColumn pImCol = v[ix];


@@ -431,9 +439,9 @@ BackprojectDiff::~BackprojectDiff()
 void
 BackprojectDiff::BackprojectView (const double* const filteredProj, const double view_angle)
 {
 void
 BackprojectDiff::BackprojectView (const double* const filteredProj, const double view_angle)
 {

-  double theta = - view_angle;  // add half PI to view angle to get perpendicular theta angle
-  double det_dx = xInc * sin (theta);
-  double det_dy = yInc * cos (theta);
+  double theta = view_angle;  // add half PI to view angle to get perpendicular theta angle
+  double det_dx = xInc * cos (theta);
+  double det_dy = yInc * sin (theta);

   double lColStart = start_r * cos (theta - start_phi);  // calculate L for first point in image
        
   for (int ix = 0; ix < nx; ix++, lColStart += det_dx) {
   double lColStart = start_r * cos (theta - start_phi);  // calculate L for first point in image
        
   for (int ix = 0; ix < nx; ix++, lColStart += det_dx) {


@@ -475,11 +483,11 @@ BackprojectDiff::BackprojectView (const double* const filteredProj, const double
 void
 BackprojectDiff2::BackprojectView (const double* const filteredProj, const double view_angle)
 {
 void
 BackprojectDiff2::BackprojectView (const double* const filteredProj, const double view_angle)
 {

-  double theta = - view_angle;  // add half PI to view angle to get perpendicular theta angle
+  double theta = view_angle;

 
   // Distance betw. detectors for an angle given in units of detectors 
 
   // Distance betw. detectors for an angle given in units of detectors 

-  double det_dx = xInc * sin (theta) / detInc;
-  double det_dy = yInc * cos (theta) / detInc;
+  double det_dx = xInc * cos (theta) / detInc;
+  double det_dy = yInc * sin (theta) / detInc;

 
   // calculate detPosition for first point in image (ix=0, iy=0) 
   double detPosColStart = start_r * cos (theta - start_phi) / detInc;
 
   // calculate detPosition for first point in image (ix=0, iy=0) 
   double detPosColStart = start_r * cos (theta - start_phi) / detInc;


@@ -524,14 +532,14 @@ BackprojectDiff2::BackprojectView (const double* const filteredProj, const doubl
 void
 BackprojectIntDiff2::BackprojectView (const double* const filteredProj, const double view_angle)
 {
 void
 BackprojectIntDiff2::BackprojectView (const double* const filteredProj, const double view_angle)
 {

-  double theta = - view_angle;  // add half PI to view angle to get perpendicular theta angle
+  double theta = view_angle;

 
   static const kint32 scale = 1 << 16;
   static const double dScale = scale;
   static const kint32 halfScale = scale / 2;
 
 
   static const kint32 scale = 1 << 16;
   static const double dScale = scale;
   static const kint32 halfScale = scale / 2;
 

-  const kint32 det_dx = nearest<kint32> (xInc * sin (theta) / detInc * scale);
-  const kint32 det_dy = nearest<kint32> (yInc * cos (theta) / detInc * scale);
+  const kint32 det_dx = nearest<kint32> (xInc * cos (theta) / detInc * scale);
+  const kint32 det_dy = nearest<kint32> (yInc * sin (theta) / detInc * scale);

 
   // calculate L for first point in image (0, 0) 
   kint32 detPosColStart = nearest<kint32> (start_r * cos (theta - start_phi) / detInc * scale);
 
   // calculate L for first point in image (0, 0) 
   kint32 detPosColStart = nearest<kint32> (start_r * cos (theta - start_phi) / detInc * scale);


@@ -576,15 +584,15 @@ BackprojectIntDiff2::BackprojectView (const double* const filteredProj, const do
 void
 BackprojectIntDiff3::BackprojectView (const double* const filteredProj, const double view_angle)
 {
 void
 BackprojectIntDiff3::BackprojectView (const double* const filteredProj, const double view_angle)
 {

-  double theta = - view_angle;  // add half PI to view angle to get perpendicular theta angle
+  double theta = view_angle;  // add half PI to view angle to get perpendicular theta angle

   static const int scaleShift = 16;
   static const kint32 scale = (1 << scaleShift);
   static const kint32 scaleBitmask = scale - 1;
   static const kint32 halfScale = scale / 2;
   static const double dInvScale = 1. / scale;
 
   static const int scaleShift = 16;
   static const kint32 scale = (1 << scaleShift);
   static const kint32 scaleBitmask = scale - 1;
   static const kint32 halfScale = scale / 2;
   static const double dInvScale = 1. / scale;
 

-  const kint32 det_dx = nearest<kint32> (xInc * sin (theta) / detInc * scale);
-  const kint32 det_dy = nearest<kint32> (yInc * cos (theta) / detInc * scale);
+  const kint32 det_dx = nearest<kint32> (xInc * cos (theta) / detInc * scale);
+  const kint32 det_dy = nearest<kint32> (yInc * sin (theta) / detInc * scale);

 
   // calculate L for first point in image (0, 0) 
   kint32 detPosColStart = nearest<kint32> ((start_r * cos (theta - start_phi) / detInc + iDetCenter) * scale);
 
   // calculate L for first point in image (0, 0) 
   kint32 detPosColStart = nearest<kint32> ((start_r * cos (theta - start_phi) / detInc + iDetCenter) * scale);


@@ -624,3 +632,82 @@ BackprojectIntDiff3::BackprojectView (const double* const filteredProj, const do
     } //end linear
   } // end for ix
 }
     } //end linear
   } // end for ix
 }

+
+
+void
+BackprojectEquiangular::BackprojectView (const double* const filteredProj, const double view_angle)
+{
+  double beta = view_angle;
+
+  for (int ix = 0; ix < nx; ix++) {
+    ImageFileColumn pImCol = v[ix];
+
+    for (int iy = 0; iy < ny; iy++) { 
+      double dAngleDiff = beta - phi[ix][iy];
+      double rcos_t = r[ix][iy] * cos (dAngleDiff);
+      double rsin_t = r[ix][iy] * sin (dAngleDiff);
+      double dFLPlusSin = m_dFocalLength + rsin_t;
+      double gamma =  atan (rcos_t / dFLPlusSin);
+      double dL2 = dFLPlusSin * dFLPlusSin + (rcos_t * rcos_t);
+
+      if (interpType == Backprojector::INTERP_NEAREST) {
+       int iDetPos =iDetCenter + nearest<int>(gamma / detInc); // calc index in the filtered raysum vector 
+
+       if (iDetPos < 0 || iDetPos >= nDet) { // check for impossible: index outside of raysum pos  
+         ; //    errorIndexOutsideDetector (ix, iy, beta, r[ix][iy], phi[ix][iy], gamma, iDetPos);
+       } else
+         pImCol[iy] += filteredProj[iDetPos] / dL2;
+      } else if (interpType == Backprojector::INTERP_LINEAR) {
+       double dPos = gamma / detInc;           // position along detector 
+       double dPosFloor = floor (dPos);
+       int iDetPos = iDetCenter + static_cast<int>(dPosFloor);
+       double frac = dPos - dPosFloor; // fraction distance from det 
+       if (iDetPos < 0 || iDetPos >= nDet - 1) {
+         ; //      errorIndexOutsideDetector (ix, iy, beta, r[ix][iy], phi[ix][iy], gamma, iDetPos);
+       } else
+         pImCol[iy] += (((1-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1])) / dL2;
+      }
+    }  // end for y 
+  }    // end for x 
+}
+
+void
+BackprojectEquilinear::BackprojectView (const double* const filteredProj, const double view_angle)
+{
+  double beta = view_angle;
+
+  for (int ix = 0; ix < nx; ix++) {
+    ImageFileColumn pImCol = v[ix];
+
+    for (int iy = 0; iy < ny; iy++) {
+      double dAngleDiff = beta - phi[ix][iy];
+      double rcos_t = r[ix][iy] * cos (dAngleDiff);
+      double rsin_t = r[ix][iy] * sin (dAngleDiff);
+
+      double dU = (m_dFocalLength + rsin_t) / m_dFocalLength;
+      double dDetPos =  rcos_t / dU;
+      // double to scale for imaginary detector that passes through origin
+      // of phantom, see Kak-Slaney Figure 3.22
+      dDetPos *= 2; 
+
+      if (interpType == Backprojector::INTERP_NEAREST) {
+       int iDetPos = iDetCenter + nearest<int>(dDetPos / detInc);      // calc index in the filtered raysum vector 
+
+       if (iDetPos < 0 || iDetPos >= nDet)     // check for impossible: index outside of raysum pos 
+         ; ///   errorIndexOutsideDetector (ix, iy, beta, r[ix][iy], phi[ix][iy], dDetPos, iDetPos);
+       else
+         pImCol[iy] += (filteredProj[iDetPos] / (dU * dU));
+      } else if (interpType == Backprojector::INTERP_LINEAR) {
+       double dPos = dDetPos / detInc;         // position along detector 
+       double dPosFloor = floor (dPos);
+       int iDetPos = iDetCenter + static_cast<int>(dPosFloor);
+       double frac = dPos - dPosFloor; // fraction distance from det 
+       if (iDetPos < 0 || iDetPos >= nDet - 1)
+         ; //      errorIndexOutsideDetector (ix, iy, beta, r[ix][iy], phi[ix][iy], dDetPos, iDetPos);
+       else
+         pImCol[iy] += (((1-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1])) / (dU * dU);
+      }
+    }  // end for y 
+  }    // end for x 
+}
+

diff --git a/libctsim/procsignal.cpp b/libctsim/procsignal.cpp
index e38fae239f0095e9b499425599526cfe0460de54..5c10c9e9ff2a44b6b9ba602c5352c52bf6c1094e 100644 (file)
--- a/libctsim/procsignal.cpp
+++ b/libctsim/procsignal.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: procsignal.cpp,v 1.4 2000/08/27 20:32:55 kevin Exp $
+**  $Id: procsignal.cpp,v 1.5 2000/08/31 08:38:58 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
 **
 **  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


@@ -77,7 +77,7 @@ const int ProcessSignal::s_iFilterGenerationCount = sizeof(s_aszFilterGeneration
 // CLASS IDENTIFICATION
 //   ProcessSignal
 //
 // CLASS IDENTIFICATION
 //   ProcessSignal
 //

-ProcessSignal::ProcessSignal (const char* szFilterName, const char* szFilterMethodName, double dBandwidth, double dSignalIncrement, int nSignalPoints, double dFilterParam, const char* szDomainName, const char* szFilterGenerationName, int iZeropad, int iPreinterpolationFactor, int iTraceLevel, int iGeometry)
+ProcessSignal::ProcessSignal (const char* szFilterName, const char* szFilterMethodName, double dBandwidth, double dSignalIncrement, int nSignalPoints, double dFilterParam, const char* szDomainName, const char* szFilterGenerationName, int iZeropad, int iPreinterpolationFactor, int iTraceLevel, int iGeometry, double dFocalLength)

     : m_adFourierCosTable(NULL), m_adFourierSinTable(NULL), m_adFilter(NULL), m_fail(false)
 {
   m_idFilterMethod = convertFilterMethodNameToID (szFilterMethodName);
     : m_adFourierCosTable(NULL), m_adFourierSinTable(NULL), m_adFilter(NULL), m_fail(false)
 {
   m_idFilterMethod = convertFilterMethodNameToID (szFilterMethodName);


@@ -109,18 +109,19 @@ ProcessSignal::ProcessSignal (const char* szFilterName, const char* szFilterMeth
     return;
   }
 
     return;
   }
 

-  init (m_idFilter, m_idFilterMethod, dBandwidth, dSignalIncrement, nSignalPoints, dFilterParam, m_idDomain, m_idFilterGeneration, iZeropad, iPreinterpolationFactor, iTraceLevel, iGeometry);
+  init (m_idFilter, m_idFilterMethod, dBandwidth, dSignalIncrement, nSignalPoints, dFilterParam, m_idDomain, m_idFilterGeneration, iZeropad, iPreinterpolationFactor, iTraceLevel, iGeometry, dFocalLength);

 }
 
 
 void
 }
 
 
 void

-ProcessSignal::init (const int idFilter, const int idFilterMethod, double dBandwidth, double dSignalIncrement, int nSignalPoints, double dFilterParam, const int idDomain, const int idFilterGeneration, const int iZeropad, const int iPreinterpolationFactor, int iTraceLevel, int iGeometry)
+ProcessSignal::init (const int idFilter, const int idFilterMethod, double dBandwidth, double dSignalIncrement, int nSignalPoints, double dFilterParam, const int idDomain, const int idFilterGeneration, const int iZeropad, const int iPreinterpolationFactor, int iTraceLevel, int iGeometry, double dFocalLength)

 {
   m_idFilter = idFilter;
   m_idDomain = idDomain;
   m_idFilterMethod = idFilterMethod;
   m_idFilterGeneration = idFilterGeneration;
   m_idGeometry = iGeometry;
 {
   m_idFilter = idFilter;
   m_idDomain = idDomain;
   m_idFilterMethod = idFilterMethod;
   m_idFilterGeneration = idFilterGeneration;
   m_idGeometry = iGeometry;

+  m_dFocalLength = dFocalLength;

 
   if (m_idFilter == SignalFilter::FILTER_INVALID || m_idDomain == SignalFilter::DOMAIN_INVALID || m_idFilterMethod == FILTER_METHOD_INVALID || m_idFilterGeneration == FILTER_GENERATION_INVALID) {
     m_fail = true;
 
   if (m_idFilter == SignalFilter::FILTER_INVALID || m_idDomain == SignalFilter::DOMAIN_INVALID || m_idFilterMethod == FILTER_METHOD_INVALID || m_idFilterGeneration == FILTER_GENERATION_INVALID) {
     m_fail = true;


@@ -136,6 +137,13 @@ ProcessSignal::init (const int idFilter, const int idFilterMethod, double dBandw
   m_iZeropad = iZeropad;
   m_iPreinterpolationFactor = iPreinterpolationFactor;
 
   m_iZeropad = iZeropad;
   m_iPreinterpolationFactor = iPreinterpolationFactor;
 

+  // scale signalInc/BW to signalInc/2 to adjust for imaginary detector
+  // through origin of phantom, see Kak-Slaney Fig 3.22, for Collinear
+  if (m_idGeometry == Scanner::GEOMETRY_EQUILINEAR) {
+    m_dSignalInc /= 2;
+    m_dBandwidth *= 2;
+  }
+

   if (m_idFilterMethod == FILTER_METHOD_FFT) {
 #if HAVE_FFTW
     m_idFilterMethod = FILTER_METHOD_RFFTW;
   if (m_idFilterMethod == FILTER_METHOD_FFT) {
 #if HAVE_FFTW
     m_idFilterMethod = FILTER_METHOD_RFFTW;


@@ -222,10 +230,24 @@ ProcessSignal::init (const int idFilter, const int idFilterMethod, double dBandw
            m_adFilter[i] /= m_dSignalInc;
        }
     }
            m_adFilter[i] /= m_dSignalInc;
        }
     }

-  } 
+    if (m_idGeometry == Scanner::GEOMETRY_EQUILINEAR) {
+       for (int i = 0; i < m_nFilterPoints; i++)
+           m_adFilter[i] *= 0.5;
+    } else if (m_idGeometry == Scanner::GEOMETRY_EQUIANGULAR) {
+       for (int i = 0; i < m_nFilterPoints; i++) {
+         int iDetFromZero = i - ((m_nFilterPoints - 1) / 2);
+         double sinScale = sin (iDetFromZero * m_dSignalInc);
+         if (fabs(sinScale) < 1E-7)
+             sinScale = 1;
+         else
+             sinScale = (iDetFromZero * m_dSignalInc) / sinScale;
+         double dScale = 0.5 * sinScale * sinScale;
+         m_adFilter[i] *= dScale;
+       }
+    } // if (geometry)
+  } // if (spatial filtering)

 
 

-  // Frequency-based filtering
-  else if (m_bFrequencyFiltering) {
+  else if (m_bFrequencyFiltering) {  // Frequency-based filtering

 
     if (m_idFilterGeneration == FILTER_GENERATION_DIRECT) {
       // calculate number of filter points with zeropadding
 
     if (m_idFilterGeneration == FILTER_GENERATION_DIRECT) {
       // calculate number of filter points with zeropadding


@@ -256,6 +278,22 @@ ProcessSignal::init (const int idFilter, const int idFilterMethod, double dBandw
       SignalFilter filter (m_idFilter, m_dFilterMin, m_dFilterMax, m_nFilterPoints, m_dBandwidth, m_dFilterParam, SignalFilter::DOMAIN_FREQUENCY);
       m_adFilter = new double [m_nFilterPoints];
       filter.copyFilterData (m_adFilter, 0, m_nFilterPoints);
       SignalFilter filter (m_idFilter, m_dFilterMin, m_dFilterMax, m_nFilterPoints, m_dBandwidth, m_dFilterParam, SignalFilter::DOMAIN_FREQUENCY);
       m_adFilter = new double [m_nFilterPoints];
       filter.copyFilterData (m_adFilter, 0, m_nFilterPoints);

+
+      if (m_idGeometry == Scanner::GEOMETRY_EQUILINEAR) {
+       for (int i = 0; i < m_nFilterPoints; i++)
+         m_adFilter[i] *= 0.5;
+      } else if (m_idGeometry == Scanner::GEOMETRY_EQUIANGULAR) {
+       for (int i = 0; i < m_nFilterPoints; i++) {
+         int iDetFromZero = i - ((m_nFilterPoints - 1) / 2);
+         double sinScale = sin (iDetFromZero * m_dSignalInc);
+         if (fabs(sinScale) < 1E-7)
+           sinScale = 1;
+         else
+           sinScale = (iDetFromZero * m_dSignalInc) / sinScale;
+         double dScale = 0.5 * sinScale * sinScale;
+         //      m_adFilter[i] *= dScale;
+       }
+      }

       if (m_traceLevel >= Trace::TRACE_PLOT) {
        SGPDriver sgpDriver ("Frequency Filter: Natural Order");
        SGP sgp (sgpDriver);
       if (m_traceLevel >= Trace::TRACE_PLOT) {
        SGPDriver sgpDriver ("Frequency Filter: Natural Order");
        SGP sgp (sgpDriver);


@@ -311,6 +349,21 @@ ProcessSignal::init (const int idFilter, const int idFilterMethod, double dBandw
        cio_put_str ("Press any key to continue");
        cio_kb_getc ();
       }
        cio_put_str ("Press any key to continue");
        cio_kb_getc ();
       }

+      if (m_idGeometry == Scanner::GEOMETRY_EQUILINEAR) {
+       for (int i = 0; i < m_nFilterPoints; i++)
+         adSpatialFilter[i] *= 0.5;
+      } else if (m_idGeometry == Scanner::GEOMETRY_EQUIANGULAR) {
+       for (int i = 0; i < m_nFilterPoints; i++) {
+         int iDetFromZero = i - ((m_nFilterPoints - 1) / 2);
+         double sinScale = sin (iDetFromZero * m_dSignalInc);
+         if (fabs(sinScale) < 1E-7)
+           sinScale = 1;
+         else
+           sinScale = (iDetFromZero * m_dSignalInc) / sinScale;
+         double dScale = 0.5 * sinScale * sinScale;
+         adSpatialFilter[i] *= dScale;
+       }
+      }

       for (int i = nSpatialPoints; i < m_nFilterPoints; i++)
        adSpatialFilter[i] = 0;
 
       for (int i = nSpatialPoints; i < m_nFilterPoints; i++)
        adSpatialFilter[i] = 0;
 


@@ -470,11 +523,26 @@ ProcessSignal::convertFilterGenerationIDToTitle (const int fgID)
 }
 
 void
 }
 
 void

-ProcessSignal::filterSignal (const float input[], double output[], int iView) const
+ProcessSignal::filterSignal (const float constInput[], double output[]) const

 {
 {

+  double input [m_nSignalPoints];
+  for (int i = 0; i < m_nSignalPoints; i++)
+    input[i] = constInput[i];
+
+  if (m_idGeometry == Scanner::GEOMETRY_EQUILINEAR) {
+    for (int i = 0; i < m_nSignalPoints; i++) {
+      int iDetFromCenter = i - (m_nSignalPoints / 2);
+      input[i] *= m_dFocalLength / sqrt (m_dFocalLength * m_dFocalLength + iDetFromCenter * iDetFromCenter * m_dSignalInc * m_dSignalInc);
+    }
+  } else if (m_idGeometry == Scanner::GEOMETRY_EQUIANGULAR) {
+    for (int i = 0; i < m_nSignalPoints; i++) {
+      int iDetFromCenter = i - (m_nSignalPoints / 2);
+      input[i] *= m_dFocalLength * cos (iDetFromCenter * m_dSignalInc);
+    }
+  }

   if (m_idFilterMethod == FILTER_METHOD_CONVOLUTION) {
   if (m_idFilterMethod == FILTER_METHOD_CONVOLUTION) {

-    for (int i = 0; i < m_nSignalPoints; i++)
-      output[i] = convolve (input, m_dSignalInc, i, m_nSignalPoints);
+      for (int i = 0; i < m_nSignalPoints; i++)
+       output[i] = convolve (input, m_dSignalInc, i, m_nSignalPoints);

   } else if (m_idFilterMethod == FILTER_METHOD_FOURIER) {
     double inputSignal[m_nFilterPoints];
     for (int i = 0; i < m_nSignalPoints; i++)
   } else if (m_idFilterMethod == FILTER_METHOD_FOURIER) {
     double inputSignal[m_nFilterPoints];
     for (int i = 0; i < m_nSignalPoints; i++)

diff --git a/libctsim/projections.cpp b/libctsim/projections.cpp
index da298f01622003e1dde47bbc1da6757a51a2cab9..be8e82941f9a47b30fe1946905fa9fba96bc7f0c 100644 (file)
--- a/libctsim/projections.cpp
+++ b/libctsim/projections.cpp
@@ -8,7 +8,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: projections.cpp,v 1.23 2000/08/27 20:32:55 kevin Exp $
+**  $Id: projections.cpp,v 1.24 2000/08/31 08:38:58 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
 **
 **  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


@@ -80,10 +80,17 @@ Projections::initFromScanner (const Scanner& scanner)
   m_phmLen = scanner.phmLen();
   m_rotInc = scanner.rotInc();
   m_detInc = scanner.detInc();
   m_phmLen = scanner.phmLen();
   m_rotInc = scanner.rotInc();
   m_detInc = scanner.detInc();

+  m_geometry = scanner.geometry();

   m_focalLength = scanner.focalLength();
   m_fieldOfView = scanner.fieldOfView();
   m_rotStart = 0;
   m_focalLength = scanner.focalLength();
   m_fieldOfView = scanner.fieldOfView();
   m_rotStart = 0;

-  m_detStart =  -(scanner.detLen() / 2) + (scanner.detInc() / 2);
+  m_detStart =  -(scanner.detLen() / 2);
+#if 0
+  if (m_geometry == Scanner::GEOMETRY_EQUILINEAR) {
+    m_detInc /= 2;
+    cout << "Kludge: detInc /= 2 in Projections::initFromScanner" << endl;
+  }
+#endif

 }
 
 void
 }
 
 void


@@ -286,6 +293,13 @@ Projections::headerRead (fnetorderstream& fs)
   return true;
 }
 
   return true;
 }
 

+bool
+Projections::read (const string& filename)
+{
+  return read (filename.c_str());
+}
+
+

 bool
 Projections::read (const char* filename)
 {
 bool
 Projections::read (const char* filename)
 {


@@ -311,6 +325,12 @@ Projections::read (const char* filename)
 }
 
 
 }
 
 

+bool
+Projections::write (const string& filename)
+{
+  return write (filename.c_str());
+}
+

 bool
 Projections::write (const char* filename)
 {
 bool
 Projections::write (const char* filename)
 {


@@ -442,14 +462,16 @@ Projections::detarrayWrite (fnetorderstream& fs, const DetectorArray& darray, co
 void
 Projections::printProjectionData (void)
 {
 void
 Projections::printProjectionData (void)
 {

-  printf("Projections Print\n\n");
+  printf("Projections Data\n\n");

   printf("Description: %s\n", m_remark.c_str());
   printf("Description: %s\n", m_remark.c_str());

-  printf("nView = %8d          nDet = %8d\n", m_nView, m_nDet);
-  printf("rotStart = %8.4f     rotInc = %8.4f\n", m_rotStart, m_rotInc);
-  printf("detStart = %8.4f     detInc = %8.4f\n", m_detStart, m_detInc);
+  printf("Geometry: %s\n", Scanner::convertGeometryIDToName (m_geometry));
+  printf("nView       = %8d           nDet = %8d\n", m_nView, m_nDet);

   printf("focalLength = %8.4f  fieldOfView = %8.4f\n", m_focalLength, m_fieldOfView);
   printf("focalLength = %8.4f  fieldOfView = %8.4f\n", m_focalLength, m_fieldOfView);

+  printf("rotStart    = %8.4f       rotInc = %8.4f\n", m_rotStart, m_rotInc);
+  printf("detStart    = %8.4f       detInc = %8.4f\n", m_detStart, m_detInc);

   if (m_projData != NULL) {
       for (int ir = 0; ir < m_nView; ir++) {
   if (m_projData != NULL) {
       for (int ir = 0; ir < m_nView; ir++) {

+         printf("View %d: angle %f\n", ir, m_projData[ir]->viewAngle());

          DetectorValue* detval = m_projData[ir]->detValues();
          for (int id = 0; id < m_projData[ir]->nDet(); id++)
              printf("%8.4f  ", detval[id]);
          DetectorValue* detval = m_projData[ir]->detValues();
          for (int id = 0; id < m_projData[ir]->nDet(); id++)
              printf("%8.4f  ", detval[id]);


@@ -459,16 +481,19 @@ Projections::printProjectionData (void)
 }
 
 void 
 }
 
 void 

-Projections::printScanInfo (void) const
+Projections::printScanInfo (ostringstream& os) const

 {
 {

-  printf ("Number of detectors: %d\n", m_nDet);
-  printf ("    Number of views: %d\n", m_nView);
-  printf ("             Remark: %s\n", m_remark.c_str());
-  printf ("             phmLen: %f\n", m_phmLen);
-  printf ("           detStart: %f\n", m_detStart);
-  printf ("             detInc: %f\n", m_detInc);
-  printf ("           rotStart: %f\n", m_rotStart);
-  printf ("             rotInc: %f\n", m_rotInc);
+  os << "Number of detectors: " << m_nDet << "\n";
+  os << "    Number of views: " << m_nView<< "\n";
+  os << "             Remark: " << m_remark.c_str()<< "\n";
+  os << "           Geometry: " << Scanner::convertGeometryIDToName (m_geometry)<< "\n";
+  os << "       Focal Length: " << m_focalLength<< "\n";
+  os << "      Field Of View: " << m_fieldOfView<< "\n";
+  os << "             phmLen: " << m_phmLen<< "\n";
+  os << "           detStart: " << m_detStart<< "\n";
+  os << "             detInc: " << m_detInc<< "\n";
+  os << "           rotStart: " << m_rotStart<< "\n";
+  os << "             rotInc: " << m_rotInc<< "\n";

 }
 
 
 }
 
 


@@ -512,7 +537,7 @@ Projections::reconstruct (ImageFile& im, const char* const filterName, double fi
 #endif
 
   double filterBW = 1. / detInc;
 #endif
 
   double filterBW = 1. / detInc;

-  ProcessSignal processSignal (filterName, filterMethodName, filterBW, m_detInc, m_nDet, filt_param, "spatial", filterGenerationName, zeropad, interpFactor, trace, m_geometry);
+  ProcessSignal processSignal (filterName, filterMethodName, filterBW, m_detInc, m_nDet, filt_param, "spatial", filterGenerationName, zeropad, interpFactor, trace, m_geometry, m_focalLength);

 
   if (processSignal.fail()) {
       sys_error (ERR_SEVERE, "%s [Projections::reconstruct]", processSignal.failMessage().c_str());
 
   if (processSignal.fail()) {
       sys_error (ERR_SEVERE, "%s [Projections::reconstruct]", processSignal.failMessage().c_str());


@@ -558,12 +583,12 @@ Projections::reconstruct (ImageFile& im, const char* const filterName, double fi
 
   for (int iview = 0; iview < m_nView; iview++)  {
     if (trace >= Trace::TRACE_CONSOLE) 
 
   for (int iview = 0; iview < m_nView; iview++)  {
     if (trace >= Trace::TRACE_CONSOLE) 

-      printf ("Reconstructing view %d (last = %d)\n",  iview, m_nView - 1);
+      cout <<"Reconstructing view " << iview << "(last = " << m_nView - 1 << ")\n";

       
     const DetectorArray& darray = getDetectorArray (iview);
     const DetectorValue* detval = darray.detValues();
 
       
     const DetectorArray& darray = getDetectorArray (iview);
     const DetectorValue* detval = darray.detValues();
 

-    processSignal.filterSignal (detval, filteredProj, iview);
+    processSignal.filterSignal (detval, filteredProj);

 
 #ifdef HAVE_BSPLINE_INTERP
     if (interp_type == I_BSPLINE) 
 
 #ifdef HAVE_BSPLINE_INTERP
     if (interp_type == I_BSPLINE) 

diff --git a/libctsim/scanner.cpp b/libctsim/scanner.cpp
index c45d4fd3234c1aa71561728d3bfb952b0840eeba..ef72c5b4cc99a7672a9d3ba2a25eac03a8d26edf 100644 (file)
--- a/libctsim/scanner.cpp
+++ b/libctsim/scanner.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: scanner.cpp,v 1.12 2000/08/27 20:32:55 kevin Exp $
+**  $Id: scanner.cpp,v 1.13 2000/08/31 08:38:58 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
 **
 **  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


@@ -30,21 +30,21 @@
 
 const int Scanner::GEOMETRY_INVALID = -1;
 const int Scanner::GEOMETRY_PARALLEL = 0;
 
 const int Scanner::GEOMETRY_INVALID = -1;
 const int Scanner::GEOMETRY_PARALLEL = 0;

-const int Scanner::GEOMETRY_EQUILINEAR = 1;
-const int Scanner::GEOMETRY_EQUIANGULAR = 2;
+const int Scanner::GEOMETRY_EQUIANGULAR = 1;
+const int Scanner::GEOMETRY_EQUILINEAR = 2;

 
 const char* Scanner::s_aszGeometryName[] = 
 {
   {"parallel"},
 
 const char* Scanner::s_aszGeometryName[] = 
 {
   {"parallel"},

-  {"equilinear"},

   {"equiangular"},
   {"equiangular"},

+  {"equilinear"},

 };
 
 const char* Scanner::s_aszGeometryTitle[] = 
 {
   {"Parallel"},
 };
 
 const char* Scanner::s_aszGeometryTitle[] = 
 {
   {"Parallel"},

-  {"Equilinear"},

   {"Equiangular"},
   {"Equiangular"},

+  {"Equilinear"},

 };
 
 const int Scanner::s_iGeometryCount = sizeof(s_aszGeometryName) / sizeof(const char*);
 };
 
 const int Scanner::s_iGeometryCount = sizeof(s_aszGeometryName) / sizeof(const char*);


@@ -101,8 +101,6 @@ Scanner::Scanner (const Phantom& phm, const char* const geometryName, int nDet,
     m_nSample = 1;
   if (nDet < 1)
     nDet = 1;
     m_nSample = 1;
   if (nDet < 1)
     nDet = 1;

-  //  if ((nDet % 2) == 0)
-  //    ++nDet;                // ensure odd number of detectors

 
   m_nDet     = nDet;
   m_nView    = nView;
 
   m_nDet     = nDet;
   m_nView    = nView;


@@ -114,23 +112,26 @@ Scanner::Scanner (const Phantom& phm, const char* const geometryName, int nDet,
 
   m_dXCenter = phm.xmin() + (phm.xmax() - phm.xmin()) / 2;
   m_dYCenter = phm.ymin() + (phm.ymax() - phm.ymin()) / 2;
 
   m_dXCenter = phm.xmin() + (phm.xmax() - phm.xmin()) / 2;
   m_dYCenter = phm.ymin() + (phm.ymax() - phm.ymin()) / 2;

+  m_rotLen  = rot_anglen;
+  m_rotInc  = m_rotLen / m_nView;

   if (m_idGeometry == GEOMETRY_PARALLEL) {
     m_detLen   = m_dFieldOfView;
     m_detInc  = m_detLen / m_nDet;
   if (m_idGeometry == GEOMETRY_PARALLEL) {
     m_detLen   = m_dFieldOfView;
     m_detInc  = m_detLen / m_nDet;

-    m_rotLen  = rot_anglen;
-    m_rotInc  = m_rotLen / m_nView;

     
     double dHalfDetLen = m_detLen / 2;
     m_initPos.xs1 = m_dXCenter - m_dFocalLength;
     
     double dHalfDetLen = m_detLen / 2;
     m_initPos.xs1 = m_dXCenter - m_dFocalLength;

-    m_initPos.ys1 = m_dYCenter - dHalfDetLen;
-    m_initPos.xs2 = m_dXCenter - m_dFocalLength;
+    m_initPos.ys1 = m_dYCenter + dHalfDetLen;
+    m_initPos.xs2 = m_dXCenter + m_dFocalLength;

     m_initPos.ys2 = m_dYCenter + dHalfDetLen;
     m_initPos.ys2 = m_dYCenter + dHalfDetLen;

-    m_initPos.xd1 = m_dXCenter + m_dFocalLength;
+    m_initPos.xd1 = m_dXCenter - m_dFocalLength;

     m_initPos.yd1 = m_dYCenter - dHalfDetLen;
     m_initPos.xd2 = m_dXCenter + m_dFocalLength;
     m_initPos.yd1 = m_dYCenter - dHalfDetLen;
     m_initPos.xd2 = m_dXCenter + m_dFocalLength;

-    m_initPos.yd2 = m_dYCenter + dHalfDetLen;
+    m_initPos.yd2 = m_dYCenter - dHalfDetLen;

     m_initPos.angle = 0.0;
   } else if (m_idGeometry == GEOMETRY_EQUILINEAR) {
     m_initPos.angle = 0.0;
   } else if (m_idGeometry == GEOMETRY_EQUILINEAR) {

+#if 0
+    double dAngle = (m_dFieldOfView / 2) / cos (asin (m_dFieldOfView / 2 / m_dFocalLength));
+#else

     double dHalfSquare = m_dFieldOfView / SQRT2 / 2;
     double dFocalPastPhm = m_dFocalLength - dHalfSquare;
     if (dFocalPastPhm <= 0.) {
     double dHalfSquare = m_dFieldOfView / SQRT2 / 2;
     double dFocalPastPhm = m_dFocalLength - dHalfSquare;
     if (dFocalPastPhm <= 0.) {


@@ -139,23 +140,26 @@ Scanner::Scanner (const Phantom& phm, const char* const geometryName, int nDet,
       return;
     }
     double dAngle = atan( dHalfSquare / dFocalPastPhm );
       return;
     }
     double dAngle = atan( dHalfSquare / dFocalPastPhm );

+#endif

     double dHalfDetLen = 2 * m_dFocalLength * tan (dAngle);
 
     m_detLen = dHalfDetLen * 2;
     m_detInc  = m_detLen / m_nDet;
     double dHalfDetLen = 2 * m_dFocalLength * tan (dAngle);
 
     m_detLen = dHalfDetLen * 2;
     m_detInc  = m_detLen / m_nDet;

-    m_rotLen  = rot_anglen;
-    m_rotInc  = m_rotLen / m_nView;

     
     

-    m_initPos.xs1 = m_dXCenter - m_dFocalLength;
-    m_initPos.ys1 = m_dYCenter;
-    m_initPos.xs2 = m_dXCenter - m_dFocalLength;
-    m_initPos.ys2 = m_dYCenter;
-    m_initPos.xd1 = m_dXCenter + m_dFocalLength;
-    m_initPos.yd1 = m_dYCenter - dHalfDetLen;
-    m_initPos.xd2 = m_dXCenter + m_dFocalLength;
-    m_initPos.yd2 = m_dYCenter + dHalfDetLen;
+    m_initPos.angle = 0.0;
+    m_initPos.xs1 = m_dXCenter;
+    m_initPos.ys1 = m_dYCenter + m_dFocalLength;
+    m_initPos.xs2 = m_dXCenter;
+    m_initPos.ys2 = m_dYCenter + m_dFocalLength;
+    m_initPos.xd1 = m_dXCenter - dHalfDetLen;
+    m_initPos.yd1 = m_dYCenter - m_dFocalLength;
+    m_initPos.xd2 = m_dXCenter + dHalfDetLen;
+    m_initPos.yd2 = m_dYCenter - m_dFocalLength;

     m_initPos.angle = 0.0;
   } else if (m_idGeometry == GEOMETRY_EQUIANGULAR) {
     m_initPos.angle = 0.0;
   } else if (m_idGeometry == GEOMETRY_EQUIANGULAR) {

+#if 0
+    double dAngle = atan ((m_dFieldOfView / 2) / m_dFocalLength);
+#else

     double dHalfSquare = m_dFieldOfView / SQRT2 / 2;
     double dFocalPastPhm = m_dFocalLength - dHalfSquare;
     if (dFocalPastPhm <= 0.) {
     double dHalfSquare = m_dFieldOfView / SQRT2 / 2;
     double dFocalPastPhm = m_dFocalLength - dHalfSquare;
     if (dFocalPastPhm <= 0.) {


@@ -163,17 +167,20 @@ Scanner::Scanner (const Phantom& phm, const char* const geometryName, int nDet,
       m_failMessage = "Focal Point inside of phantom";
       return;
     }
       m_failMessage = "Focal Point inside of phantom";
       return;
     }

-    double dAngle = 2 * atan( dHalfSquare / dFocalPastPhm );
+    double dAngle =  atan ( dHalfSquare / dFocalPastPhm ); 
+#endif

     m_detLen = 2 * dAngle;
     m_detInc = m_detLen / m_nDet;
     m_detLen = 2 * dAngle;
     m_detInc = m_detLen / m_nDet;

-    m_rotLen = rot_anglen;
-    m_rotInc = m_rotLen / m_nView;

 
 

-    m_initPos.xs1 = m_dXCenter - m_dFocalLength;
-    m_initPos.ys1 = m_dYCenter;
-    m_initPos.xs2 = m_dXCenter - m_dFocalLength;
-    m_initPos.ys2 = m_dYCenter;
-    m_initPos.angle = -dAngle;
+    m_dAngularDetIncrement = m_detInc * 2; // Angular Position 2x gamma angle
+    m_dAngularDetLen = m_detLen * 2;
+    m_initPos.dAngularDet = -m_dAngularDetLen / 2;
+
+    m_initPos.angle = 0;
+    m_initPos.xs1 = m_dXCenter;
+    m_initPos.ys1 = m_dYCenter + m_dFocalLength;;
+    m_initPos.xs2 = m_dXCenter;
+    m_initPos.ys2 = m_dYCenter + m_dFocalLength;

   }
 
   // Calculate incrementatal rotation matrix 
   }
 
   // Calculate incrementatal rotation matrix 


@@ -260,10 +267,8 @@ Scanner::collectProjections (Projections& proj, const Phantom& phm, const int iS
   xlat_mtx2 (temp, m_dXCenter, m_dYCenter);
   mult_mtx2 (rotmtx_initial, temp, rotmtx_initial);
 
   xlat_mtx2 (temp, m_dXCenter, m_dYCenter);
   mult_mtx2 (rotmtx_initial, temp, rotmtx_initial);
 

-  double xd1=0, yd1=0, xd2=0, yd2=0, dDetAngle=0;
-  if (m_idGeometry == GEOMETRY_EQUIANGULAR)
-      dDetAngle = m_initPos.angle + start_angle;
-  else {
+  double xd1=0, yd1=0, xd2=0, yd2=0;
+  if (m_idGeometry != GEOMETRY_EQUIANGULAR) {

       xd1 = m_initPos.xd1;
       yd1 = m_initPos.yd1;
       xd2 = m_initPos.xd2;
       xd1 = m_initPos.xd1;
       yd1 = m_initPos.yd1;
       xd2 = m_initPos.xd2;


@@ -303,12 +308,12 @@ Scanner::collectProjections (Projections& proj, const Phantom& phm, const int iS
     m_pSGP->eraseWindow ();
     m_pSGP->setWindow (m_dXMinWin, m_dYMinWin, m_dXMaxWin, m_dYMaxWin);
     m_pSGP->setRasterOp (RO_COPY);
     m_pSGP->eraseWindow ();
     m_pSGP->setWindow (m_dXMinWin, m_dYMinWin, m_dXMaxWin, m_dYMaxWin);
     m_pSGP->setRasterOp (RO_COPY);

-    m_pSGP->setColor (C_BLUE);
+    m_pSGP->setColor (C_RED);

     m_pSGP->moveAbs (0., 0.);
     m_pSGP->drawRect (m_dXCenter - dHalfPhmLen, m_dYCenter - dHalfPhmLen, m_dXCenter + dHalfPhmLen, m_dYCenter + dHalfPhmLen);
     m_pSGP->moveAbs (0., 0.);
     m_pSGP->drawCircle (m_dFocalLength);
     m_pSGP->moveAbs (0., 0.);
     m_pSGP->drawRect (m_dXCenter - dHalfPhmLen, m_dYCenter - dHalfPhmLen, m_dXCenter + dHalfPhmLen, m_dYCenter + dHalfPhmLen);
     m_pSGP->moveAbs (0., 0.);
     m_pSGP->drawCircle (m_dFocalLength);

-    m_pSGP->setColor (C_RED);
+    m_pSGP->setColor (C_BLUE);

     phm.draw (*m_pSGP);
     m_dTextHeight = m_pSGP->getCharHeight ();
 
     phm.draw (*m_pSGP);
     m_dTextHeight = m_pSGP->getCharHeight ();
 


@@ -348,15 +353,15 @@ Scanner::collectProjections (Projections& proj, const Phantom& phm, const int iS
        m_pSGP->lineAbs (xs2, ys2);
        m_pSGP->setPenWidth (2);
        m_pSGP->moveAbs (0., 0.);
        m_pSGP->lineAbs (xs2, ys2);
        m_pSGP->setPenWidth (2);
        m_pSGP->moveAbs (0., 0.);

-       m_pSGP->drawArc (m_dFocalLength, start_angle + m_initPos.angle, start_angle - m_initPos.angle);
+       m_pSGP->drawArc (m_dFocalLength, viewAngle + 3 * HALFPI - (m_dAngularDetLen/2), viewAngle + 3 * HALFPI + (m_dAngularDetLen/2));

       }
       m_pSGP->setPenWidth (1);
     }
     if (m_trace >= Trace::TRACE_CONSOLE)
       }
       m_pSGP->setPenWidth (1);
     }
     if (m_trace >= Trace::TRACE_CONSOLE)

-      traceShowParam ("Current View:", "%d (%.0f%%)", PROJECTION_TRACE_ROW_CURR_VIEW, C_RED, iView + iStartView, (iView + iStartView) / m_nView * 100.);
+      traceShowParam ("Current View:", "%d (%.0f%%)", PROJECTION_TRACE_ROW_CURR_VIEW, C_RED, iView + iStartView, (iView + iStartView) / static_cast<double>(m_nView) * 100.);

 #endif
            
 #endif
            

-    projectSingleView (phm, detArray, xd1, yd1, xd2, yd2, xs1, ys1, xs2, ys2, dDetAngle);
+    projectSingleView (phm, detArray, xd1, yd1, xd2, yd2, xs1, ys1, xs2, ys2, viewAngle + 3 * HALFPI);

     detArray.setViewAngle (viewAngle);
       
 #ifdef HAVE_SGP
     detArray.setViewAngle (viewAngle);
       
 #ifdef HAVE_SGP


@@ -366,9 +371,7 @@ Scanner::collectProjections (Projections& proj, const Phantom& phm, const int iS
 #endif
     xform_mtx2 (m_rotmtxIncrement, xs1, ys1);
     xform_mtx2 (m_rotmtxIncrement, xs2, ys2);
 #endif
     xform_mtx2 (m_rotmtxIncrement, xs1, ys1);
     xform_mtx2 (m_rotmtxIncrement, xs2, ys2);

-    if (m_idGeometry == GEOMETRY_EQUIANGULAR)
-       dDetAngle += m_detInc;
-    else {
+    if (m_idGeometry != GEOMETRY_EQUIANGULAR) {

        xform_mtx2 (m_rotmtxIncrement, xd1, yd1);  // rotate detector endpoints 
        xform_mtx2 (m_rotmtxIncrement, xd2, yd2);
     }
        xform_mtx2 (m_rotmtxIncrement, xd1, yd1);  // rotate detector endpoints 
        xform_mtx2 (m_rotmtxIncrement, xd2, yd2);
     }


@@ -411,10 +414,10 @@ Scanner::projectSingleView (const Phantom& phm, DetectorArray& detArray, const d
   double ddx=0, ddy=0, ddx2=0, ddy2=0, ddx2_ofs=0, ddy2_ofs=0, xd_maj=0, yd_maj=0;
   double dAngleInc=0, dAngleSampleInc=0, dAngleSampleOffset=0, dAngleMajor=0;
   if (m_idGeometry == GEOMETRY_EQUIANGULAR) {
   double ddx=0, ddy=0, ddx2=0, ddy2=0, ddx2_ofs=0, ddy2_ofs=0, xd_maj=0, yd_maj=0;
   double dAngleInc=0, dAngleSampleInc=0, dAngleSampleOffset=0, dAngleMajor=0;
   if (m_idGeometry == GEOMETRY_EQUIANGULAR) {

-      dAngleInc = m_detInc;
+      dAngleInc = m_dAngularDetIncrement;

       dAngleSampleInc = dAngleInc / m_nSample;
       dAngleSampleOffset = dAngleSampleInc / 2;
       dAngleSampleInc = dAngleInc / m_nSample;
       dAngleSampleOffset = dAngleSampleInc / 2;

-      dAngleMajor = dDetAngle + dAngleSampleOffset;
+      dAngleMajor = dDetAngle - (m_dAngularDetLen/2) + dAngleSampleOffset;

   } else {
       ddx = (xd2 - xd1) / detArray.nDet();  // change in coords 
       ddy = (yd2 - yd1) / detArray.nDet();  // between detectors
   } else {
       ddx = (xd2 - xd1) / detArray.nDet();  // change in coords 
       ddy = (yd2 - yd1) / detArray.nDet();  // between detectors


@@ -456,10 +459,9 @@ Scanner::projectSingleView (const Phantom& phm, DetectorArray& detArray, const d
 #ifdef HAVE_SGP
        if (m_pSGP && m_trace >= Trace::TRACE_PROJECTIONS) {
          m_pSGP->setColor (C_YELLOW);
 #ifdef HAVE_SGP
        if (m_pSGP && m_trace >= Trace::TRACE_PROJECTIONS) {
          m_pSGP->setColor (C_YELLOW);

-         m_pSGP->setRasterOp (RO_OR_REVERSE);
+         m_pSGP->setRasterOp (RO_AND);

          m_pSGP->moveAbs (xs, ys);
          m_pSGP->lineAbs (xd, yd);
          m_pSGP->moveAbs (xs, ys);
          m_pSGP->lineAbs (xd, yd);

-         m_pSGP->setRasterOp (RO_SET);

        }
 #endif
 
        }
 #endif
 


@@ -470,14 +472,6 @@ Scanner::projectSingleView (const Phantom& phm, DetectorArray& detArray, const d
          traceShowParam ("Attenuation:", "%s", PROJECTION_TRACE_ROW_ATTEN, C_LTMAGENTA, "        ");
          traceShowParam ("Attenuation:", "%.3f", PROJECTION_TRACE_ROW_ATTEN, C_LTMAGENTA, sum);
        }
          traceShowParam ("Attenuation:", "%s", PROJECTION_TRACE_ROW_ATTEN, C_LTMAGENTA, "        ");
          traceShowParam ("Attenuation:", "%.3f", PROJECTION_TRACE_ROW_ATTEN, C_LTMAGENTA, sum);
        }

-
-       //      if (m_pSGP && m_trace >= Trace::TRACE_PROJECTIONS) {
-       //        m_pSGP->setColor (C_YELLOW);
-       //        m_pSGP->setRasterOp (RO_XOR);
-       //        m_pSGP->moveAbs (xs, ys);
-       //        m_pSGP->lineAbs (xd, yd);
-       //        m_pSGP->setRasterOp (RO_SET);
-       //      }

 #endif
        if (m_idGeometry == GEOMETRY_EQUIANGULAR)
            dAngle += dAngleSampleInc;
 #endif
        if (m_idGeometry == GEOMETRY_EQUIANGULAR)
            dAngle += dAngleSampleInc;

diff --git a/libctsupport/mathfuncs.cpp b/libctsupport/mathfuncs.cpp
index bc3b38308bfde55d9b6b2f7002d1801b15661b98..1773c8b7f6cbc9db54bf7aacb15f297bc776235f 100644 (file)
--- a/libctsupport/mathfuncs.cpp
+++ b/libctsupport/mathfuncs.cpp
@@ -2,7 +2,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: mathfuncs.cpp,v 1.1 2000/06/22 10:17:28 kevin Exp $
+**  $Id: mathfuncs.cpp,v 1.2 2000/08/31 08:38:58 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
 **
 **  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


@@ -63,10 +63,10 @@ integrateSimpson (const double xmin, const double xmax, const double *y, const i
 
 
 /* NAME
 
 
 /* NAME

- *    norm_angle       Normalize angle to 0 to 2 * PI range
+ *    normalizeAngle       Normalize angle to 0 to 2 * PI range

  *
  * SYNOPSIS
  *
  * SYNOPSIS

- *    t = norm_angle (theta)
+ *    t = normalizeAngle (theta)

  *    double t        Normalized angle
  *    double theta     Input angle
  */
  *    double t        Normalized angle
  *    double theta     Input angle
  */

diff --git a/libctsupport/syserror.cpp b/libctsupport/syserror.cpp
index e59449f0161c5695a4efa4b471f5af34acadbc12..939f024bee0b29a6a6bfcb99a4b9ee7fe7bb276f 100644 (file)
--- a/libctsupport/syserror.cpp
+++ b/libctsupport/syserror.cpp
@@ -2,7 +2,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: syserror.cpp,v 1.4 2000/08/25 15:59:13 kevin Exp $
+**  $Id: syserror.cpp,v 1.5 2000/08/31 08:38:58 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
 **
 **  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


@@ -66,7 +66,7 @@ void sys_verror (int severity, const char *msg, va_list arg)
       return;
     else if (nErrorCount == MAX_ERROR_COUNT) {
       cout << "*****************************************************************" << endl;
       return;
     else if (nErrorCount == MAX_ERROR_COUNT) {
       cout << "*****************************************************************" << endl;

-      cout << "***     M A X I M U M  E R R O R  C O U N T  R E A C H E D    ***" << endl;
+      cout << "***   M A X I M U M   E R R O R   C O U N T   R E A C H E D   ***" << endl;

       cout << "***                                                           ***" << endl;
       cout << "***           No further errors will be reported              ***" << endl;
       cout << "*****************************************************************" << endl;
       cout << "***                                                           ***" << endl;
       cout << "***           No further errors will be reported              ***" << endl;
       cout << "*****************************************************************" << endl;

diff --git a/src/ctsim.cpp b/src/ctsim.cpp
index bb3b841dba227b70c364ad6df08bd13d749ddeb1..385fbee15445f120343b24bf7830efb62d649474 100644 (file)
--- a/src/ctsim.cpp
+++ b/src/ctsim.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: ctsim.cpp,v 1.9 2000/07/23 01:49:03 kevin Exp $
+**  $Id: ctsim.cpp,v 1.10 2000/08/31 08:38:58 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
 **
 **  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


@@ -53,7 +53,7 @@
 #include <getopt.h>
 #endif
 
 #include <getopt.h>
 #endif
 

-static const char* rcsindent = "$Id";
+static const char* rcsindent = "$Id: ctsim.cpp,v 1.10 2000/08/31 08:38:58 kevin Exp $";

 
 class CTSimApp* theApp = NULL;
 
 
 class CTSimApp* theApp = NULL;
 

diff --git a/src/dialogs.cpp b/src/dialogs.cpp
index a3c27f748a43c760b6a533d137158693d181c8f0..57c95f8b1bbda0d5aee2cbe57b8c37ff5c4828f2 100644 (file)
--- a/src/dialogs.cpp
+++ b/src/dialogs.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: dialogs.cpp,v 1.10 2000/08/27 20:32:55 kevin Exp $
+**  $Id: dialogs.cpp,v 1.11 2000/08/31 08:38:58 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
 **
 **  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


@@ -630,7 +630,7 @@ DialogGetReconstructionParameters::getFilterGenerationName (void)
 }
 
 
 }
 
 

-DialogAutoScaleParameters::DialogAutoScaleParameters (wxFrame *pParent, const ImageFile& rIF)
+DialogAutoScaleParameters::DialogAutoScaleParameters (wxFrame *pParent, const ImageFile& rIF, double dDefaultScaleFactor)

   : wxDialog (pParent, -1, "Auto Scale Parameters", wxDefaultPosition, wxDefaultSize, wxDEFAULT_DIALOG_STYLE | wxCAPTION), m_rImageFile(rIF)
 {
   wxBoxSizer* pTopSizer = new wxBoxSizer (wxVERTICAL);
   : wxDialog (pParent, -1, "Auto Scale Parameters", wxDefaultPosition, wxDefaultSize, wxDEFAULT_DIALOG_STYLE | wxCAPTION), m_rImageFile(rIF)
 {
   wxBoxSizer* pTopSizer = new wxBoxSizer (wxVERTICAL);


@@ -650,7 +650,9 @@ DialogAutoScaleParameters::DialogAutoScaleParameters (wxFrame *pParent, const Im
 
   wxGridSizer *pGridSizer = new wxGridSizer (2);
   pGridSizer->Add (new wxStaticText (this, -1, "Standard Deviation Factor"), 0, wxALIGN_RIGHT | wxALIGN_CENTER_VERTICAL);
 
   wxGridSizer *pGridSizer = new wxGridSizer (2);
   pGridSizer->Add (new wxStaticText (this, -1, "Standard Deviation Factor"), 0, wxALIGN_RIGHT | wxALIGN_CENTER_VERTICAL);

-  m_pTextCtrlStdDevFactor = new wxTextCtrl (this, -1, "1.0", wxDefaultPosition, wxSize(100, 25), 0);
+  ostringstream osDefaultFactor;
+  osDefaultFactor << dDefaultScaleFactor;
+  m_pTextCtrlStdDevFactor = new wxTextCtrl (this, -1, osDefaultFactor.str().c_str(), wxDefaultPosition, wxSize(100, 25), 0);

   pGridSizer->Add (m_pTextCtrlStdDevFactor, 0, wxALIGN_CENTER_VERTICAL);
   pTopSizer->Add (pGridSizer, 1, wxALL, 10);
 
   pGridSizer->Add (m_pTextCtrlStdDevFactor, 0, wxALIGN_CENTER_VERTICAL);
   pTopSizer->Add (pGridSizer, 1, wxALL, 10);
 


@@ -694,3 +696,15 @@ DialogAutoScaleParameters::getMinMax (double* pMin, double* pMax)
   *pMax = dCenter + dHalfWidth;
   *theApp->getLog() << "Setting minimum to " << *pMin << " and maximum to " << *pMax << "\n";
 }
   *pMax = dCenter + dHalfWidth;
   *theApp->getLog() << "Setting minimum to " << *pMin << " and maximum to " << *pMax << "\n";
 }

+
+double
+DialogAutoScaleParameters::getAutoScaleFactor ()
+{
+  wxString sStddevFactor = m_pTextCtrlStdDevFactor->GetValue();
+  double dValue = 1.;
+  if (! sStddevFactor.ToDouble (&dValue)) {
+    *theApp->getLog() << "Error: Non-numeric Standard Deviation Factor of " << sStddevFactor << "\n";
+  }
+
+  return dValue;
+}

diff --git a/src/dialogs.h b/src/dialogs.h
index 77b7cbf8e4cb6722b21b67b7c46cd2e773f81150..7271339bd9a8c789e254157e8e22b38ffef5d948 100644 (file)
--- a/src/dialogs.h
+++ b/src/dialogs.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: dialogs.h,v 1.11 2000/08/27 20:32:55 kevin Exp $
+**  $Id: dialogs.h,v 1.12 2000/08/31 08:38:58 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
 **
 **  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


@@ -185,10 +185,11 @@ class DialogGetReconstructionParameters : public wxDialog
 class DialogAutoScaleParameters : public wxDialog
 {
  public:
 class DialogAutoScaleParameters : public wxDialog
 {
  public:

-    DialogAutoScaleParameters (wxFrame* pParent, const ImageFile& rImageFile);
+    DialogAutoScaleParameters (wxFrame* pParent, const ImageFile& rImageFile, double dDefaultScaleFactor = 1.);

     virtual ~DialogAutoScaleParameters() {}
 
     void getMinMax (double* pMin, double* pMax);
     virtual ~DialogAutoScaleParameters() {}
 
     void getMinMax (double* pMin, double* pMax);

+    double getAutoScaleFactor ();

 
  private:
     const ImageFile& m_rImageFile;
 
  private:
     const ImageFile& m_rImageFile;

diff --git a/src/dlgprojections.cpp b/src/dlgprojections.cpp
index 8f66dc3caa474f53f474231be32a580c4c2a9c28..4f793342a37e3e2ad031cb284247a79242077dda 100644 (file)
--- a/src/dlgprojections.cpp
+++ b/src/dlgprojections.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: dlgprojections.cpp,v 1.1 2000/08/27 20:32:55 kevin Exp $
+**  $Id: dlgprojections.cpp,v 1.2 2000/08/31 08:38:58 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
 **
 **  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


@@ -58,17 +58,20 @@ static const int LAYOUT_Y_MARGIN = 4;
 
 BEGIN_EVENT_TABLE(ProjectionsDialog, wxDialog)
    EVT_BUTTON(wxID_CANCEL, ProjectionsDialog::OnCancel)
 
 BEGIN_EVENT_TABLE(ProjectionsDialog, wxDialog)
    EVT_BUTTON(wxID_CANCEL, ProjectionsDialog::OnCancel)

+   EVT_BUTTON(ID_BTN_PAUSE, ProjectionsDialog::OnPause)
+   EVT_BUTTON(ID_BTN_STEP, ProjectionsDialog::OnStep)

    EVT_CLOSE(ProjectionsDialog::OnClose)
    EVT_CLOSE(ProjectionsDialog::OnClose)

+   EVT_PAINT(ProjectionsDialog::OnPaint)

 END_EVENT_TABLE()
 
 IMPLEMENT_CLASS(ProjectionsDialog, wxDialog)
 
 
 ProjectionsDialog::ProjectionsDialog (Scanner& rScanner, Projections& rProj, const Phantom& rPhantom, const int iTrace, wxWindow *parent)
 END_EVENT_TABLE()
 
 IMPLEMENT_CLASS(ProjectionsDialog, wxDialog)
 
 
 ProjectionsDialog::ProjectionsDialog (Scanner& rScanner, Projections& rProj, const Phantom& rPhantom, const int iTrace, wxWindow *parent)

-  : wxDialog(parent, -1, "Collect Projections"), m_rScanner(rScanner), m_rProjections(rProj), m_rPhantom(rPhantom), m_pSGPDriver(NULL), m_pSGP(NULL), m_iTrace(iTrace), m_pDC(NULL)
+    : wxDialog(parent, -1, "Collect Projections"), m_rScanner(rScanner), m_rProjections(rProj), m_rPhantom(rPhantom), m_pSGPDriver(NULL), m_pSGP(NULL), m_iTrace(iTrace), m_pDC(NULL), m_btnAbort(0), m_btnPause(0), m_btnStep(0)

 {
     m_state = Continue;
 {
     m_state = Continue;

-
+    m_iLastView = -1;

     m_parentTop = parent;
     while ( m_parentTop && m_parentTop->GetParent() )
         m_parentTop = m_parentTop->GetParent();
     m_parentTop = parent;
     while ( m_parentTop && m_parentTop->GetParent() )
         m_parentTop = m_parentTop->GetParent();


@@ -84,6 +87,22 @@ ProjectionsDialog::ProjectionsDialog (Scanner& rScanner, Projections& rProj, con
 
     m_btnAbort->SetConstraints(c);
 
 
     m_btnAbort->SetConstraints(c);
 

+    m_btnPause = new wxButton (this, ID_BTN_PAUSE, wxString("Pause"));
+    wxLayoutConstraints* cPause = new wxLayoutConstraints;
+    cPause->right.SameAs(this, wxRight, 3*LAYOUT_X_MARGIN + sizeBtn.x);
+    cPause->bottom.SameAs(this, wxBottom, 2*LAYOUT_Y_MARGIN);
+    cPause->width.Absolute(sizeBtn.x);
+    cPause->height.Absolute(sizeBtn.y);
+    m_btnPause->SetConstraints(cPause);
+
+    m_btnStep = new wxButton (this, ID_BTN_STEP, wxString("Step"));
+    wxLayoutConstraints* cStep = new wxLayoutConstraints;
+    cStep->right.SameAs(this, wxRight, 5*LAYOUT_X_MARGIN + sizeBtn.x * 2);
+    cStep->bottom.SameAs(this, wxBottom, 2*LAYOUT_Y_MARGIN);
+    cStep->width.Absolute(sizeBtn.x);
+    cStep->height.Absolute(sizeBtn.y);
+    m_btnStep->SetConstraints(cStep);
+

     SetAutoLayout(TRUE);
     Layout();
 
     SetAutoLayout(TRUE);
     Layout();
 


@@ -92,7 +111,13 @@ ProjectionsDialog::ProjectionsDialog (Scanner& rScanner, Projections& rProj, con
       sizeDlg.x = max(sizeDlg.x,sizeDlg.y);
       sizeDlg.y = max(sizeDlg.x,sizeDlg.y);
     }
       sizeDlg.x = max(sizeDlg.x,sizeDlg.y);
       sizeDlg.y = max(sizeDlg.x,sizeDlg.y);
     }

+    if (m_iTrace >= Trace::TRACE_PLOT)
+      sizeDlg.x += 200;
+
+    m_iClientX = sizeDlg.x;
+    m_iClientY = sizeDlg.y;

     SetClientSize(sizeDlg);
     SetClientSize(sizeDlg);

+    m_bitmap.Create (m_iClientX, m_iClientY); // save a copy of screen

 
     Centre(wxCENTER_FRAME | wxBOTH);
 
 
     Centre(wxCENTER_FRAME | wxBOTH);
 


@@ -110,34 +135,64 @@ ProjectionsDialog::ProjectionsDialog (Scanner& rScanner, Projections& rProj, con
 
     wxYield();     // Update the display
 
 
     wxYield();     // Update the display
 

+    m_pSGPDriver->idWX()->SetFont(*wxSWISS_FONT);

 #ifdef __WXMAC__
     MacUpdateImmediately();
 #endif
 }
 
 #ifdef __WXMAC__
     MacUpdateImmediately();
 #endif
 }
 

+void
+ProjectionsDialog::showView (int iViewNumber)
+{
+    if ( iViewNumber < m_rProjections.nView() ) {
+       wxYield();        // update the display
+       m_iLastView = iViewNumber;
+       if (m_iTrace >= Trace::TRACE_PLOT)
+           m_pSGP->setViewport (0, 0, 0.75, 1);
+       m_rScanner.collectProjections (m_rProjections, m_rPhantom, iViewNumber, 1, true, m_iTrace, m_pSGP);
+       if (m_iTrace >= Trace::TRACE_PLOT) {
+           const DetectorArray& detArray = m_rProjections.getDetectorArray (iViewNumber);
+           const DetectorValue* detValues = detArray.detValues();
+           double detPos [detArray.nDet()];
+           for (int i = 0; i < detArray.nDet(); i++)
+               detPos[i] = i;
+           EZPlot ezplot (*m_pSGP);
+           ezplot.ezset("xporigin 0.75");
+           ezplot.ezset("yporigin 0.10");
+           ezplot.ezset("xlength  0.25");
+           ezplot.ezset("ylength  0.90");
+           ezplot.ezset("grid");
+           ezplot.ezset("box");
+           ezplot.addCurve (detValues, detPos, detArray.nDet());
+           ezplot.plot();
+       }
+    }
+}
+

 bool
 ProjectionsDialog::projectView (int iViewNumber)
 {
 bool
 ProjectionsDialog::projectView (int iViewNumber)
 {

-    if ( iViewNumber < m_rProjections.nView() ) {
-      m_rScanner.collectProjections (m_rProjections, m_rPhantom, iViewNumber, 1, true, m_iTrace, m_pSGP);
-      wxYield();        // update the display
-    } else {
-      m_state = Finished;    // so that we return TRUE below and 
-                             // that [Cancel] handler knew what to do
+    showView (iViewNumber);
+    wxYield();        // update the display
+    if (m_iTrace >= Trace::TRACE_PLOT)
+       sleep(1);
+    else {
+       m_state = Finished;    // so that we return TRUE below and 
+       // that [Cancel] handler knew what to do

 #if 0
 #if 0

-      if ( m_btnAbort ) 
-       m_btnAbort->SetLabel(_("Close")); // tell the user what he should do...
-      wxYield();
-      
-      (void)ShowModal();
+       if ( m_btnAbort ) 
+           m_btnAbort->SetLabel(_("Close")); // tell the user what he should do...
+       wxYield();
+       
+       (void)ShowModal();

 #endif
     }
 #endif
     }

-      
+    

 #ifdef __WXMAC__
     MacUpdateImmediately();
 #endif
 
 #ifdef __WXMAC__
     MacUpdateImmediately();
 #endif
 

-    return m_state != Canceled;
+    return m_state != Cancelled;

 }
 
 
 }
 
 


@@ -152,7 +207,7 @@ void ProjectionsDialog::OnCancel (wxCommandEvent& event)
   } else {
     // request to cancel was received, the next time Update() is called we
     // will handle it
   } else {
     // request to cancel was received, the next time Update() is called we
     // will handle it

-    m_state = Canceled;
+    m_state = Cancelled;

 
     // update the button state immediately so that the user knows that the
     // request has been noticed
 
     // update the button state immediately so that the user knows that the
     // request has been noticed


@@ -160,14 +215,71 @@ void ProjectionsDialog::OnCancel (wxCommandEvent& event)
   }
 }
 
   }
 }
 

+
+void 
+ProjectionsDialog::OnPause (wxCommandEvent& event)
+{
+  if ( m_state == Finished ) {
+    // this means that the count down is already finished and we're being
+    // shown as a modal dialog - so just let the default handler do the job
+    event.Skip();
+  } else {
+      if (m_state == Continue) {
+         m_memoryDC.SelectObject (m_bitmap);       // in memoryDC
+         m_pSGP->setDC (&m_memoryDC);
+         m_memoryDC.SetFont (*wxSWISS_FONT);
+         showView (m_iLastView);
+         m_pSGP->setDC (m_pDC);
+         m_memoryDC.SelectObject(wxNullBitmap);
+         m_state = Paused;
+         m_btnPause->SetLabel (wxString("Resume"));
+      } else if (m_state == Paused) {
+         m_state = Continue;
+         m_btnPause->SetLabel (wxString("Pause"));
+      }
+  }
+}
+
+void 
+ProjectionsDialog::OnStep (wxCommandEvent& event)
+{
+  if ( m_state == Finished ) {
+    // this means that the count down is already finished and we're being
+    // shown as a modal dialog - so just let the default handler do the job
+    event.Skip();
+  } else {
+      if (m_state == Continue) {
+         m_memoryDC.SelectObject (m_bitmap);       // in memoryDC
+         m_pSGP->setDC (&m_memoryDC);
+         m_memoryDC.SetFont (*wxSWISS_FONT);
+         m_rScanner.collectProjections (m_rProjections, m_rPhantom, m_iLastView, 1, true, m_iTrace, m_pSGP);
+         m_pSGP->setDC (m_pDC);
+         m_memoryDC.SelectObject(wxNullBitmap);
+         m_state = Paused;
+         m_btnPause->SetLabel (wxString("Resume"));
+      } else if (m_state == Paused) {
+         projectView (m_iLastView + 1);
+      }
+  }
+}
+

 void ProjectionsDialog::OnClose(wxCloseEvent& event)
 {
 void ProjectionsDialog::OnClose(wxCloseEvent& event)
 {

-    if ( m_state == Uncancelable )
+    if ( m_state == Uncancellable )

       event.Veto(TRUE);    // can't close this dialog
     else if ( m_state == Finished )
       event.Skip(); // let the default handler close the window as we already terminated
     else
       event.Veto(TRUE);    // can't close this dialog
     else if ( m_state == Finished )
       event.Skip(); // let the default handler close the window as we already terminated
     else

-      m_state = Canceled;          // next Update() will notice it
+      m_state = Cancelled;          // next Update() will notice it
+}
+
+void
+ProjectionsDialog::OnPaint (wxPaintEvent& event)
+{
+  wxPaintDC paintDC (this);
+  if (m_state == Paused) {
+    paintDC.DrawBitmap(m_bitmap, 0, 0, false);
+  }

 }
 
 
 }
 
 

diff --git a/src/dlgprojections.h b/src/dlgprojections.h
index 440d2ff2ceaf744759017ca387757b46060dd43c..b8798bb66fd486f65a46cc389b51322473614e60 100644 (file)
--- a/src/dlgprojections.h
+++ b/src/dlgprojections.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: dlgprojections.h,v 1.1 2000/08/27 20:32:55 kevin Exp $
+**  $Id: dlgprojections.h,v 1.2 2000/08/31 08:38:58 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
 **
 **  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


@@ -34,6 +34,7 @@
 
 #include "wx/setup.h"
 #include "wx/dialog.h"
 
 #include "wx/setup.h"
 #include "wx/dialog.h"

+#include "wx/dcmemory.h"

 
 class wxButton;
 class wxStaticText;
 
 class wxButton;
 class wxStaticText;


@@ -69,30 +70,51 @@ public:
    void OnCancel(wxCommandEvent& event);
        // callback to disable "hard" window closing
    void OnClose(wxCloseEvent& event);
    void OnCancel(wxCommandEvent& event);
        // callback to disable "hard" window closing
    void OnClose(wxCloseEvent& event);

+   void OnPaint(wxPaintEvent& event);
+
+   void OnPause(wxCommandEvent& event);
+   void OnStep(wxCommandEvent& event);
+
+   bool isPaused() const {return m_state == Paused;}
+
+   bool isCancelled() const {return m_state == Cancelled;}

 
 private:
    // parent top level window (may be NULL)
    wxWindow *m_parentTop;
 
 private:
    // parent top level window (may be NULL)
    wxWindow *m_parentTop;

+   int m_iLastView;
+   int m_iClientX;   // size of client window
+   int m_iClientY;
+
+   Scanner& m_rScanner;
+   Projections& m_rProjections;
+   const Phantom& m_rPhantom;
+   SGPDriver* m_pSGPDriver;
+   SGP* m_pSGP;
+   const int m_iTrace;
+   wxDC* m_pDC;
+
+   wxButton *m_btnAbort;    // the abort button (or NULL if none)
+   wxButton *m_btnPause; 
+   wxButton *m_btnStep;
+
+   wxMemoryDC m_memoryDC;  // for restoring image on OnPaint
+   wxBitmap m_bitmap;

 
    // continue processing or not (return value for Update())
    enum
    {
 
    // continue processing or not (return value for Update())
    enum
    {

-      Uncancelable = -1,   // dialog can't be canceled
-      Canceled,            // can be cancelled and, in fact, was
+      Uncancellable = -1,   // dialog can't be canceled
+      Paused,
+      Cancelled,            // can be cancelled and, in fact, was

       Continue,            // can be cancelled but wasn't
       Finished             // finished, waiting to be removed from screen
    } m_state;
 
       Continue,            // can be cancelled but wasn't
       Finished             // finished, waiting to be removed from screen
    } m_state;
 

-   // the abort button (or NULL if none)
-   wxButton *m_btnAbort;
+   const static int ID_BTN_PAUSE = 19998;
+   const static int ID_BTN_STEP = 19999;

 
 

-   Scanner& m_rScanner;
-   Projections& m_rProjections;
-   const Phantom& m_rPhantom;
-   SGPDriver* m_pSGPDriver;
-   SGP* m_pSGP;
-   const int m_iTrace;
-   wxDC* m_pDC;
+   void showView (int iViewNumber);

 
    DECLARE_EVENT_TABLE()
 };
 
    DECLARE_EVENT_TABLE()
 };

diff --git a/src/views.cpp b/src/views.cpp
index edc35e6a00bdb417e313eed63798fd412bc2251c..3df51fa5e3c94227185b8acb602024602de88d15 100644 (file)
--- a/src/views.cpp
+++ b/src/views.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: views.cpp,v 1.17 2000/08/27 20:32:55 kevin Exp $
+**  $Id: views.cpp,v 1.18 2000/08/31 08:38:58 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
 **
 **  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


@@ -152,7 +152,7 @@ void
 ImageFileView::OnScaleAuto (wxCommandEvent& event)
 {
     const ImageFile& rIF = GetDocument()->getImageFile();
 ImageFileView::OnScaleAuto (wxCommandEvent& event)
 {
     const ImageFile& rIF = GetDocument()->getImageFile();

-    DialogAutoScaleParameters dialogAutoScale (m_frame, rIF);
+    DialogAutoScaleParameters dialogAutoScale (m_frame, rIF, m_dAutoScaleFactor);

     int iRetVal = dialogAutoScale.ShowModal();
     if (iRetVal == wxID_OK) {
       m_bMinSpecified = true;
     int iRetVal = dialogAutoScale.ShowModal();
     if (iRetVal == wxID_OK) {
       m_bMinSpecified = true;


@@ -161,6 +161,7 @@ ImageFileView::OnScaleAuto (wxCommandEvent& event)
       dialogAutoScale.getMinMax (&dMin, &dMax);
       m_dMinPixel = dMin;
       m_dMaxPixel = dMax;
       dialogAutoScale.getMinMax (&dMin, &dMax);
       m_dMinPixel = dMin;
       m_dMaxPixel = dMax;

+      m_dAutoScaleFactor = dialogAutoScale.getAutoScaleFactor();

       OnUpdate (this, NULL);
     }
 }
       OnUpdate (this, NULL);
     }
 }


@@ -256,6 +257,7 @@ ImageFileView::OnCreate(wxDocument *doc, long WXUNUSED(flags) )
     
     m_bMinSpecified = false;
     m_bMaxSpecified = false;
     
     m_bMinSpecified = false;
     m_bMaxSpecified = false;

+    m_dAutoScaleFactor = 1.;

 
     int width, height;
     m_frame->GetClientSize(&width, &height);
 
     int width, height;
     m_frame->GetClientSize(&width, &height);


@@ -430,14 +432,21 @@ PhantomView::OnProjections (wxCommandEvent& event)
       if (iTrace > Trace::TRACE_CONSOLE) {
        ProjectionsDialog dialogProjections (theScanner, rProj, rPhantom, iTrace, dynamic_cast<wxWindow*>(m_frame));
        for (int iView = 0; iView < rProj.nView(); iView++) {
       if (iTrace > Trace::TRACE_CONSOLE) {
        ProjectionsDialog dialogProjections (theScanner, rProj, rPhantom, iTrace, dynamic_cast<wxWindow*>(m_frame));
        for (int iView = 0; iView < rProj.nView(); iView++) {

-         if (! dialogProjections.projectView (iView)) {
+         ::wxYield();
+         ::wxYield();
+         if (dialogProjections.isCancelled() || ! dialogProjections.projectView (iView)) {

            pProjectionDoc->DeleteAllViews();
            return;
          }
          ::wxYield();
            pProjectionDoc->DeleteAllViews();
            return;
          }
          ::wxYield();

+         ::wxYield();
+         while (dialogProjections.isPaused()) {
+             ::wxYield();
+             ::wxUsleep(50);
+         }

        }
       } else
        }
       } else

-       theScanner.collectProjections (rProj, rPhantom);
+       theScanner.collectProjections (rProj, rPhantom, iTrace);

 
       pProjectionDoc->Modify(true);
       pProjectionDoc->UpdateAllViews(this);
 
       pProjectionDoc->Modify(true);
       pProjectionDoc->UpdateAllViews(this);


@@ -447,6 +456,7 @@ PhantomView::OnProjections (wxCommandEvent& event)
        }
       ostringstream os;
       os << "Projections for " << rPhantom.name() << ": nDet=" << nDet << ", nView=" << nView << ", nSamples=" << nSamples << ", RotAngle=" << dRotAngle << ", FocalLengthRatio=" << dFocalLengthRatio << ", FieldOfViewRatio=" << dFieldOfViewRatio << ", Geometry=" << sGeometry.c_str() << "\n";
        }
       ostringstream os;
       os << "Projections for " << rPhantom.name() << ": nDet=" << nDet << ", nView=" << nView << ", nSamples=" << nSamples << ", RotAngle=" << dRotAngle << ", FocalLengthRatio=" << dFocalLengthRatio << ", FieldOfViewRatio=" << dFieldOfViewRatio << ", Geometry=" << sGeometry.c_str() << "\n";

+      rProj.setRemark (os.str());

       *theApp->getLog() << os.str().c_str();
     }
   }
       *theApp->getLog() << os.str().c_str();
     }
   }


@@ -642,7 +652,7 @@ ProjectionFileView::OnProperties (wxCommandEvent& event)
 {
   const Projections& rProj = GetDocument()->getProjections();
   ostringstream os;
 {
   const Projections& rProj = GetDocument()->getProjections();
   ostringstream os;

-  os << "ProjectionFile " << rProj.getFilename() << ":  Number of Detectors = " << rProj.nDet() << ", Number of Views = " << rProj.nView() << "\n";
+  rProj.printScanInfo(os);

   *theApp->getLog() << os.str().c_str();
 }
 
   *theApp->getLog() << os.str().c_str();
 }
 

diff --git a/src/views.h b/src/views.h
index a11cfaf8d8df918d3dbc7899b0f8298ded9ef1d4..71efce20081a48f79fc7bc5ce942fce6825697c4 100644 (file)
--- a/src/views.h
+++ b/src/views.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: views.h,v 1.7 2000/08/25 15:59:13 kevin Exp $
+**  $Id: views.h,v 1.8 2000/08/31 08:38:58 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
 **
 **  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


@@ -54,6 +54,7 @@ private:
     bool m_bMaxSpecified;
     double m_dMinPixel;
     double m_dMaxPixel;
     bool m_bMaxSpecified;
     double m_dMinPixel;
     double m_dMaxPixel;

+    double m_dAutoScaleFactor;

 
 public:
     ImageFileView(void);
 
 public:
     ImageFileView(void);

diff --git a/tools/Makefile.am b/tools/Makefile.am
index a79c7e211b9b15c348d93b6b6cc07edbd345c81d..ff514447abde20be325180d114b4f98f87357076 100644 (file)
--- a/tools/Makefile.am
+++ b/tools/Makefile.am
@@ -1,4 +1,4 @@
-bin_PROGRAMS = pjrec phm2pj pj2if @lamprograms@  ifinfo phm2if if-1 if-2 if2img 
+bin_PROGRAMS = pjrec phm2pj pj2if @lamprograms@  ifinfo phm2if if-1 if-2 if2img pjinfo

 bin_SCRIPTS = sample-ctsim.sh
 EXTRA_PROGRAMS = pjrec-lam phm2pj-lam phm2if-lam
 INCLUDES=@my_includes@
 bin_SCRIPTS = sample-ctsim.sh
 EXTRA_PROGRAMS = pjrec-lam phm2pj-lam phm2if-lam
 INCLUDES=@my_includes@


@@ -25,6 +25,9 @@ if2img_DEPENDENCIES=$(SOURCE_DEPEND)
 pj2if_SOURCES = pj2if.cpp
 pj2if_LDADD=@ctlibs@
 pj2if_DEPENDENCIES=$(SOURCE_DEPEND)
 pj2if_SOURCES = pj2if.cpp
 pj2if_LDADD=@ctlibs@
 pj2if_DEPENDENCIES=$(SOURCE_DEPEND)

+pjinfo_SOURCES = pjinfo.cpp
+pjinfo_LDADD=@ctlibs@
+pjinfo_DEPENDENCIES=$(SOURCE_DEPEND)

 if_1_SOURCES=if-1.cpp 
 if_1_LDADD=@ctlibs@
 if_1_DEPENDENCIES=$(SOURCE_DEPEND)
 if_1_SOURCES=if-1.cpp 
 if_1_LDADD=@ctlibs@
 if_1_DEPENDENCIES=$(SOURCE_DEPEND)

diff --git a/tools/phm2pj.cpp b/tools/phm2pj.cpp
index 0459f53f7b10bdfa87c61656a8963bcdb112232d..c42f0af80618451f26cc973dd08355c42a5be2c2 100644 (file)
--- a/tools/phm2pj.cpp
+++ b/tools/phm2pj.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: phm2pj.cpp,v 1.11 2000/08/27 20:32:55 kevin Exp $
+**  $Id: phm2pj.cpp,v 1.12 2000/08/31 08:38:58 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
 **
 **  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


@@ -50,7 +50,7 @@ static struct option phm2pj_options[] =
   {0, 0, 0, 0}
 };
 
   {0, 0, 0, 0}
 };
 

-static const char* g_szIdStr = "$Id: phm2pj.cpp,v 1.11 2000/08/27 20:32:55 kevin Exp $";
+static const char* g_szIdStr = "$Id: phm2pj.cpp,v 1.12 2000/08/31 08:38:58 kevin Exp $";

 
 
 void 
 
 
 void 


@@ -235,7 +235,7 @@ phm2pj_main (int argc, char* argv[])
     }
 
     ostringstream desc;
     }
 
     ostringstream desc;

-    desc << "phm2pj: NDet=" << opt_ndet << ", Nview=" << opt_nview << ", NRay=" << opt_nray << ", RotAngle=" << opt_rotangle << ", ";
+    desc << "phm2pj: NDet=" << opt_ndet << ", Nview=" << opt_nview << ", NRay=" << opt_nray << ", RotAngle=" << opt_rotangle << ", Geometry=" << optGeometryName << ", ";

     if (optPhmFileName.length()) {
       desc << "PhantomFile=" << optPhmFileName;
     } else if (optPhmName != "") {
     if (optPhmFileName.length()) {
       desc << "PhantomFile=" << optPhmFileName;
     } else if (optPhmName != "") {


@@ -296,8 +296,11 @@ phm2pj_main (int argc, char* argv[])
   if (mpiWorld.getRank() == 0)
     pjGlobal.initFromScanner (scanner);
   
   if (mpiWorld.getRank() == 0)
     pjGlobal.initFromScanner (scanner);
   

-  if (opt_verbose)
-    pjGlobal.printScanInfo();
+  if (opt_verbose) {
+    ostringstream os;
+    pjGlobal.printScanInfo(os);
+    cout << os;
+  }

 
   Projections pjLocal (scanner);
   pjLocal.setNView (mpiWorld.getMyLocalWorkUnits());
 
   Projections pjLocal (scanner);
   pjLocal.setNView (mpiWorld.getMyLocalWorkUnits());


@@ -348,8 +351,9 @@ phm2pj_main (int argc, char* argv[])
       if (opt_verbose) {
        phm.print();
        cout << endl;
       if (opt_verbose) {
        phm.print();
        cout << endl;

-       pjGlobal.printScanInfo();
-       cout << endl;
+       ostringstream os;
+       pjGlobal.printScanInfo (os);
+       cout << os.str() << endl;

        cout << "  Remark: " << pjGlobal.remark() << endl;
        cout << "Run time: " << pjGlobal.calcTime() << " seconds\n";
       }
        cout << "  Remark: " << pjGlobal.remark() << endl;
        cout << "Run time: " << pjGlobal.calcTime() << " seconds\n";
       }

diff --git a/tools/pj2if.cpp b/tools/pj2if.cpp
index f65d55d12e9d6aa57e285e29ab2251bc54bec0b1..b9bd2e9e6987c4d8c48b77454f715b05207feb63 100644 (file)
--- a/tools/pj2if.cpp
+++ b/tools/pj2if.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: pj2if.cpp,v 1.3 2000/08/03 09:57:29 kevin Exp $
+**  $Id: pj2if.cpp,v 1.4 2000/08/31 08:38:58 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
 **
 **  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


@@ -36,17 +36,18 @@
 #include "timer.h"
 
 
 #include "timer.h"
 
 

-enum { O_VERBOSE, O_HELP, O_VERSION };
+enum { O_VERBOSE, O_DUMP, O_HELP, O_VERSION };

 
 static struct option my_options[] =
 {
   {"verbose", 0, 0, O_VERBOSE},
 
 static struct option my_options[] =
 {
   {"verbose", 0, 0, O_VERBOSE},

+  {"dump", 0, 0, O_DUMP},

   {"help", 0, 0, O_HELP},
   {"version", 0, 0, O_VERSION},
   {0, 0, 0, 0}
 };
 
   {"help", 0, 0, O_HELP},
   {"version", 0, 0, O_VERSION},
   {0, 0, 0, 0}
 };
 

-static const char* g_szIdStr = "$Id: pj2if.cpp,v 1.3 2000/08/03 09:57:29 kevin Exp $";
+static const char* g_szIdStr = "$Id: pj2if.cpp,v 1.4 2000/08/31 08:38:58 kevin Exp $";

 
 void 
 pj2if_usage (const char *program)
 
 void 
 pj2if_usage (const char *program)


@@ -55,6 +56,7 @@ pj2if_usage (const char *program)
   cout << "Converts a projection file to a IF file" << endl;
   cout << endl;
   cout << "   --verbose   Verbose mode" << endl;
   cout << "Converts a projection file to a IF file" << endl;
   cout << endl;
   cout << "   --verbose   Verbose mode" << endl;

+  cout << "   --dump      Dump all scan data" << endl;

   cout << "   --version   Print version" << endl;
   cout << "   --help      Print this help message" << endl;
 }
   cout << "   --version   Print version" << endl;
   cout << "   --help      Print this help message" << endl;
 }


@@ -66,6 +68,7 @@ pj2if_main (const int argc, char *const argv[])
 {
   char *pj_name, *im_name;
   bool optVerbose = false;
 {
   char *pj_name, *im_name;
   bool optVerbose = false;

+  bool optDump = false;

   extern int optind;
   Timer timerProgram;
 
   extern int optind;
   Timer timerProgram;
 


@@ -80,6 +83,9 @@ pj2if_main (const int argc, char *const argv[])
        case O_VERBOSE:
          optVerbose = true;
          break;
        case O_VERBOSE:
          optVerbose = true;
          break;

+       case O_DUMP:
+         optDump = true;
+         break;

         case O_VERSION:
 #ifdef VERSION
          cout << "Version " << VERSION << endl << g_szIdStr << endl;
         case O_VERSION:
 #ifdef VERSION
          cout << "Version " << VERSION << endl << g_szIdStr << endl;


@@ -111,8 +117,13 @@ pj2if_main (const int argc, char *const argv[])
     return (1);
   }
 
     return (1);
   }
 

-  if (optVerbose)
-      pj.printScanInfo();
+  if (optDump)
+    pj.printProjectionData();
+  else if (optVerbose) {
+    ostringstream os;
+    pj.printScanInfo (os);
+    cout << os.str();
+  }

   
   ImageFile im (pj.nDet(), pj.nView());
   
   
   ImageFile im (pj.nDet(), pj.nView());
   

diff --git a/tools/pjrec.cpp b/tools/pjrec.cpp
index 157d6584c1cdd28e8c06d3041817f0ed8a1503c7..9fcb8903a6de987e3a130f7ae720eeb28883dc45 100644 (file)
--- a/tools/pjrec.cpp
+++ b/tools/pjrec.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **

-**  $Id: pjrec.cpp,v 1.15 2000/08/27 20:32:55 kevin Exp $
+**  $Id: pjrec.cpp,v 1.16 2000/08/31 08:38:58 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
 **
 **  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


@@ -49,7 +49,7 @@ static struct option my_options[] =
   {0, 0, 0, 0}
 };
 
   {0, 0, 0, 0}
 };
 

-static const char* g_szIdStr = "$Id: pjrec.cpp,v 1.15 2000/08/27 20:32:55 kevin Exp $";
+static const char* g_szIdStr = "$Id: pjrec.cpp,v 1.16 2000/08/31 08:38:58 kevin Exp $";

 
 void 
 pjrec_usage (const char *program)
 
 void 
 pjrec_usage (const char *program)


@@ -260,8 +260,11 @@ pjrec_main (int argc, char * argv[])
 #ifdef HAVE_MPI
   if (mpiWorld.getRank() == 0) {
     projGlobal.read (pszFilenameProj);
 #ifdef HAVE_MPI
   if (mpiWorld.getRank() == 0) {
     projGlobal.read (pszFilenameProj);

-    if (bOptVerbose)
-      projGlobal.printScanInfo();
+    if (bOptVerbose) {
+      ostringstream os;
+      projGlobal.printScanInfo (os);
+      cout << os.str();
+    }

 
     mpi_ndet = projGlobal.nDet();
     mpi_nview = projGlobal.nView();
 
     mpi_ndet = projGlobal.nDet();
     mpi_nview = projGlobal.nView();


@@ -310,8 +313,11 @@ pjrec_main (int argc, char * argv[])
   imLocal = new ImageFile (nx, ny);
 #else
   projGlobal.read (pszFilenameProj);
   imLocal = new ImageFile (nx, ny);
 #else
   projGlobal.read (pszFilenameProj);

-  if (bOptVerbose)
-    projGlobal.printScanInfo();
+  if (bOptVerbose) {
+    ostringstream os;
+    projGlobal.printScanInfo(os);
+    cout << os.str();
+  }

 
   imGlobal = new ImageFile (nx, ny);
 #endif
 
   imGlobal = new ImageFile (nx, ny);
 #endif