** 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.15 2000/09/07 14:29:05 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 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"},
- {"equilinear"},
{"equiangular"},
+ {"equilinear"},
};
const char* Scanner::s_aszGeometryTitle[] =
{
{"Parallel"},
- {"Equilinear"},
{"Equiangular"},
+ {"Equilinear"},
};
const int Scanner::s_iGeometryCount = sizeof(s_aszGeometryName) / sizeof(const char*);
return;
}
- if (nView < 1)
- nView = 1;
+ if (nView < 1 || nDet < 1) {
+ m_fail = true;
+ m_failMessage = "nView & nDet must be greater than 0";
+ return;
+ }
if (nSample < 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_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;
- m_rotLen = rot_anglen;
- m_rotInc = m_rotLen / m_nView;
+ if (m_nDet % 2 == 0) // Adjust for Even number of detectors
+ m_detInc = m_detLen / (m_nDet - 1); // center detector = (nDet/2)-1
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.ys2 = m_dYCenter + dHalfDetLen;
- 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.xs1 = m_dXCenter - dHalfDetLen;
+ m_initPos.ys1 = m_dYCenter + m_dFocalLength;
+ m_initPos.xs2 = m_dXCenter + dHalfDetLen;
+ 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_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.) {
return;
}
double dAngle = atan( dHalfSquare / dFocalPastPhm );
+#endif
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;
+ if (m_nDet % 2 == 0) // Adjust for Even number of detectors
+ m_detInc = m_detLen / (m_nDet - 1); // center detector = (nDet/2)-1
- 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) {
+#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.) {
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_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;
+ if (m_nDet % 2 == 0) // Adjust for Even number of detectors
+ m_detInc = m_detLen / (m_nDet - 1); // center detector = (nDet/2)-1
+ 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
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;
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->setColor (C_RED);
+ m_pSGP->setColor (C_BLUE);
phm.draw (*m_pSGP);
m_dTextHeight = m_pSGP->getCharHeight ();
- traceShowParam ("Projection Collector", "%s", PROJECTION_TRACE_ROW_TITLE, C_BLACK, " ");
- traceShowParam ("________________", "%s", PROJECTION_TRACE_ROW_TITLE2, C_LTGRAY, " ");
traceShowParam ("Phantom:", "%s", PROJECTION_TRACE_ROW_PHANT_ID, C_BLACK, phm.name().c_str());
traceShowParam ("Geometry:", "%s", PROJECTION_TRACE_ROW_GEOMETRY, C_BLUE, convertGeometryIDToName(m_idGeometry));
traceShowParam ("Focal Length Ratio:", "%.2f", PROJECTION_TRACE_ROW_FOCAL_LENGTH, C_BLUE, m_dFocalLengthRatio);
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)
- 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
- 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
#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);
}
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;
- 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
#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->setRasterOp (RO_SET);
}
#endif
sum += projectSingleLine (phm, xd, yd, xs, ys);
#ifdef HAVE_SGP
- if (m_trace >= Trace::TRACE_CLIPPING) {
- 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);
+ // if (m_trace >= Trace::TRACE_CLIPPING) {
+ // traceShowParam ("Attenuation:", "%s", PROJECTION_TRACE_ROW_ATTEN, C_LTMAGENTA, " ");
+ // traceShowParam ("Attenuation:", "%.3f", PROJECTION_TRACE_ROW_ATTEN, C_LTMAGENTA, sum);
// }
#endif
if (m_idGeometry == GEOMETRY_EQUIANGULAR)
{
va_list arg;
va_start(arg, color);
+#ifdef HAVE_SGP
traceShowParamRasterOp (RO_COPY, szLabel, fmt, row, color, arg);
+#else
+ traceShowParamRasterOp (0, szLabel, fmt, row, color, arg);
+#endif
va_end(arg);
}
{
va_list arg;
va_start(arg, color);
+#ifdef HAVE_SGP
traceShowParamRasterOp (RO_XOR, szLabel, fmt, row, color, arg);
+#else
+ traceShowParamRasterOp (0, szLabel, fmt, row, color, arg);
+#endif
va_end(arg);
}
// cio_set_text_clr (color - 8, 0);
// cio_set_text_clr (color, 0);
+#ifdef HAVE_SGP
if (m_pSGP) {
m_pSGP->setRasterOp (iRasterOp);
double dYPos = m_dYMaxWin - (row * m_dTextHeight);
m_pSGP->moveAbs (m_dXMinWin, dYPos);
m_pSGP->setTextColor (color, -1);
m_pSGP->drawText (szLabel);
- double dValueOffset = (m_dXMaxWin - m_dXMinWin) / 5;
+ double dValueOffset = (m_dXMaxWin - m_dXMinWin) / 4;
m_pSGP->moveAbs (m_dXMinWin + dValueOffset, dYPos);
m_pSGP->drawText (szValue);
- } else {
- cio_put_str (szLabel);
- cio_put_str (szValue);
- cio_put_str ("\n");
- }
+ } else
+#endif
+ {
+ cio_put_str (szLabel);
+ cio_put_str (szValue);
+ cio_put_str ("\n");
+ }
}