projects
/
ctsim.git
/ blobdiff
commit
grep
author
committer
pickaxe
?
search:
re
summary
|
shortlog
|
log
|
commit
|
commitdiff
|
tree
raw
|
inline
| side by side
r616: *** empty log message ***
[ctsim.git]
/
libctsim
/
backprojectors.cpp
diff --git
a/libctsim/backprojectors.cpp
b/libctsim/backprojectors.cpp
index c7a8f27f47f7bd34b1523af9d852760ff967cf35..bfd594e8697ef0075b98ce6bce5b325261b4b082 100644
(file)
--- a/
libctsim/backprojectors.cpp
+++ b/
libctsim/backprojectors.cpp
@@
-8,7
+8,7
@@
** This is part of the CTSim program
** Copyright (c) 1983-2001 Kevin Rosenberg
**
** This is part of the CTSim program
** Copyright (c) 1983-2001 Kevin Rosenberg
**
-** $Id: backprojectors.cpp,v 1.
26 2001/02/11 04:56:37
kevin Exp $
+** $Id: backprojectors.cpp,v 1.
30 2001/03/07 19:02:44
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
@@
-32,7
+32,7
@@
const int Backprojector::BPROJ_TABLE = 1;
const int Backprojector::BPROJ_DIFF = 2;
const int Backprojector::BPROJ_IDIFF = 3;
const int Backprojector::BPROJ_DIFF = 2;
const int Backprojector::BPROJ_IDIFF = 3;
-const char* Backprojector::s_aszBackprojectName[] =
+const char*
const
Backprojector::s_aszBackprojectName[] =
{
{"trig"},
{"table"},
{
{"trig"},
{"table"},
@@
-40,7
+40,7
@@
const char* Backprojector::s_aszBackprojectName[] =
{"idiff"},
};
{"idiff"},
};
-const char* Backprojector::s_aszBackprojectTitle[] =
+const char*
const
Backprojector::s_aszBackprojectTitle[] =
{
{"Direct Trigometric"},
{"Trigometric Table"},
{
{"Direct Trigometric"},
{"Trigometric Table"},
@@
-62,7
+62,7
@@
const int Backprojector::INTERP_2BSPLINE = 6;
const int Backprojector::INTERP_3BSPLINE = 7;
#endif
const int Backprojector::INTERP_3BSPLINE = 7;
#endif
-const char* Backprojector::s_aszInterpName[] =
+const char*
const
Backprojector::s_aszInterpName[] =
{
{"nearest"},
{"linear"},
{
{"nearest"},
{"linear"},
@@
-78,7
+78,7
@@
const char* Backprojector::s_aszInterpName[] =
#endif
};
#endif
};
-const char* Backprojector::s_aszInterpTitle[] =
+const char*
const
Backprojector::s_aszInterpTitle[] =
{
{"Nearest"},
{"Linear"},
{
{"Nearest"},
{"Linear"},
@@
-113,6
+113,13
@@
Backprojector::BackprojectView (const double* const viewData, const double viewA
m_pBackprojectImplem->BackprojectView (viewData, viewAngle);
}
m_pBackprojectImplem->BackprojectView (viewData, viewAngle);
}
+void
+Backprojector::PostProcessing()
+{
+ if (m_pBackprojectImplem != NULL)
+ m_pBackprojectImplem->PostProcessing();
+}
+
Backprojector::~Backprojector ()
{
delete m_pBackprojectImplem;
Backprojector::~Backprojector ()
{
delete m_pBackprojectImplem;
@@
-254,7
+261,7
@@
Backprojector::convertInterpIDToTitle (const int interpID)
// Pure virtual base class for all backprojectors.
Backproject::Backproject (const Projections& proj, ImageFile& im, int interpType, const int interpFactor)
// Pure virtual base class for all backprojectors.
Backproject::Backproject (const Projections& proj, ImageFile& im, int interpType, const int interpFactor)
-: proj(proj), im(im), interpType(interpType), m_interpFactor(interpFactor)
+: proj(proj), im(im), interpType(interpType), m_interpFactor(interpFactor)
, m_bPostProcessingDone(false)
{
detInc = proj.detInc();
nDet = proj.nDet();
{
detInc = proj.detInc();
nDet = proj.nDet();
@@
-282,11
+289,18
@@
Backproject::Backproject (const Projections& proj, ImageFile& im, int interpType
yInc = (yMax - yMin) / ny;
m_dFocalLength = proj.focalLength();
yInc = (yMax - yMin) / ny;
m_dFocalLength = proj.focalLength();
+ m_dSourceDetectorLength = proj.sourceDetectorLength();
}
Backproject::~Backproject ()
{}
}
Backproject::~Backproject ()
{}
+void
+Backproject::PostProcessing()
+{
+ m_bPostProcessingDone = true;
+}
+
void
Backproject::ScaleImageByRotIncrement ()
{
void
Backproject::ScaleImageByRotIncrement ()
{
@@
-351,7
+365,7
@@
BackprojectTrig::BackprojectView (const double* const filteredProj, const double
double frac = p - pFloor; // fraction distance from det
if (iDetPos >= 0 && iDetPos < nDet - 1)
v[ix][iy] += rotScale * ((1-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1]);
double frac = p - pFloor; // fraction distance from det
if (iDetPos >= 0 && iDetPos < nDet - 1)
v[ix][iy] += rotScale * ((1-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1]);
- } else if (interpType = Backprojector::INTERP_CUBIC) {
+ } else if (interpType =
=
Backprojector::INTERP_CUBIC) {
double p = iDetCenter + (L / detInc); // position along detector
if (p >= 0 && p < nDet)
v[ix][iy] += rotScale * pCubicInterp->interpolate (p);
double p = iDetCenter + (L / detInc); // position along detector
if (p >= 0 && p < nDet)
v[ix][iy] += rotScale * pCubicInterp->interpolate (p);
@@
-389,7
+403,15
@@
BackprojectTable::BackprojectTable (const Projections& proj, ImageFile& im, int
BackprojectTable::~BackprojectTable ()
{
BackprojectTable::~BackprojectTable ()
{
- ScaleImageByRotIncrement();
+}
+
+void
+BackprojectTable::PostProcessing()
+{
+ if (! m_bPostProcessingDone) {
+ ScaleImageByRotIncrement();
+ m_bPostProcessingDone = true;
+ }
}
void
}
void
@@
-419,7
+441,7
@@
BackprojectTable::BackprojectView (const double* const filteredProj, const doubl
double frac = dPos - dPosFloor; // fraction distance from det
if (iDetPos >= 0 && iDetPos < nDet - 1)
pImCol[iy] += ((1-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1]);
double frac = dPos - dPosFloor; // fraction distance from det
if (iDetPos >= 0 && iDetPos < nDet - 1)
pImCol[iy] += ((1-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1]);
- } else if (interpType = Backprojector::INTERP_CUBIC) {
+ } else if (interpType =
=
Backprojector::INTERP_CUBIC) {
double p = iDetCenter + (L / detInc); // position along detector
if (p >= 0 && p < nDet)
pImCol[iy] += pCubicInterp->interpolate (p);
double p = iDetCenter + (L / detInc); // position along detector
if (p >= 0 && p < nDet)
pImCol[iy] += pCubicInterp->interpolate (p);
@@
-451,11
+473,18
@@
BackprojectDiff::BackprojectDiff (const Projections& proj, ImageFile& im, int in
im.arrayDataClear();
}
im.arrayDataClear();
}
-BackprojectDiff::~BackprojectDiff()
+BackprojectDiff::~BackprojectDiff
()
{
{
- ScaleImageByRotIncrement();
}
}
+void
+BackprojectDiff::PostProcessing()
+{
+ if (! m_bPostProcessingDone) {
+ ScaleImageByRotIncrement();
+ m_bPostProcessingDone = true;
+ }
+}
void
BackprojectDiff::BackprojectView (const double* const filteredProj, const double view_angle)
void
BackprojectDiff::BackprojectView (const double* const filteredProj, const double view_angle)
@@
-489,7
+518,7
@@
BackprojectDiff::BackprojectView (const double* const filteredProj, const double
double frac = curDetPos - detPosFloor; // fraction distance from det
if (iDetPos > 0 && iDetPos < nDet - 1)
*pImCol++ += filteredProj[iDetPos] + (frac * (filteredProj[iDetPos+1] - filteredProj[iDetPos]));
double frac = curDetPos - detPosFloor; // fraction distance from det
if (iDetPos > 0 && iDetPos < nDet - 1)
*pImCol++ += filteredProj[iDetPos] + (frac * (filteredProj[iDetPos+1] - filteredProj[iDetPos]));
- } else if (interpType = Backprojector::INTERP_CUBIC) {
+ } else if (interpType =
=
Backprojector::INTERP_CUBIC) {
double p = iDetCenter + curDetPos; // position along detector
if (p >= 0 && p < nDet)
*pImCol++ += pCubicInterp->interpolate (p);
double p = iDetCenter + curDetPos; // position along detector
if (p >= 0 && p < nDet)
*pImCol++ += pCubicInterp->interpolate (p);
@@
-560,7
+589,7
@@
BackprojectIntDiff::BackprojectView (const double* const filteredProj, const dou
if (iDetPos >= 0 && iDetPos <= iLastDet)
*pImCol++ += filteredProj[iDetPos] + (detRemainder * deltaFilteredProj[iDetPos]);
} // end for iy
if (iDetPos >= 0 && iDetPos <= iLastDet)
*pImCol++ += filteredProj[iDetPos] + (detRemainder * deltaFilteredProj[iDetPos]);
} // end for iy
- } else if (interpType = Backprojector::INTERP_CUBIC) {
+ } else if (interpType =
=
Backprojector::INTERP_CUBIC) {
for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
*pImCol++ += pCubicInterp->interpolate (static_cast<double>(curDetPos) / 65536);
}
for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
*pImCol++ += pCubicInterp->interpolate (static_cast<double>(curDetPos) / 65536);
}
@@
-636,10
+665,8
@@
BackprojectEquilinear::BackprojectView (const double* const filteredProj, const
double dU = (m_dFocalLength + rsin_t) / m_dFocalLength;
double dDetPos = rcos_t / dU;
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. This assumes that the detector is also
- // located focal-length away from the origin.
- dDetPos *= 2;
+ // Scale for imaginary detector that passes through origin of phantom, see Kak-Slaney Figure 3.22.
+ dDetPos *= m_dSourceDetectorLength / m_dFocalLength;
double dPos = dDetPos / detInc; // position along detector array
if (interpType == Backprojector::INTERP_NEAREST) {
double dPos = dDetPos / detInc; // position along detector array
if (interpType == Backprojector::INTERP_NEAREST) {