Convert pol to C++ object.
-Add loading / saving of arbitrary phantom objections in ASCII format.
-
Adding animation of X-ray projection collection
Add animation of filtering & backprojections
+FFT filtering requires FFTW library. Add a public domain FFT routine.
+
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: filter.cpp,v 1.26 2000/08/22 07:02:48 kevin Exp $
+** $Id: filter.cpp,v 1.27 2000/08/22 16:49:56 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 SignalFilter::FILTER_INVALID = -1 ;
const int SignalFilter::FILTER_ABS_BANDLIMIT = 0; // filter times |x = |
-const int SignalFilter::FILTER_ABS_SINC = 1;
-const int SignalFilter::FILTER_ABS_G_HAMMING = 2;
-const int SignalFilter::FILTER_ABS_COSINE = 3;
+const int SignalFilter::FILTER_ABS_G_HAMMING = 1;
+const int SignalFilter::FILTER_ABS_COSINE = 2;
+const int SignalFilter::FILTER_ABS_SINC = 3;
const int SignalFilter::FILTER_SHEPP = 4;
const int SignalFilter::FILTER_BANDLIMIT = 5;
const int SignalFilter::FILTER_SINC = 6;
const char* SignalFilter::s_aszFilterName[] = {
{"abs_bandlimit"},
- {"abs_sinc"},
{"abs_hamming"},
{"abs_cosine"},
+ {"abs_sinc"},
{"shepp"},
{"bandlimit"},
{"sinc"},
const char* SignalFilter::s_aszFilterTitle[] = {
{"Abs(w) * Bandlimit"},
- {"Abs(w) * Sinc"},
{"Abs(w) * Hamming"},
{"Abs(w) * Cosine"},
+ {"Abs(w) * Sinc"},
{"Shepp"},
{"Bandlimit"},
{"Sinc"},
switch (filterID) {
case FILTER_BANDLIMIT:
- if (fabs(au) > fabs(bw / 2) + F_EPSILON)
+ if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
q = 0.;
else
q = 1;
break;
case FILTER_ABS_BANDLIMIT:
- if (fabs(au) > fabs(bw / 2) + F_EPSILON)
+ if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
q = 0.;
else
q = au;
break;
case FILTER_TRIANGLE:
- if (fabs(au) > fabs(bw / 2) + F_EPSILON)
+ if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
q = 0;
else
q = 1 - au / bw;
break;
case FILTER_COSINE:
- if (fabs(au) > fabs(bw / 2) + F_EPSILON)
+ if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
q = 0;
else
q = cos(PI * u / bw);
break;
case FILTER_ABS_COSINE:
- if (fabs(au) > fabs(bw / 2) + F_EPSILON)
+ if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
q = 0;
else
q = au * cos(PI * u / bw);
q = au * bw * sinc (PI * bw * u, 1.);
break;
case FILTER_G_HAMMING:
- if (fabs(au) > fabs(bw / 2) + F_EPSILON)
+ if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
q = 0;
else
q = param + (1 - param) * cos (TWOPI * u / bw);
break;
case FILTER_ABS_G_HAMMING:
- if (fabs(au) > fabs(bw / 2) + F_EPSILON)
+ if (fabs(au) >= fabs(bw / 2) + F_EPSILON)
q = 0;
else
q = au * (param + (1 - param) * cos(TWOPI * u / bw));
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: imagefile.cpp,v 1.9 2000/07/15 08:36:13 kevin Exp $
+** $Id: imagefile.cpp,v 1.10 2000/08/22 16:49:56 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
}
mean /= (nx * ny);
- static const int nbin = 256;
+ static const int nbin = 1024;
int hist[ nbin ] = {0};
double spread = max - min;
mode = 0;
}
mode = (max_binindex * spread / (nbin - 1)) + min;
-
- median = 0.;
+
+ int nPixels = nx * ny;
+ slist<double> vecImage;
+ for (int ix = 0; ix < nx; ix++)
+ for (int iy = 0; iy < ny; iy++)
+ vecImage.push_front (v[ix][iy]);
+ vecImage.sort();
+ slist<double>::const_iterator iter = vecImage.begin();
+ for (int i = 0; i < (nPixels / 2) - 1; i++)
+ iter++; // Advance iterator to (nPixels / 2) - 1;
+
+ if (nPixels % 2) { // Odd
+ iter++;
+ median = *iter;
+ } else // Even
+ median = (*iter++ + *iter) / 2;
}
** This is part of the CTSim program
** Copyright (C) 1983-2000 Kevin Rosenberg
**
-** $Id: views.cpp,v 1.14 2000/08/22 07:02:48 kevin Exp $
+** $Id: views.cpp,v 1.15 2000/08/22 16:49:56 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
ProjectionFileView::OnReconstruct (wxCommandEvent& event)
{
- DialogGetReconstructionParameters dialogReconstruction (m_frame, 256, 256, SignalFilter::FILTER_ABS_BANDLIMIT, 1., ProcessSignal::FILTER_METHOD_CONVOLUTION, ProcessSignal::FILTER_GENERATION_INVALID, 3, Backprojector::INTERP_LINEAR, 1, Backprojector::BPROJ_IDIFF3);
+#if HAVE_FFTW
+ DialogGetReconstructionParameters dialogReconstruction (m_frame, 256, 256, SignalFilter::FILTER_ABS_BANDLIMIT, 1., ProcessSignal::FILTER_METHOD_RFFTW, ProcessSignal::FILTER_GENERATION_INVERSE_FOURIER, 1, Backprojector::INTERP_LINEAR, 1, Backprojector::BPROJ_IDIFF3);
+#else
+ DialogGetReconstructionParameters dialogReconstruction (m_frame, 256, 256, SignalFilter::FILTER_ABS_BANDLIMIT, 1., ProcessSignal::FILTER_METHOD_CONVOLUTION, ProcessSignal::FILTER_GENERATION_DIRECT, 1, Backprojector::INTERP_LINEAR, 1, Backprojector::BPROJ_IDIFF3);
+#endif
int retVal = dialogReconstruction.ShowModal();
if (retVal == wxID_OK) {
int xSize = dialogReconstruction.getXSize();
pReconDoc->Modify(true);
pReconDoc->UpdateAllViews(this);
ostringstream os;
- os << "Reconstruct " << rProj.getFilename() << ": xSize=" << xSize << ", ySize=" << ySize << ", Filter=" << optFilterName.c_str() << ", FilterParam=" << optFilterParam << ", FilterMethod=" << optFilterMethodName.c_str() << ", Zeropad=" << optZeropad << ", Interpolation=" << optInterpName.c_str() << ", InterpolationParam=" << optInterpParam << ", Backprojection=" << optBackprojectName.c_str() << "\n";
+ os << "Reconstruct " << rProj.getFilename() << ": xSize=" << xSize << ", ySize=" << ySize << ", Filter=" << optFilterName.c_str() << ", FilterParam=" << optFilterParam << ", FilterMethod=" << optFilterMethodName.c_str() << ", FilterGeneration=" << optFilterGenerationName.c_str() << ", Zeropad=" << optZeropad << ", Interpolation=" << optInterpName.c_str() << ", InterpolationParam=" << optInterpParam << ", Backprojection=" << optBackprojectName.c_str() << "\n";
*theApp->getLog() << os.str().c_str();
}
}