+2.0.0-b10 - 8/31/00
+ ctsim: Added Auto Scaling for image windows
+ ctsim: Change title from "windowing" to "display scaling"
+ Added FieldOfView and FocalLength ratio parameters to projection collection
+ Changed name of Rowland Phantom to correct name of Shepp-Logan
+ Added data collection for equilinear scanner geometry
+ Fixed bug in backprojection selection
+ Improved projection collection animation
+
2.0.0-b9 - 8/22/00
Added RCS Id strings to executable files
Added RPM Spec file for RPM package creation
PACKAGE=ctsim
-VERSION=2.0.0b9
+VERSION=2.0.0b10
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; }
dnl CDPATH=
AC_INIT(src/ctsim.cpp)
-AM_INIT_AUTOMAKE(ctsim,2.0.0b9)
+AM_INIT_AUTOMAKE(ctsim,2.0.0b10)
AM_CONFIG_HEADER(config.h)
dnl Checks for programs.
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: backprojectors.h,v 1.12 2000/07/22 15:45:33 kevin Exp $
+** $Id: backprojectors.h,v 1.13 2000/08/25 15:59:13 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
{
public:
BackprojectTrig (const Projections& proj, ImageFile& im, int interpID, const int interpFactor)
- : Backproject::Backproject (proj, im, interpType, interpFactor)
+ : Backproject::Backproject (proj, im, interpID, interpFactor)
{}
void BackprojectView (const double* const t, double view_angle);
{
public:
BackprojectDiff2 (const Projections& proj, ImageFile& im, int interpID, const int interpFactor)
- : BackprojectDiff::BackprojectDiff (proj, im, interpType, interpFactor)
+ : BackprojectDiff::BackprojectDiff (proj, im, interpID, interpFactor)
{}
void BackprojectView (const double* const t, double view_angle);
{
public:
BackprojectIntDiff2 (const Projections& proj, ImageFile& im, int interpID, const int interpFactor)
- : BackprojectDiff::BackprojectDiff (proj, im, interpType, interpFactor)
+ : BackprojectDiff::BackprojectDiff (proj, im, interpID, interpFactor)
{}
void BackprojectView (const double* const t, double view_angle);
class BackprojectIntDiff3 : public BackprojectDiff
{
public:
- BackprojectIntDiff3 (const Projections& proj, ImageFile& im, int interpType, const int interpFactor)
- : BackprojectDiff::BackprojectDiff (proj, im, interpType, interpFactor)
+ BackprojectIntDiff3 (const Projections& proj, ImageFile& im, int interpID, const int interpFactor)
+ : BackprojectDiff::BackprojectDiff (proj, im, interpID, interpFactor)
{}
void BackprojectView (const double* const t, double view_angle);
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: phantom.h,v 1.12 2000/07/31 14:48:35 kevin Exp $
+** $Id: phantom.h,v 1.13 2000/08/25 15:59:13 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
public:
static const int PHM_INVALID;
static const int PHM_HERMAN;
- static const int PHM_BHERMAN;
- static const int PHM_ROWLAND;
- static const int PHM_BROWLAND;
+ static const int PHM_B_HERMAN;
+ static const int PHM_SHEPP_LOGAN;
+ static const int PHM_B_SHEPP_LOGAN;
static const int PHM_UNITPULSE;
Phantom ();
void addStdHerman ();
void addStdHermanBordered ();
- void addStdRowland ();
- void addStdRowlandBordered ();
+ void addStdSheppLogan ();
+ void addStdSheppLoganBordered ();
void print () const;
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: scanner.h,v 1.11 2000/08/02 18:06:00 kevin Exp $
+** $Id: scanner.h,v 1.12 2000/08/25 15:59:13 kevin Exp $
**
**
** This program is free software; you can redistribute it and/or modify
static const int Scanner::GEOMETRY_EQUIANGULAR;
- Scanner (const Phantom& phm, const char* const geometryName, int nDet, int nView, int nSample, const double rot_anglen);
+ Scanner (const Phantom& phm, const char* const geometryName, int nDet, int nView, int nSample, const double rot_anglen, double dFieldOfView, double dFocalLength);
~Scanner();
- void collectProjections (Projections& proj, const Phantom& phm, const int start_view = 0, const int trace = TRACE_NONE, SGP* pSGP = NULL);
+ void collectProjections (Projections& proj, const Phantom& phm, const int trace = TRACE_NONE, SGP* pSGP = NULL);
+
+ void collectProjections (Projections& proj, const Phantom& phm, const int iStartView, const int iNumViews, bool bStoreAtViewPosition, const int trace = TRACE_NONE, SGP* pSGP = NULL);
void setNView (int nView);
const double phmLen(void) const {return m_phmLen;}
const double rotInc(void) const {return m_rotInc;}
const double detInc(void) const {return m_detInc;}
- const double radius(void) const {return m_radius;}
+ const double detLen(void) const {return m_detLen;}
static const int getGeometryCount() {return s_iGeometryCount;}
static const char** getGeometryNameArray() {return s_aszGeometryName;}
unsigned int m_nDet; /* Number of detectors in array */
unsigned int m_nView; /* Number of rotated views */
unsigned int m_nSample; /* Number of rays per detector */
- double m_detLen; /* Total length of detector array */
- double m_rotLen; /* Rotation angle length in radians (norm 2PI) */
- double m_detInc; /* Increment between centers of detectors */
- double m_rotInc; /* Increment in rotation angle between views */
- double m_radius; /* Radius of rotation. Distance from center of phm to center of det */
- double m_phmLen; /* Maximum Length of phantom or area of interest */
+ double m_dFieldOfView; // Field of View
+ double m_dFocalLength; // Focal Length
+ double m_detLen; // Total length of detector array
+ double m_rotLen; // Rotation angle length in radians (norm 2PI)
+ double m_detInc; // Increment between centers of detectors
+ double m_rotInc; // Increment in rotation angle between views
+ double m_phmLen; // Maximum Length of phantom or area of interest
int m_trace;
struct {
double xd1,yd1,xd2,yd2; /* Coordinates of detector endpoints */
static const char* s_aszGeometryName[];
static const char* s_aszGeometryTitle[];
static const int s_iGeometryCount;
- static const int N_EXTRA_DETECTORS=0; // Number of extra detectors widths when calculating detlen
void projectSingleView (const Phantom& phm, DetectorArray& darray, const double xd1, const double yd1, const double xd2, const double yd2, const double xs1, const double ys1, const double xs2, const double ys2, SGP* pSGP);
double projectLineAgainstPElem (const PhantomElement& pelem, const double x1, const double y1, const double x2, const double y2, SGP* pSGP);
- void traceShowParam (const char *label, const char *fmt, int row, int color, ...);
+ void traceShowParam (SGP* pSGP, const char *label, const char *fmt, int row, int color, ...);
};
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: sgp.h,v 1.13 2000/08/02 18:06:00 kevin Exp $
+** $Id: sgp.h,v 1.14 2000/08/25 15:59:13 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
static RGBColor s_aRGBColor[];
static int s_iRGBColorCount;
+#if HAVE_WXWINDOWS
+ wxPen* m_pPen;
+#endif
+
public:
enum { // linestyles
LS_NOLINE = 0,
#ifndef TRACE_H
#define TRACE_H
-enum TraceID {
- TRACE_INVALID=-1,
- TRACE_NONE=0, /* No tracing */
- TRACE_TEXT, /* Minimal status */
- TRACE_PHM, /* Show phantom */
- TRACE_RAYS, /* Show all rays */
- TRACE_PLOT, /* Plot raysums */
- TRACE_CLIPPING /* Plot clipping */
-};
+static const int TRACE_INVALID = 0xFFFF;
+static const int TRACE_NONE = 0x0000;
+static const int TRACE_TEXT = 0x0001;
+static const int TRACE_PHM = 0x0002;
+static const int TRACE_RAYS = 0x0004;
+static const int TRACE_PLOT = 0x0008;
+static const int TRACE_CLIPPING = 0x0010;
class TraceLevel
{
int getTraceLevel(void) const { return m_traceLevel; }
- static TraceID convertTraceNameToID (const char* traceName);
+ static int convertTraceNameToID (const char* traceName);
private:
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: sgp.cpp,v 1.9 2000/08/02 18:06:00 kevin Exp $
+** $Id: sgp.cpp,v 1.10 2000/08/25 15:59:13 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
SGP::SGP (const SGPDriver& driver)
: m_driver (driver)
+#if HAVE_WXWINDOWS
+ , m_pPen(NULL)
+#endif
{
m_iPhysicalXSize = m_driver.getPhysicalXSize();
m_iPhysicalYSize = m_driver.getPhysicalYSize();
#if HAVE_WXWINDOWS
if (m_driver.isWX()) {
wxColor colour (s_aRGBColor[icol].getRed(), s_aRGBColor[icol].getGreen(), s_aRGBColor[icol].getBlue());
- wxPen pen (colour, 1, wxSOLID);
- m_driver.idWX()->SetPen (pen);
+ delete m_pPen;
+ m_pPen = new wxPen (colour, 1, wxSOLID);
+ m_driver.idWX()->SetPen (*m_pPen);
}
#endif
}
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: backprojectors.cpp,v 1.11 2000/07/23 01:49:03 kevin Exp $
+** $Id: backprojectors.cpp,v 1.12 2000/08/25 15:59:13 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
void Backproject::errorIndexOutsideDetector (int ix, int iy, double theta, double r, double phi, double L, int iDetPos)
{
- printf ("r=%f, phi=%f\n", r, phi);
+ sys_error (ERR_WARNING, "r=%f, phi=%f", r, phi);
errorIndexOutsideDetector (ix, iy, theta, L, iDetPos);
}
void Backproject::errorIndexOutsideDetector (int ix, int iy, double theta, double L, int iDetPos)
{
- printf ("ix=%d, iy=%d\n", ix, iy);
- printf ("theta=%f, L=%f, detInc=%f\n", theta, L, detInc);
- printf ("proj.ndet=%d, proj.detInc=%.4f, iDetCenter=%d\n", nDet, detInc, iDetCenter);
- printf ("xMin=%15.8f, xMax=%15.8f, xInc=%15.8f\n", xMin, xMax, xInc);
- printf ("yMin=%15.8f, yMax=%15.8f, yInc=%15.8f\n", yMin, yMax, yInc);
- sys_error (ERR_WARNING, "iDetPos index outside bounds: %d [backprojector]", iDetPos);
+ ostringstream os;
+ os << "ix=" << ix << ", iy=" << iy << ", theta=" << theta << ", L=" << L << ", detinc=" << detInc << "\n";
+ os << "ndet=" << nDet << ", detInc=" << detInc << ", iDetCenter=" << iDetCenter << "\n";
+ os << "xMin=" << xMin << ", xMax=" << xMax << ", xInc=" << xInc << "\n";
+ os << "yMin=" << yMin << ", yMax=" << yMax << ", yInc=" << yInc << "\n";
+ os << "iDetPos index outside bounds: " << iDetPos << " [backprojector]";;
+
+ sys_error (ERR_WARNING, os.str().c_str());
}
int detPosNearest = (curDetPos >= 0 ? ((curDetPos + halfScale) / scale) : ((curDetPos - halfScale) / scale));
int iDetPos = iDetCenter + detPosNearest; // calc index in the filtered raysum vector
- if (iDetPos < 0 || iDetPos >= nDet) // check for impossible: index outside of raysum pos
+ if (iDetPos < 0 || iDetPos >= nDet) // check for index outside of raysum pos
errorIndexOutsideDetector (ix, iy, theta, curDetPos, iDetPos);
else
*pImCol++ += filteredProj[iDetPos];
kint32 detPosColStart = nearest<kint32> ((start_r * cos (theta - start_phi) / detInc + iDetCenter) * scale);
// precalculate scaled difference for linear interpolation
- double deltaFilteredProj [nDet - 1];
+ double deltaFilteredProj [nDet];
if (interpType == Backprojector::INTERP_LINEAR) {
for (int i = 0; i < nDet - 1; i++)
deltaFilteredProj[i] = (filteredProj[i+1] - filteredProj[i]) * dInvScale;
}
+ deltaFilteredProj[nDet - 1] = 0; // last detector
+ int iLastDet = nDet - 1;
for (int ix = 0; ix < nx; ix++, detPosColStart += det_dx) {
kint32 curDetPos = detPosColStart;
ImageFileColumn pImCol = v[ix];
if (interpType == Backprojector::INTERP_NEAREST) {
for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
const int iDetPos = (curDetPos + halfScale) >> 16;
- assert(iDetPos >= 0 && iDetPos < nDet);
- *pImCol++ += filteredProj[iDetPos];
+ if (iDetPos >= 0 && iDetPos <= iLastDet)
+ *pImCol++ += filteredProj[iDetPos];
} // end for iy
} else if (interpType == Backprojector::INTERP_FREQ_PREINTERPOLATION) {
for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
const int iDetPos = ((curDetPos + halfScale) >> 16) * m_interpFactor;
- assert(iDetPos >= 0 && iDetPos < nDet);
+ if (iDetPos >= 0 && iDetPos <= iLastDet)
*pImCol++ += filteredProj[iDetPos];
} // end for iy
} else if (interpType == Backprojector::INTERP_LINEAR) {
for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
const kint32 iDetPos = curDetPos >> scaleShift;
const kint32 detRemainder = curDetPos & scaleBitmask;
- assert(iDetPos >= 0 && iDetPos < nDet - 1);
- *pImCol++ += filteredProj[iDetPos] + (detRemainder * deltaFilteredProj[iDetPos]);
+ if (iDetPos >= 0 && iDetPos <= iLastDet)
+ *pImCol++ += filteredProj[iDetPos] + (detRemainder * deltaFilteredProj[iDetPos]);
} // end for iy
} //end linear
} // end for ix
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: imagefile.cpp,v 1.10 2000/08/22 16:49:56 kevin Exp $
+** $Id: imagefile.cpp,v 1.11 2000/08/25 15:59:13 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
return 0;
#if HAVE_G2_H
- int pens [nx * ny * scale * scale ];
+ int* pPens = new int [nx * ny * scale * scale ];
double view_scale = 255 / (pmax - pmin);
int id_X11 = g2_open_X11 (nx * scale, ny * scale);
grayscale[i] = g2_ink (id_X11, cval, cval, cval);
}
- for (int i= 0, iy = ny - 1; iy >= 0; iy--) {
- for (int ix = 0; ix < nx; ix++) {
- int cval = static_cast<int>((v[ix][iy] - pmin) * view_scale);
- if (cval < 0)
- cval = 0;
- else if (cval > 255)
- cval = 255;
- pens[i++] = grayscale[cval];
+ for (int iy = ny - 1; iy >= 0; iy--) {
+ int iRowPos = ((ny - 1 - iy) * scale) * (nx * scale);
+ for (int ix = 0; ix < nx; ix++) {
+ int cval = static_cast<int>((v[ix][iy] - pmin) * view_scale);
+ if (cval < 0)
+ cval = 0;
+ else if (cval > 255)
+ cval = 255;
+ for (int sy = 0; sy < scale; sy++)
+ for (int sx = 0; sx < scale; sx++)
+ pPens[iRowPos+(sy * nx * scale)+(sx + (ix * scale))] = grayscale[cval];
}
}
- g2_image (id_X11, 0., 0., nx, ny, pens);
+ g2_image (id_X11, 0., 0., nx * scale, ny * scale, pPens);
+ delete pPens;
return (id_X11);
#endif
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: phantom.cpp,v 1.15 2000/08/02 18:06:00 kevin Exp $
+** $Id: phantom.cpp,v 1.16 2000/08/25 15:59:13 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
const int Phantom::PHM_INVALID = -1;
const int Phantom::PHM_HERMAN = 0;
-const int Phantom::PHM_BHERMAN = 1;
-const int Phantom::PHM_ROWLAND = 2;
-const int Phantom::PHM_BROWLAND = 3;
+const int Phantom::PHM_B_HERMAN = 1;
+const int Phantom::PHM_SHEPP_LOGAN = 2;
+const int Phantom::PHM_B_SHEPP_LOGAN = 3;
const int Phantom::PHM_UNITPULSE = 4;
const char* Phantom::s_aszPhantomName[] =
{
{"herman"},
- {"bherman"},
- {"rowland"},
- {"browland"},
+ {"herman-b"},
+ {"shepp-logan"},
+ {"shepp-logan-b"},
{"unitpulse"},
};
const char* Phantom::s_aszPhantomTitle[] =
{
{"Herman Head"},
- {"Herman Head Bordered"},
- {"Rowland Head"},
- {"Rowland Head Bordered"},
+ {"Herman Head (Bordered)"},
+ {"Shepp-Logan"},
+ {"Shepp-Logan (Bordered)"},
{"Unit Pulse"},
};
case PHM_HERMAN:
addStdHerman();
break;
- case PHM_BHERMAN:
+ case PHM_B_HERMAN:
addStdHermanBordered();
break;
- case PHM_ROWLAND:
- addStdRowland();
+ case PHM_SHEPP_LOGAN:
+ addStdSheppLogan();
break;
- case PHM_BROWLAND:
- addStdRowlandBordered ();
+ case PHM_B_SHEPP_LOGAN:
+ addStdSheppLoganBordered();
break;
case PHM_UNITPULSE:
m_composition = P_UNIT_PULSE;
/* NAME
- * addStdRowland Make head phantom of S.W. Rowland
+ * addStdSheppLogan Make head phantom of Shepp-Logan
*
* REFERENCES
* S. W. Rowland, "Computer Implementation of Image Reconstruction
*/
void
-Phantom::addStdRowland ()
+Phantom::addStdSheppLogan ()
{
addPElem ("ellipse", 0.0000, 0.0000, 0.6900, 0.9200, 0.0, 1.00);
addPElem ("ellipse", 0.0000, -0.0184, 0.6624, 0.8740, 0.0, -0.98);
}
void
-Phantom::addStdRowlandBordered ()
+Phantom::addStdSheppLoganBordered ()
{
- addStdRowland ();
+ addStdSheppLogan ();
addPElem ("rectangle", 0.000, 0.0000, 0.7500, 1.000, 0.0, 0.00);
}
Phantom::addStdHermanBordered ()
{
addStdHerman();
- addPElem ("rectangle", 0.000, 0.000, 8.650, 8.650, 0.00, 0.000);
+ addPElem ("rectangle", 0.000, 0.00, 9.000, 6.800, 0.00, 0.000);
}
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: procsignal.cpp,v 1.2 2000/08/22 07:02:48 kevin Exp $
+** $Id: procsignal.cpp,v 1.3 2000/08/25 15:59:13 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
// 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 = 0, int iPreinterpolationFactor = 1, int iTraceLevel = TRACE_NONE)
+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)
: m_adFourierCosTable(NULL), m_adFourierSinTable(NULL), m_adFilter(NULL), m_fail(false)
{
m_idFilterMethod = convertFilterMethodNameToID (szFilterMethodName);
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: projections.cpp,v 1.21 2000/08/22 07:02:48 kevin Exp $
+** $Id: projections.cpp,v 1.22 2000/08/25 15:59:13 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
m_rotInc = scanner.rotInc();
m_detInc = scanner.detInc();
m_rotStart = 0;
- m_detStart = -scanner.radius() + (scanner.detInc() / 2);
- m_phmLen = scanner.phmLen();
+ m_detStart = -(scanner.detLen() / 2) + (scanner.detInc() / 2);
}
void
processSignal.filterSignal (detval, filteredProj);
-
-
#ifdef HAVE_BSPLINE_INTERP
if (interp_type == I_BSPLINE)
bspline (m_nDet, zoom_factor, spline_order, filteredProj, filteredProj);
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: scanner.cpp,v 1.10 2000/08/03 09:21:12 kevin Exp $
+** $Id: scanner.cpp,v 1.11 2000/08/25 15:59:13 kevin Exp $
**
** This program is free software; you can redistribute it and/or modify
** it under the terms of the GNU General Public License (version 2) as
* int nSample Number of rays per detector
*/
-Scanner::Scanner (const Phantom& phm, const char* const geometryName, int nDet, int nView, int nSample, const double rot_anglen)
+Scanner::Scanner (const Phantom& phm, const char* const geometryName, int nDet, int nView, int nSample, const double rot_anglen, const double dFocalLengthRatio, const double dFieldOfViewRatio)
{
m_phmLen = phm.maxAxisLength(); // maximal length along an axis
m_nSample = 1;
if (nDet < 1)
nDet = 1;
- if ((nDet % 2) == 0)
- ++nDet; // ensure odd number of detectors
+ // if ((nDet % 2) == 0)
+ // ++nDet; // ensure odd number of detectors
m_nDet = nDet;
m_nView = nView;
m_nSample = nSample;
- m_detLen = SQRT2 * m_phmLen * ((m_nDet + N_EXTRA_DETECTORS) / static_cast<double>(m_nDet));
-
- m_rotLen = rot_anglen;
-
- m_radius = m_detLen / 2;
- m_detInc = m_detLen / m_nDet;
- m_rotInc = m_rotLen / m_nView;
-
- m_initPos.xd1 = m_detLen / 2;
- m_initPos.yd1 = -m_detLen / 2;
- m_initPos.xd2 = m_detLen / 2;
- m_initPos.yd2 = m_detLen / 2;
- m_initPos.xs1 = -m_detLen / 2;
- m_initPos.ys1 = -m_detLen / 2;
- m_initPos.xs2 = -m_detLen / 2;
- m_initPos.ys2 = m_detLen / 2;
- m_initPos.angle = 0.0;
+ m_dFocalLength = (m_phmLen * SQRT2 / 2) * dFocalLengthRatio;
+ m_dFieldOfView = m_phmLen * SQRT2 * dFieldOfViewRatio;
+
+ 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_dFocalLength;
+ m_initPos.ys1 = -dHalfDetLen;
+ m_initPos.xs2 = -m_dFocalLength;
+ m_initPos.ys2 = dHalfDetLen;
+ m_initPos.xd1 = m_dFocalLength;
+ m_initPos.yd1 = -dHalfDetLen;
+ m_initPos.xd2 = m_dFocalLength;
+ m_initPos.yd2 = dHalfDetLen;
+ m_initPos.angle = 0.0;
+ } else if (m_idGeometry == GEOMETRY_EQUILINEAR) {
+ double dHalfSquare = m_phmLen / 2;
+ double dFocalPastPhm = m_dFocalLength - dHalfSquare;
+ if (dFocalPastPhm <= 0.) {
+ m_fail = true;
+ m_failMessage = "Focal Point inside of phantom";
+ return;
+ }
+ double dAngle = atan( dHalfSquare / dFocalPastPhm );
+ 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_dFocalLength;
+ m_initPos.ys1 = 0;
+ m_initPos.xs2 = -m_dFocalLength;
+ m_initPos.ys2 = 0;
+ m_initPos.xd1 = m_dFocalLength;
+ m_initPos.yd1 = -dHalfDetLen;
+ m_initPos.xd2 = m_dFocalLength;
+ m_initPos.yd2 = dHalfDetLen;
+ m_initPos.angle = 0.0;
+ } else if (m_idGeometry == GEOMETRY_EQUIANGULAR) {
+ double dHalfSquare = m_phmLen / 2;
+ double dFocalPastPhm = m_dFocalLength - dHalfSquare;
+ if (dFocalPastPhm <= 0.) {
+ m_fail = true;
+ m_failMessage = "Focal Point inside of phantom";
+ return;
+ }
+ double dAngle = atan( dHalfSquare / dFocalPastPhm );
+ m_detLen = dAngle * 2;
+ }
}
Scanner::~Scanner (void)
/* NAME
- * raysum_collect Calculate ray sums for a Phantom
+ * collectProjections Calculate projections for a Phantom
*
* SYNOPSIS
- * rs = raysum_collect (det, phm, start_view, trace, unit_pulse)
- * Scanner& det Scanner specifications**
- * RAYSUM *rs Calculated ray sum matrix
- * Phantom& phm Phantom we are taking ray sums of
- * int trace Boolean flag to signal ray sum tracing
+ * collectProjections (proj, phm, start_view, nView, bStoreViewPos, trace)
+ * Projectrions& proj Projection storage
+ * Phantom& phm Phantom for which we collect projections
+ * bool bStoreViewPos TRUE then storage proj at normal view position
+ * int trace Trace level
*/
+
void
-Scanner::collectProjections (Projections& proj, const Phantom& phm, const int start_view = 0, const int trace = TRACE_NONE, SGP* pSGP = NULL)
+Scanner::collectProjections (Projections& proj, const Phantom& phm, const int trace, SGP* pSGP)
+{
+ collectProjections (proj, phm, 0, proj.nView(), true, trace, pSGP);
+}
+
+void
+Scanner::collectProjections (Projections& proj, const Phantom& phm, const int iStartView, const int iNumViews, bool bStoreAtViewPosition, const int trace, SGP* pSGP)
{
GRFMTX_2D rotmtx_initial, temp;
GRFMTX_2D rotmtx_incr;
- double start_angle = start_view * proj.rotInc();
+ double start_angle = iStartView * proj.rotInc();
double xcent = phm.xmin() + (phm.xmax() - phm.xmin()) / 2;
double ycent = phm.ymin() + (phm.ymax() - phm.ymin()) / 2;
m_trace = trace;
-#ifdef HAVE_SGP
- if (pSGP && m_trace >= TRACE_PHM) {
- double halfPhmLen = m_phmLen / 2;
- double wsize = SQRT2 * halfPhmLen;
-
- pSGP->setRasterOp (RO_SET);
- pSGP->eraseWindow ();
- pSGP->setColor (C_LTBLUE);
- pSGP->setWindow (xcent - wsize, ycent - wsize, xcent + wsize, ycent + wsize);
- pSGP->setColor (C_BROWN);
- pSGP->drawRect (xcent - halfPhmLen, ycent - halfPhmLen, xcent + halfPhmLen, ycent + halfPhmLen);
- pSGP->setColor (C_BROWN);
- pSGP->moveAbs (0., 0.);
- pSGP->drawCircle (wsize);
-
- traceShowParam ("X-Ray Simulator", "%s", RAYSUM_TRACE_ROW_TITLE, C_BLACK, " ");
- traceShowParam ("---------------", "%s", RAYSUM_TRACE_ROW_TITLE2, C_BLACK, " ");
- traceShowParam ("Phantom:", "%s", RAYSUM_TRACE_ROW_PHANT_ID, C_YELLOW, " Herman");
- traceShowParam ("Chomaticity :", "%s", RAYSUM_TRACE_ROW_CHROMATIC, C_LTGREEN, "Mono");
- traceShowParam ("Scatter :", "%5.1f", RAYSUM_TRACE_ROW_SCATTER, C_LTGREEN, 0.);
- traceShowParam ("Photon Uncert:", "%5.1f", RAYSUM_TRACE_ROW_PHOT_STAT, C_LTGREEN, 0.);
- traceShowParam ("Num Detectors:", "%5d", RAYSUM_TRACE_ROW_NDET, C_LTRED, proj.nDet());
- traceShowParam ("Num Views :", "%5d", RAYSUM_TRACE_ROW_NVIEW, C_LTRED, proj.nView());
- traceShowParam ("Samples / Ray:", "%5d", RAYSUM_TRACE_ROW_SAMPLES, C_LTRED, m_nSample);
-
- pSGP->setColor (C_LTGREEN);
- phm.draw (*pSGP);
-
- pSGP->setMarker (SGP::MARK_BDIAMOND, C_LTGREEN);
- }
-#endif
-
/* Calculate initial rotation matrix */
xlat_mtx2 (rotmtx_initial, -xcent, -ycent);
rot_mtx2 (temp, start_angle);
int iview;
double viewAngle;
- for (iview = 0, viewAngle = start_angle; iview < proj.nView(); iview++, viewAngle += proj.rotInc()) {
- DetectorArray& detArray = proj.getDetectorArray( iview );
+ for (iview = 0, viewAngle = start_angle; iview < iNumViews; iview++, viewAngle += proj.rotInc()) {
+ int iStoragePosition = iview;
+ if (bStoreAtViewPosition)
+ iStoragePosition += iStartView;
+
+ DetectorArray& detArray = proj.getDetectorArray( iStoragePosition );
+
+#ifdef HAVE_SGP
+ if (pSGP && m_trace >= TRACE_PHM) {
+ pSGP->eraseWindow();
+ double dWindowSize = max(m_detLen, m_dFocalLength * 2) * SQRT2;
+ double dHalfWindowSize = dWindowSize / 2;
+ double dHalfPhmLen = m_phmLen / 2;
+
+ pSGP->setRasterOp (RO_SET);
+ pSGP->eraseWindow ();
+ pSGP->setWindow (xcent - dHalfWindowSize, ycent - dHalfWindowSize, xcent + dHalfWindowSize, ycent + dHalfWindowSize);
+ pSGP->setColor (C_BROWN);
+ pSGP->moveAbs (0., 0.);
+ pSGP->drawRect (xcent - dHalfPhmLen, ycent - dHalfPhmLen, xcent + dHalfPhmLen, ycent + dHalfPhmLen);
+
+#if 0
+ traceShowParam (pSGP, "X-Ray Simulator", "%s", RAYSUM_TRACE_ROW_TITLE, C_BLACK, " ");
+ traceShowParam (pSGP, "---------------", "%s", RAYSUM_TRACE_ROW_TITLE2, C_BLACK, " ");
+ traceShowParam (pSGP, "Phantom:", "%s", RAYSUM_TRACE_ROW_PHANT_ID, C_YELLOW, " Herman");
+ traceShowParam (pSGP, "Chomaticity :", "%s", RAYSUM_TRACE_ROW_CHROMATIC, C_LTGREEN, "Mono");
+ traceShowParam (pSGP, "Scatter :", "%5.1f", RAYSUM_TRACE_ROW_SCATTER, C_LTGREEN, 0.);
+ traceShowParam (pSGP, "Photon Uncert:", "%5.1f", RAYSUM_TRACE_ROW_PHOT_STAT, C_LTGREEN, 0.);
+ traceShowParam (pSGP, "Num Detectors:", "%5d", RAYSUM_TRACE_ROW_NDET, C_LTRED, proj.nDet());
+ traceShowParam (pSGP, "Num Views :", "%5d", RAYSUM_TRACE_ROW_NVIEW, C_LTRED, proj.nViews());
+ traceShowParam (pSGP, "Samples / Ray:", "%5d", RAYSUM_TRACE_ROW_SAMPLES, C_LTRED, m_nSample);
+#endif
+
+ pSGP->setColor (C_RED);
+ phm.draw (*pSGP);
+
+ pSGP->setMarker (SGP::MARK_BDIAMOND, C_LTGREEN);
+ }
+#endif
#ifdef HAVE_SGP
if (pSGP && m_trace >= TRACE_PHM) {
pSGP->lineAbs (xd2, yd2);
pSGP->moveAbs (xs1, ys1);
pSGP->lineAbs (xs2, ys2);
+ pSGP->setRasterOp (RO_SET);
}
- if (m_trace)
- traceShowParam ("Current View :", "%5d", RAYSUM_TRACE_ROW_CURR_VIEW, C_LTMAGENTA, iview);
+ if (m_trace >= TRACE_TEXT)
+ traceShowParam (pSGP, "Current View :", "%5d", RAYSUM_TRACE_ROW_CURR_VIEW, C_LTMAGENTA, iview);
#endif
projectSingleView (phm, detArray, xd1, yd1, xd2, yd2, xs1, ys1, xs2, ys2, pSGP);
if (pSGP && m_trace >= TRACE_PHM) {
// rs_plot (detArray, xd1, yd1, xcent, ycent, theta);
pSGP->setColor (C_RED);
+ pSGP->setRasterOp (RO_XOR);
pSGP->moveAbs (xd1, yd1);
pSGP->lineAbs (xd2, yd2);
pSGP->moveAbs (xs1, ys1);
pSGP->lineAbs (xs2, ys2);
+ pSGP->setRasterOp (RO_SET);
}
#endif
xform_mtx2 (rotmtx_incr, xd1, yd1); // rotate detector endpoints
#ifdef HAVE_SGP
if (pSGP && m_trace >= TRACE_RAYS) {
pSGP->setColor (C_LTBLUE);
+ pSGP->setRasterOp (RO_XOR);
pSGP->moveAbs (xs, ys);
pSGP->lineAbs (xd, yd);
+ pSGP->setRasterOp (RO_SET);
}
#endif
sum += projectSingleLine (phm, xd, yd, xs, ys, pSGP);
#ifdef HAVE_SGP
- if (m_trace >= TRACE_RAYS)
- traceShowParam ("Attenuation :", "%5.2f", RAYSUM_TRACE_ROW_ATTEN, C_LTMAGENTA, "sum");
+ if (m_trace >= TRACE_RAYS) {
+ traceShowParam (pSGP, "Attenuation :", "%5.2f", RAYSUM_TRACE_ROW_ATTEN, C_LTMAGENTA, sum);
+ }
if (pSGP && m_trace >= TRACE_RAYS) {
pSGP->setColor (C_LTBLUE);
+ pSGP->setRasterOp (RO_XOR);
pSGP->moveAbs (xs, ys);
pSGP->lineAbs (xd, yd);
+ pSGP->setRasterOp (RO_SET);
}
#endif
xd += ddx2;
void
-Scanner::traceShowParam (const char *label, const char *fmt, int row, int color, ...)
+Scanner::traceShowParam (SGP* pSGP, const char *label, const char *fmt, int row, int color, ...)
{
- char s[80];
+ char s[256];
va_list arg;
va_start(arg, color);
- // cio_set_cpos (raysum_trace_menu_column, row);
snprintf (s, sizeof(s), label, "%s");
- // cio_set_text_clr (color - 8, 0);
- cio_put_str (s);
+ string strOut(s);
vsnprintf (s, sizeof(s), fmt, arg);
+ strOut += s;
+
+ // cio_set_cpos (raysum_trace_menu_column, row);
+ // cio_set_text_clr (color - 8, 0);
// cio_set_text_clr (color, 0);
- cio_put_str (s);
- cout << "\n";
+
+ if (pSGP) {
+ pSGP->moveAbs (0., row * 0.04);
+ pSGP->setTextColor (color, -1);
+ pSGP->drawText (strOut);
+ } else {
+ cio_put_str (strOut.c_str());
+ cio_put_str ("\n");
+ }
+
va_end(arg);
}
#ifdef HAVE_SGP
if (pSGP && m_trace == TRACE_CLIPPING) {
+ pSGP->setRasterOp (RO_XOR);
pSGP->moveAbs (x1, y1);
pSGP->lineAbs (x2, y2);
cio_tone (8000., 0.05);
pSGP->moveAbs (x1, y1);
pSGP->lineAbs (x2, y2);
+ pSGP->setRasterOp (RO_SET);
}
#endif
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: trace.cpp,v 1.1 2000/07/20 11:17:31 kevin Exp $
+** $Id: trace.cpp,v 1.2 2000/08/25 15:59:13 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
const char TraceLevel::TRACE_CLIPPING_STR[]= "clipping";
-TraceID
+int
TraceLevel::convertTraceNameToID (const char *traceString)
{
- TraceID traceID = TRACE_INVALID;
+ int traceID = TRACE_INVALID;
if (strcasecmp (traceString, TRACE_NONE_STR) == 0)
traceID = TRACE_NONE;
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: syserror.cpp,v 1.3 2000/06/19 19:04:05 kevin Exp $
+** $Id: syserror.cpp,v 1.4 2000/08/25 15:59:13 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
void sys_verror (int severity, const char *msg, va_list arg)
{
if (severity < errorlevel)
- return; /* ignore error if less than max level */
+ return; // ignore error if less than max level
nErrorCount++;
if (severity != ERR_FATAL) {
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: ctsim.h,v 1.4 2000/07/23 01:49:03 kevin Exp $
+** $Id: ctsim.h,v 1.5 2000/08/25 15:59:13 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
PJMENU_FILE_PROPERTIES,
PHMMENU_FILE_PROPERTIES,
PJMENU_PROCESS_RECONSTRUCT,
- IFMENU_VIEW_WINDOW_AUTO,
- IFMENU_VIEW_WINDOW_MINMAX,
+ IFMENU_VIEW_SCALE_AUTO,
+ IFMENU_VIEW_SCALE_MINMAX,
PHMMENU_PROCESS_RASTERIZE,
PHMMENU_PROCESS_PROJECTIONS,
};
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: dialogs.cpp,v 1.8 2000/08/22 07:02:48 kevin Exp $
+** $Id: dialogs.cpp,v 1.9 2000/08/25 15:59:13 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
/////////////////////////////////////////////////////////////////////
-DialogGetProjectionParameters::DialogGetProjectionParameters (wxFrame* pParent, int iDefaultNDet = 0, int iDefaultNView = 0, int iDefaultNSamples = 1, double dDefaultRotAngle = 1., int iDefaultGeometry = Scanner::GEOMETRY_PARALLEL)
- : wxDialog (pParent, -1, "Set Projection Parameters", wxDefaultPosition, wxDefaultSize, wxDEFAULT_DIALOG_STYLE | wxCAPTION)
+DialogGetProjectionParameters::DialogGetProjectionParameters (wxFrame* pParent, int iDefaultNDet = 0, int iDefaultNView = 0, int iDefaultNSamples = 1, double dDefaultRotAngle = 1., double dDefaultFocalLength = 1, double dDefaultFieldOfView = 1., int iDefaultGeometry = Scanner::GEOMETRY_PARALLEL)
+ : wxDialog (pParent, -1, "Set Projection Parameters", wxDefaultPosition, wxDefaultSize, wxDEFAULT_DIALOG_STYLE | wxCAPTION)
{
wxBoxSizer* pTopSizer = new wxBoxSizer (wxVERTICAL);
+ m_dDefaultRotAngle = dDefaultRotAngle;
+ m_dDefaultFocalLength = dDefaultFocalLength;
+ m_dDefaultFieldOfView = dDefaultFieldOfView;
+ m_iDefaultNSamples = iDefaultNSamples;
+ m_iDefaultNView = iDefaultNView;
+ m_iDefaultNDet = iDefaultNDet;
+
pTopSizer->Add (new wxStaticText (this, -1, "Set Projection Parameters"), 0, wxALIGN_CENTER | wxTOP | wxLEFT | wxRIGHT, 5);
pTopSizer->Add (new wxStaticLine (this, -1, wxDefaultPosition, wxSize(3,3), wxHORIZONTAL), 0, wxEXPAND | wxALL, 5);
ostringstream osRotAngle;
osRotAngle << dDefaultRotAngle;
m_pTextCtrlRotAngle = new wxTextCtrl (this, -1, osRotAngle.str().c_str(), wxDefaultPosition, wxSize(100, 25), 0);
+ ostringstream osFocalLength;
+ osFocalLength << dDefaultFocalLength;
+ m_pTextCtrlFocalLength = new wxTextCtrl (this, -1, osFocalLength.str().c_str(), wxDefaultPosition, wxSize(100, 25), 0);
+ ostringstream osFieldOfView;
+ osFieldOfView << dDefaultFieldOfView;
+ m_pTextCtrlFieldOfView = new wxTextCtrl (this, -1, osFieldOfView.str().c_str(), wxDefaultPosition, wxSize(100, 25), 0);
wxGridSizer *pGridSizer = new wxGridSizer (2, 4, 5);
pGridSizer->Add (new wxStaticText (this, -1, "Detectors"), 0, wxALIGN_RIGHT | wxALIGN_CENTER_VERTICAL);
pGridSizer->Add (m_pTextCtrlNSamples, 0, wxALIGN_CENTER_VERTICAL);
pGridSizer->Add (new wxStaticText (this, -1, "Rotation Angle (PI units)"), 0, wxALIGN_RIGHT | wxALIGN_CENTER_VERTICAL);
pGridSizer->Add (m_pTextCtrlRotAngle, 0, wxALIGN_CENTER_VERTICAL);
+ pGridSizer->Add (new wxStaticText (this, -1, "Focal Length Ratio (phantom radius units)"), 0, wxALIGN_RIGHT | wxALIGN_CENTER_VERTICAL);
+ pGridSizer->Add (m_pTextCtrlFocalLength, 0, wxALIGN_CENTER_VERTICAL);
+ pGridSizer->Add (new wxStaticText (this, -1, "Field of View (phantom diameter units)"), 0, wxALIGN_RIGHT | wxALIGN_CENTER_VERTICAL);
+ pGridSizer->Add (m_pTextCtrlFieldOfView, 0, wxALIGN_CENTER_VERTICAL);
pTopSizer->Add (pGridSizer, 1, wxALL, 10);
pTopSizer->Add (new wxStaticLine (this, -1, wxDefaultPosition, wxSize(3,3), wxHORIZONTAL), 0, wxEXPAND | wxALL, 5);
return (m_dDefaultRotAngle);
}
+double
+DialogGetProjectionParameters::getFocalLengthRatio (void)
+{
+ wxString strCtrl = m_pTextCtrlFocalLength->GetValue();
+ double dValue;
+ if (strCtrl.ToDouble (&dValue))
+ return (dValue);
+ else
+ return (m_dDefaultFocalLength);
+}
+
+double
+DialogGetProjectionParameters::getFieldOfViewRatio (void)
+{
+ wxString strCtrl = m_pTextCtrlFieldOfView->GetValue();
+ double dValue;
+ if (strCtrl.ToDouble (&dValue))
+ return (dValue);
+ else
+ return (m_dDefaultFieldOfView);
+}
+
const char*
DialogGetProjectionParameters::getGeometry (void)
{
{
return m_pListBoxFilterGeneration->getSelectionStringValue();
}
+
+
+DialogAutoScaleParameters::DialogAutoScaleParameters (wxFrame *pParent, const ImageFile& rIF)
+ : wxDialog (pParent, -1, "Auto Scale Parameters", wxDefaultPosition, wxDefaultSize, wxDEFAULT_DIALOG_STYLE | wxCAPTION), m_rImageFile(rIF)
+{
+ wxBoxSizer* pTopSizer = new wxBoxSizer (wxVERTICAL);
+
+ pTopSizer->Add (new wxStaticText (this, -1, "Auto Scale Parameters"), 0, wxALIGN_CENTER | wxTOP | wxLEFT | wxRIGHT, 5);
+
+ pTopSizer->Add (new wxStaticLine (this, -1, wxDefaultPosition, wxSize(3,3), wxHORIZONTAL), 0, wxEXPAND | wxALL, 5);
+
+ wxString asTitle[3];
+ asTitle[0] = "Median";
+ asTitle[1] = "Mode";
+ asTitle[2] = "Mean";
+
+ m_pListBoxCenter = new wxListBox (this, -1, wxDefaultPosition, wxDefaultSize, 3, asTitle, wxLB_SINGLE | wxLB_NEEDED_SB);
+ m_pListBoxCenter->SetSelection (0);
+ pTopSizer->Add (m_pListBoxCenter, 0, wxALL | wxALIGN_CENTER | wxEXPAND);
+
+ wxGridSizer *pGridSizer = new wxGridSizer (2, 2, 5);
+ 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);
+ pGridSizer->Add (m_pTextCtrlStdDevFactor, 0, wxALIGN_CENTER_VERTICAL);
+ pTopSizer->Add (pGridSizer, 1, wxALL, 10);
+
+ pTopSizer->Add (new wxStaticLine (this, -1, wxDefaultPosition, wxSize(3,3), wxHORIZONTAL), 0, wxEXPAND | wxALL, 5);
+
+ wxBoxSizer* pButtonSizer = new wxBoxSizer (wxHORIZONTAL);
+ wxButton* pButtonOk = new wxButton (this, wxID_OK, "Okay");
+ wxButton* pButtonCancel = new wxButton (this, wxID_CANCEL, "Cancel");
+ pButtonSizer->Add (pButtonOk, 0, wxEXPAND | wxALL, 10);
+ pButtonSizer->Add (pButtonCancel, 0, wxEXPAND | wxALL, 10);
+
+ pTopSizer->Add (pButtonSizer, 0, wxALIGN_CENTER);
+
+ SetAutoLayout (true);
+ SetSizer (pTopSizer);
+ pTopSizer->Fit (this);
+ pTopSizer->SetSizeHints (this);
+}
+
+void
+DialogAutoScaleParameters::getMinMax (double* pMin, double* pMax)
+{
+ int iCenter = m_pListBoxCenter->GetSelection();
+ double min, max, mean, mode, median, stddev;
+ m_rImageFile.statistics (min, max, mean, mode, median, stddev);
+ double dCenter = median;
+ if (iCenter == 1)
+ dCenter = mode;
+ else if (iCenter == 2)
+ dCenter = mean;
+
+ wxString sStddevFactor = m_pTextCtrlStdDevFactor->GetValue();
+ double dValue;
+ if (! sStddevFactor.ToDouble (&dValue)) {
+ *theApp->getLog() << "Error: Non-numeric Standard Deviation Factor of " << sStddevFactor << "\n";
+ *pMin = min;
+ *pMax = max;
+ }
+ double dHalfWidth = dValue * stddev / 2;
+ *pMin = dCenter - dHalfWidth;
+ *pMax = dCenter + dHalfWidth;
+ *theApp->getLog() << "Setting minimum to " << *pMin << " and maximum to " << *pMax << "\n";
+}
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: dialogs.h,v 1.9 2000/08/22 07:02:48 kevin Exp $
+** $Id: dialogs.h,v 1.10 2000/08/25 15:59:13 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
class DialogGetProjectionParameters : public wxDialog
{
public:
- DialogGetProjectionParameters (wxFrame* pParent, int iDefaultNDet = 0, int iDefaultNView = 0, int iDefaultNSamples = 1, double dDefaultRotAngle = 1., int iDefaultGeometry = Scanner::GEOMETRY_PARALLEL);
+ DialogGetProjectionParameters (wxFrame* pParent, int iDefaultNDet = 0, int iDefaultNView = 0, int iDefaultNSamples = 1, double dDefaultRotAngle = 1., double dDefaultFocalLength = 1, double dDefaultFieldOfView = 1., int iDefaultGeometry = Scanner::GEOMETRY_PARALLEL);
~DialogGetProjectionParameters (void);
unsigned int getNDet (void);
unsigned int getNView (void);
unsigned int getNSamples (void);
double getRotAngle (void);
+ double getFieldOfViewRatio (void);
+ double getFocalLengthRatio (void);
const char* getGeometry(void);
private:
wxTextCtrl* m_pTextCtrlNView;
wxTextCtrl* m_pTextCtrlNSamples;
wxTextCtrl* m_pTextCtrlRotAngle;
+ wxTextCtrl* m_pTextCtrlFocalLength;
+ wxTextCtrl* m_pTextCtrlFieldOfView;
StringValueAndTitleListBox* m_pListBoxGeometry;
int m_iDefaultNView;
int m_iDefaultNSamples;
double m_dDefaultRotAngle;
+ double m_dDefaultFocalLength;
+ double m_dDefaultFieldOfView;
};
int m_iDefaultInterpParam;
};
+class DialogAutoScaleParameters : public wxDialog
+{
+ public:
+ DialogAutoScaleParameters (wxFrame* pParent, const ImageFile& rImageFile);
+ virtual ~DialogAutoScaleParameters() {}
+
+ void getMinMax (double* pMin, double* pMax);
+
+ private:
+ const ImageFile& m_rImageFile;
+
+ wxTextCtrl* m_pTextCtrlStdDevFactor;
+ wxListBox* m_pListBoxCenter;
+};
+
#endif
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: views.cpp,v 1.15 2000/08/22 16:49:56 kevin Exp $
+** $Id: views.cpp,v 1.16 2000/08/25 15:59:13 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
BEGIN_EVENT_TABLE(ImageFileView, wxView)
EVT_MENU(IFMENU_FILE_PROPERTIES, ImageFileView::OnProperties)
- EVT_MENU(IFMENU_VIEW_WINDOW_MINMAX, ImageFileView::OnWindowMinMax)
- EVT_MENU(IFMENU_VIEW_WINDOW_AUTO, ImageFileView::OnWindowAuto)
+ EVT_MENU(IFMENU_VIEW_SCALE_MINMAX, ImageFileView::OnScaleMinMax)
+ EVT_MENU(IFMENU_VIEW_SCALE_AUTO, ImageFileView::OnScaleAuto)
END_EVENT_TABLE()
ImageFileView::ImageFileView(void)
}
void
-ImageFileView::OnWindowAuto (wxCommandEvent& event)
+ImageFileView::OnScaleAuto (wxCommandEvent& event)
{
- *theApp->getLog() << "ImageFile: Window Auto\n";
+ const ImageFile& rIF = GetDocument()->getImageFile();
+ DialogAutoScaleParameters dialogAutoScale (m_frame, rIF);
+ int iRetVal = dialogAutoScale.ShowModal();
+ if (iRetVal == wxID_OK) {
+ m_bMinSpecified = true;
+ m_bMaxSpecified = true;
+ double dMin, dMax;
+ dialogAutoScale.getMinMax (&dMin, &dMax);
+ m_dMinPixel = dMin;
+ m_dMaxPixel = dMax;
+ OnUpdate (this, NULL);
+ }
}
void
-ImageFileView::OnWindowMinMax (wxCommandEvent& event)
+ImageFileView::OnScaleMinMax (wxCommandEvent& event)
{
const ImageFile& rIF = GetDocument()->getImageFile();
double min, max;
file_menu->Append(wxID_PREVIEW, "Print Pre&view");
wxMenu *view_menu = new wxMenu;
- view_menu->Append(IFMENU_VIEW_WINDOW_MINMAX, "&Set Window Min/Max");
- view_menu->Append(IFMENU_VIEW_WINDOW_AUTO, "&Auto Window...");
+ view_menu->Append(IFMENU_VIEW_SCALE_MINMAX, "Display Scale &Set...");
+ view_menu->Append(IFMENU_VIEW_SCALE_AUTO, "Display Scale &Auto...");
wxMenu *help_menu = new wxMenu;
help_menu->Append(MAINMENU_HELP_ABOUT, "&About");
void
PhantomView::OnProjections (wxCommandEvent& event)
{
- DialogGetProjectionParameters dialogProjection (m_frame, 367, 320, 1, 1., Scanner::GEOMETRY_PARALLEL);
+ DialogGetProjectionParameters dialogProjection (m_frame, 367, 320, 1, 1., 1., 1., Scanner::GEOMETRY_PARALLEL);
int retVal = dialogProjection.ShowModal();
if (retVal == wxID_OK) {
int nDet = dialogProjection.getNDet();
int nView = dialogProjection.getNView();
int nSamples = dialogProjection.getNSamples();
double dRotAngle = dialogProjection.getRotAngle();
+ double dFocalLengthRatio = dialogProjection.getFocalLengthRatio();
+ double dFieldOfViewRatio = dialogProjection.getFieldOfViewRatio();
+
wxString sGeometry = dialogProjection.getGeometry();
if (nDet > 0 && nView > 0 && sGeometry != "") {
const Phantom& rPhantom = GetDocument()->getPhantom();
ProjectionFileDocument* pProjectionDoc = dynamic_cast<ProjectionFileDocument*>(theApp->getDocManager()->CreateDocument("untitled.pj", wxDOC_SILENT));
Projections& rProj = pProjectionDoc->getProjections();
- Scanner theScanner (rPhantom, sGeometry.c_str(), nDet, nView, nSamples, dRotAngle);
+ Scanner theScanner (rPhantom, sGeometry.c_str(), nDet, nView, nSamples, dRotAngle, dFocalLengthRatio, dFieldOfViewRatio);
+ if (theScanner.fail()) {
+ *theApp->getLog() << "Failed making scanner: " << theScanner.failMessage().c_str() << "\n";
+ return;
+ }
rProj.initFromScanner (theScanner);
#if 1
frame.GetClientSize(&x, &y);
SGPDriver driver (dynamic_cast<wxDC*>(&dc), x, y);
SGP sgp (driver);
- theScanner.collectProjections (rProj, rPhantom, 0, TRACE_PHM, &sgp);
+ for (int iView = 0; iView < rProj.nView(); iView++) {
+ theScanner.collectProjections (rProj, rPhantom, iView, 1, true, TRACE_RAYS, &sgp);
+ }
#else
theScanner.collectProjections (rProj, rPhantom);
#endif
pProjectionDoc->Modify(true);
pProjectionDoc->UpdateAllViews(this);
ostringstream os;
- os << "Projections for " << rPhantom.name() << ": nDet=" << nDet << ", nView=" << nView << ", nSamples=" << nSamples << ", RotAngle=" << dRotAngle << ", Geometry=" << sGeometry.c_str() << "\n";
+ os << "Projections for " << rPhantom.name() << ": nDet=" << nDet << ", nView=" << nView << ", nSamples=" << nSamples << ", RotAngle=" << dRotAngle << ", FocalLengthRatio=" << dFocalLengthRatio << ", FieldOfViewRatio=" << dFieldOfViewRatio << ", Geometry=" << sGeometry.c_str() << "\n";
*theApp->getLog() << os.str().c_str();
}
}
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: views.h,v 1.6 2000/07/28 08:28:08 kevin Exp $
+** $Id: views.h,v 1.7 2000/08/25 15:59:13 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
void OnUpdate(wxView *sender, wxObject *hint = NULL);
bool OnClose (bool deleteWindow = true);
void OnProperties (wxCommandEvent& event);
- void OnWindowAuto (wxCommandEvent& event);
- void OnWindowMinMax (wxCommandEvent& event);
+ void OnScaleAuto (wxCommandEvent& event);
+ void OnScaleMinMax (wxCommandEvent& event);
ImageFileDocument* GetDocument(void)
{ return dynamic_cast<ImageFileDocument*>(wxView::GetDocument()); }
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: if2img.cpp,v 1.5 2000/08/22 17:03:54 kevin Exp $
+** $Id: if2img.cpp,v 1.6 2000/08/25 15:59:13 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
{0, 0, 0, 0}
};
-static const char* g_szIdStr = "$Id: if2img.cpp,v 1.5 2000/08/22 17:03:54 kevin Exp $";
+static const char* g_szIdStr = "$Id: if2img.cpp,v 1.6 2000/08/25 15:59:13 kevin Exp $";
enum { O_AUTO_FULL, O_AUTO_STD0_1, O_AUTO_STD0_5, O_AUTO_STD1, O_AUTO_STD2, O_AUTO_STD3 };
static const char O_AUTO_FULL_STR[]="full";
}
}
+#if HAVE_SGP
+ if (optind + 1 == argc)
+ opt_format = O_FORMAT_DISP; // display image if no options
+#endif
+
if ((opt_format == O_FORMAT_DISP && optind + 1 != argc)
- || (opt_format != O_FORMAT_DISP && optind + 2 != argc))
- {
- if2img_usage(argv[0]);
- return (1);
- }
+ || (opt_format != O_FORMAT_DISP && optind + 2 != argc)) {
+ if2img_usage(argv[0]);
+ return (1);
+ }
in_file = argv[optind];
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: phm2if.cpp,v 1.9 2000/08/03 09:57:29 kevin Exp $
+** $Id: phm2if.cpp,v 1.10 2000/08/25 15:59:13 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
{0, 0, 0, 0}
};
-static const char* g_szIdStr = "$Id: phm2if.cpp,v 1.9 2000/08/03 09:57:29 kevin Exp $";
+static const char* g_szIdStr = "$Id: phm2if.cpp,v 1.10 2000/08/25 15:59:13 kevin Exp $";
void
phm2if_usage (const char *program)
cout << " ny Number of pixels Y-axis" << endl;
cout << " --phantom Phantom to use for projection" << endl;
cout << " herman Herman head phantom" << endl;
- cout << " bherman Bordered Herman head phantom" << endl;
- cout << " rowland Rowland head phantom" << endl;
- cout << " browland Bordered Rowland head phantom" << endl;
+ cout << " herman-b Herman head phantom (Bordered)" << endl;
+ cout << " shepp-logan Shepp-Logan head phantom" << endl;
+ cout << " shepp-logan-b Shepp-Logan head phantom (Bordered)" << endl;
cout << " unitpulse Unit pulse phantom" << endl;
cout << " --phmfile Generate Phantom from phantom file" << endl;
cout << " --filter Generate Phantom from a filter function" << endl;
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: phm2pj.cpp,v 1.9 2000/08/03 09:57:29 kevin Exp $
+** $Id: phm2pj.cpp,v 1.10 2000/08/25 15:59:13 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
#include "timer.h"
-enum { O_PHANTOM, O_DESC, O_NRAY, O_ROTANGLE, O_PHMFILE,
+enum { O_PHANTOM, O_DESC, O_NRAY, O_ROTANGLE, O_PHMFILE, O_GEOMETRY, O_FOCAL_LENGTH, O_FIELD_OF_VIEW,
O_TRACE, O_VERBOSE, O_HELP, O_DEBUG, O_VERSION };
static struct option phm2pj_options[] =
{"desc", 1, 0, O_DESC},
{"nray", 1, 0, O_NRAY},
{"rotangle", 1, 0, O_ROTANGLE},
+ {"geometry", 1, 0, O_GEOMETRY},
+ {"focal-length", 1, 0, O_FOCAL_LENGTH},
+ {"field-of-view", 1, 0, O_FIELD_OF_VIEW},
{"trace", 1, 0, O_TRACE},
{"verbose", 0, 0, O_VERBOSE},
{"help", 0, 0, O_HELP},
{0, 0, 0, 0}
};
-static const char* g_szIdStr = "$Id: phm2pj.cpp,v 1.9 2000/08/03 09:57:29 kevin Exp $";
+static const char* g_szIdStr = "$Id: phm2pj.cpp,v 1.10 2000/08/25 15:59:13 kevin Exp $";
void
phm2pj_usage (const char *program)
{
- cout << "usage: " << fileBasename(program) << " outfile ndet nview [--phantom phantom-name] [--phmfile filename] [OPTIONS]" << endl;
- cout << "Calculate (projections) through phantom object, either a predefined --phantom or a --phmfile" << endl;
- cout << "" << endl;
- cout << " outfile Name of output file for raysums" << endl;
- cout << " ndet Number of detectors" << endl;
- cout << " nview Number of rotated views" << endl;
- cout << " --phantom Phantom to use for projection" << endl;
- cout << " herman Herman head phantom" << endl;
- cout << " bherman Bordered herman head phantom" << endl;
- cout << " rowland Rowland head phantom" << endl;
- cout << " browland Bordered Rowland head phantom" << endl;
- cout << " unitpulse Unit pulse phantom" << endl;
- cout << " --phmfile Get Phantom from phantom file" << endl;
- cout << " --desc Description of raysum" << endl;
- cout << " --nray Number of rays per detector (default = 1)" << endl;
- cout << " --rotangle Degrees to rotate view through, multiple of PI (default = 1)" << endl;
- cout << " --trace Trace level to use" << endl;
- cout << " none No tracing (default)" << endl;
- cout << " text Trace text level" << endl;
- cout << " phm Trace phantom image" << endl;
- cout << " rays Trace rays" << endl;
- cout << " plot Trace plot" << endl;
- cout << " clipping Trace clipping" << endl;
- cout << " --verbose Verbose mode" << endl;
- cout << " --debug Debug mode" << endl;
- cout << " --version Print version" << endl;
- cout << " --help Print this help message" << endl;
+ cout << "usage: " << fileBasename(program) << " outfile ndet nview [--phantom phantom-name] [--phmfile filename] [OPTIONS]\n";
+ cout << "Calculate (projections) through phantom object, either a predefined --phantom or a --phmfile\n\n";
+ cout << " outfile Name of output file for raysums\n";
+ cout << " ndet Number of detectors\n";
+ cout << " nview Number of rotated views\n";
+ cout << " --phantom Phantom to use for projection\n";
+ cout << " herman Herman head phantom\n";
+ cout << " herman-b Herman head phantom (Bordered)\n";
+ cout << " shepp-logan Shepp-Logan head phantom\n";
+ cout << " shepp-logan-b Shepp-Logan head phantom (Bordered)\n";
+ cout << " unitpulse Unit pulse phantom\n";
+ cout << " --phmfile Get Phantom from phantom file\n";
+ cout << " --desc Description of raysum\n";
+ cout << " --nray Number of rays per detector (default = 1)\n";
+ cout << " --rotangle Degrees to rotate view through (multiple of PI)\n";
+ cout << " (default = 1)\n";
+ cout << " --focal-length Focal length ratio (ratio to radius of phantom)\n";
+ cout << " (default = 1)\n";
+ cout << " --field-of-view Field of view (ratio to diameter of phantom square)\n";
+ cout << " (default = 1)\n";
+ cout << " --trace Trace level to use\n";
+ cout << " none No tracing (default)\n";
+ cout << " text Trace text level\n";
+ cout << " phm Trace phantom image\n";
+ cout << " rays Trace rays\n";
+ cout << " plot Trace plot\n";
+ cout << " clipping Trace clipping\n";
+ cout << " --verbose Verbose mode\n";
+ cout << " --debug Debug mode\n";
+ cout << " --version Print version\n";
+ cout << " --help Print this help message\n";
}
#ifdef HAVE_MPI
int opt_ndet;
int opt_nview;
int opt_nray = 1;
+ double dOptFocalLength = 1.;
+ double dOptFieldOfView = 1.;
int opt_trace = 0;
string optPhmName (Phantom::convertPhantomIDToName(Phantom::PHM_HERMAN));
int opt_verbose = 0;
return (1);
}
break;
+ case O_GEOMETRY:
+ optGeometryName = optarg;
+ break;
+ case O_FOCAL_LENGTH:
+ dOptFocalLength = strtod(optarg, &endptr);
+ endstr = optarg + strlen(optarg);
+ if (endptr != endstr) {
+ cerr << "Error setting --focal-length to " << optarg << endl;
+ phm2pj_usage(argv[0]);
+ return (1);
+ }
+ break;
+ case O_FIELD_OF_VIEW:
+ dOptFieldOfView = strtod(optarg, &endptr);
+ endstr = optarg + strlen(optarg);
+ if (endptr != endstr) {
+ cerr << "Error setting --field-of-view to " << optarg << endl;
+ phm2pj_usage(argv[0]);
+ return (1);
+ }
+ break;
case O_NRAY:
opt_nray = strtol(optarg, &endptr, 10);
endstr = optarg + strlen(optarg);
#ifdef VERSION
cout << "Version: " << VERSION << endl << g_szIdStr << endl;
#else
- cout << "Unknown version number" << endl;
+ cout << "Unknown version number\n";
#endif
return (0);
case O_HELP:
}
if (optPhmName == "" && optPhmFileName == "") {
- cerr << "No phantom defined" << endl << endl;
+ cerr << "No phantom defined\n" << endl;
phm2pj_usage(argv[0]);
return (1);
}
if (optPhmFileName != "") {
#ifdef HAVE_MPI
- cerr << "Can not read phantom from file in MPI mode" << endl;
+ cerr << "Can not read phantom from file in MPI mode\n";
return (1);
#endif
phm.createFromFile (optPhmFileName.c_str());
#endif
opt_rotangle *= PI;
- Scanner scanner (phm, optGeometryName.c_str(), opt_ndet, opt_nview, opt_nray, opt_rotangle);
+ Scanner scanner (phm, optGeometryName.c_str(), opt_ndet, opt_nview, opt_nray, opt_rotangle, dOptFocalLength, dOptFieldOfView);
if (scanner.fail()) {
cout << "Scanner Creation Error: " << scanner.failMessage() << endl;
return (1);
pjLocal.setNView (mpiWorld.getMyLocalWorkUnits());
if (opt_debug)
- cout << "pjLocal->nview = " << pjLocal.nView() << " (process " << mpiWorld.getRank() << ")" << endl;;
+ cout << "pjLocal->nview = " << pjLocal.nView() << " (process " << mpiWorld.getRank() << ")\n";;
TimerCollectiveMPI timerProject (mpiWorld.getComm());
- scanner.collectProjections (pjLocal, phm, mpiWorld.getMyStartWorkUnit(), opt_trace);
+ scanner.collectProjections (pjLocal, phm, mpiWorld.getMyStartWorkUnit(), mpiWorld.getMyLocalWorkUnits(), false, opt_trace);
if (opt_verbose)
timerProject.timerEndAndReport ("Time to collect projections");
pSGPDriver = new SGPDriver ("phm2pj", 600, 600);
pSGP = new SGP (*pSGPDriver);
}
- scanner.collectProjections (pjGlobal, phm, 0, opt_trace, pSGP);
+ scanner.collectProjections (pjGlobal, phm, opt_trace, pSGP);
if (opt_trace >= TRACE_PHM) {
cout << "Press enter to continue\n";
cio_kb_getc();
}
#else
- scanner.collectProjections (pjGlobal, phm, 0, opt_trace);
+ scanner.collectProjections (pjGlobal, phm, opt_trace);
#endif
#endif
pjGlobal.printScanInfo();
cout << endl;
cout << " Remark: " << pjGlobal.remark() << endl;
- cout << "Run time: " << pjGlobal.calcTime() << " seconds" << endl;
+ cout << "Run time: " << pjGlobal.calcTime() << " seconds\n";
}
}
} catch (exception e) {
cerr << "Exception: " << e.what() << endl;
} catch (...) {
- cerr << "Unknown exception" << endl;
+ cerr << "Unknown exception\n";
}
return (retval);
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: pjrec.cpp,v 1.13 2000/08/22 07:02:48 kevin Exp $
+** $Id: pjrec.cpp,v 1.14 2000/08/25 15:59:13 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
{0, 0, 0, 0}
};
-static const char* g_szIdStr = "$Id: pjrec.cpp,v 1.13 2000/08/22 07:02:48 kevin Exp $";
+static const char* g_szIdStr = "$Id: pjrec.cpp,v 1.14 2000/08/25 15:59:13 kevin Exp $";
void
pjrec_usage (const char *program)
cout << " --zeropad n Set zeropad level (default = 0)\n";
cout << " set n to number of powers to two to pad\n";
cout << " --filter-generation Filter Generation mode\n";
- cout << " direct Use direct filter in spatial or frequency domain\n";
+ cout << " direct Use direct filter in spatial or frequency domain (default)\n";
cout << " inverse_fourier Use inverse fourier transform of inverse filter\n";
cout << " --backproj Backprojection Method" << endl;
cout << " trig Trigometric functions at every point" << endl;
string sRemark;
bool bOptVerbose = false;
bool bOptDebug = 1;
- int iOptZeropad = 0;
+ int iOptZeropad = 1;
int optTrace = TRACE_NONE;
double dOptFilterParam = -1;
string sOptFilterName (SignalFilter::convertFilterIDToName (SignalFilter::FILTER_ABS_BANDLIMIT));
string sOptFilterMethodName (ProcessSignal::convertFilterMethodIDToName (ProcessSignal::FILTER_METHOD_CONVOLUTION));
- string sOptFilterGenerationName (ProcessSignal::convertFilterGenerationIDToName (ProcessSignal::FILTER_GENERATION_INVERSE_FOURIER));
+ string sOptFilterGenerationName (ProcessSignal::convertFilterGenerationIDToName (ProcessSignal::FILTER_GENERATION_DIRECT));
string sOptInterpName (Backprojector::convertInterpIDToName (Backprojector::INTERP_LINEAR));
string sOptBackprojectName (Backprojector::convertBackprojectIDToName (Backprojector::BPROJ_IDIFF3));
int iOptPreinterpolationFactor = 1;
projects / ctsim.git / commitdiff