projects / ctsim.git / blobdiff
? search:
summary | shortlog | log | commit | commitdiff | tree
raw | inline | side by side
r385: no message
[ctsim.git] / src / views.cpp
diff --git a/src/views.cpp b/src/views.cpp
index f52b6116e37c99685b3c5f11cc87af252a789d33..f8efe44f6949a45e748e8dc7c00c611734c8846c 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
 **
-**  $Id: views.cpp,v 1.53 2001/01/07 22:53:36 kevin Exp $
+**  $Id: views.cpp,v 1.60 2001/01/12 21:53:27 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
@@ -1056,7 +1056,8 @@ ImageFileView::OnPlotRow (wxCommandEvent& event)
   if (v != NULL && yCursor < ny) {
     double* pX = new double [nx];
     double* pYReal = new double [nx];
-    double *pYImag, *pYMag;
+    double *pYImag = NULL;
+    double *pYMag = NULL;
     if (rIF.isComplex()) {
       pYImag = new double [nx];
       pYMag = new double [nx];
@@ -1103,7 +1104,7 @@ ImageFileView::OnPlotRow (wxCommandEvent& event)
         rPlotFile.addColumn (2, pYImag);
         rPlotFile.addColumn (3, pYMag);
       }
-      for (int iL = 0; iL < rIF.nLabels(); iL++)
+      for (unsigned int iL = 0; iL < rIF.nLabels(); iL++)
         rPlotFile.addDescription (rIF.labelGet(iL).getLabelString().c_str());
       os << " Plot of " << GetDocument()->GetFirstView()->GetFrame()->GetTitle().c_str();
       *theApp->getLog() << os.str().c_str() << "\n";
@@ -1138,11 +1139,12 @@ ImageFileView::OnPlotCol (wxCommandEvent& event)
   
   if (v != NULL && xCursor < nx) {
     double* pX = new double [ny];
-    double* pYReal = new double [nx];
-    double *pYImag, *pYMag;
+    double* pYReal = new double [ny];
+    double* pYImag = NULL;
+    double* pYMag = NULL;
     if (rIF.isComplex()) {
-      pYImag = new double [nx];
-      pYMag = new double [nx];
+      pYImag = new double [ny];
+      pYMag = new double [ny];
     }
     for (int i = 0; i < ny; i++) {
       pX[i] = i;
@@ -1177,16 +1179,16 @@ ImageFileView::OnPlotCol (wxCommandEvent& event)
         rPlotFile.addEzsetCommand ("curve 3");
         rPlotFile.addEzsetCommand ("color 0");
         rPlotFile.addEzsetCommand ("solid");
-        rPlotFile.setCurveSize (4, nx);
+        rPlotFile.setCurveSize (4, ny);
       } else
-        rPlotFile.setCurveSize (2, nx);
+        rPlotFile.setCurveSize (2, ny);
       rPlotFile.addColumn (0, pX);
       rPlotFile.addColumn (1, pYReal); 
       if (rIF.isComplex()) {
         rPlotFile.addColumn (2, pYImag);
         rPlotFile.addColumn (3, pYMag);
       }
-      for (int iL = 0; iL < rIF.nLabels(); iL++)
+      for (unsigned int iL = 0; iL < rIF.nLabels(); iL++)
         rPlotFile.addDescription (rIF.labelGet(iL).getLabelString().c_str());
       os << " Plot of " << GetDocument()->GetFirstView()->GetFrame()->GetTitle().c_str();
       *theApp->getLog() << os.str().c_str() << "\n";
@@ -1272,7 +1274,6 @@ ImageFileView::OnPlotFFTRow (wxCommandEvent& event)
         rPlotFile.addEzsetCommand ("curve 3");
         rPlotFile.addEzsetCommand ("color 0");
         rPlotFile.addEzsetCommand ("solid");
-        rPlotFile.setCurveSize (4, nx);
        rPlotFile.addEzsetCommand ("box");
       rPlotFile.addEzsetCommand ("grid");
       rPlotFile.setCurveSize (4, nx);
@@ -1309,16 +1310,46 @@ ImageFileView::OnPlotFFTCol (wxCommandEvent& event)
   
   const ImageFile& rIF = dynamic_cast<ImageFileDocument*>(GetDocument())->getImageFile();
   ImageFileArrayConst v = rIF.getArray();
+  ImageFileArrayConst vImag = rIF.getImaginaryArray();
   int nx = rIF.nx();
   int ny = rIF.ny();
   
   if (v != NULL && xCursor < nx) {
+    fftw_complex* pcIn = new fftw_complex [ny];
+    double *pdTemp = new double [ny];
+
+    int i;
+    for (i = 0; i < ny; i++)
+      pdTemp[i] = v[xCursor][i];
+    Fourier::shuffleNaturalToFourierOrder (pdTemp, ny);
+    for (i = 0; i < ny; i++) 
+      pcIn[i].re = pdTemp[i];
+
+    for (i = 0; i < ny; i++) {
+      if (rIF.isComplex())
+        pdTemp[i] = vImag[xCursor][i];
+      else
+      pdTemp[i] = 0;
+    }
+    Fourier::shuffleNaturalToFourierOrder (pdTemp, ny);
+    for (i = 0; i < ny; i++)
+      pcIn[i].im = pdTemp[i];
+
+    fftw_plan plan = fftw_create_plan (ny, FFTW_BACKWARD, FFTW_IN_PLACE);
+    fftw_one (plan, pcIn, NULL);
+    fftw_destroy_plan (plan);
+
     double* pX = new double [ny];
-    double* pY = new double [ny];
-    for (int i = 0; i < ny; i++) {
+    double* pYReal = new double [ny];
+    double* pYImag = new double [ny];
+    double* pYMag = new double [ny];
+    for (i = 0; i < ny; i++) {
       pX[i] = i;
-      pY[i] = v[xCursor][i];
+      pYReal[i] = pcIn[i].re;
+      pYImag[i] = pcIn[i].im;
+      pYMag[i] = ::sqrt (pcIn[i].re * pcIn[i].re + pcIn[i].im * pcIn[i].im);
     }
+
     PlotFileDocument* pPlotDoc = dynamic_cast<PlotFileDocument*>(theApp->getDocManager()->CreateDocument("untitled.plt", wxDOC_SILENT));
     if (! pPlotDoc) {
       sys_error (ERR_SEVERE, "Internal error: unable to create Plot file");
@@ -1329,17 +1360,38 @@ ImageFileView::OnPlotFFTCol (wxCommandEvent& event)
       std::string title("title ");
       title += os.str();
       rPlotFile.addEzsetCommand (title.c_str());
-      rPlotFile.addEzsetCommand ("xlabel Row");
+      rPlotFile.addEzsetCommand ("xlabel Column");
       rPlotFile.addEzsetCommand ("ylabel Pixel Value");
       rPlotFile.addEzsetCommand ("lxfrac 0");
-      rPlotFile.addEzsetCommand ("box");
+      rPlotFile.addEzsetCommand ("curve 1");
+      rPlotFile.addEzsetCommand ("color 1");
+       rPlotFile.addEzsetCommand ("dash 1");
+        rPlotFile.addEzsetCommand ("curve 2");
+        rPlotFile.addEzsetCommand ("color 4");
+        rPlotFile.addEzsetCommand ("dash 3");
+        rPlotFile.addEzsetCommand ("curve 3");
+        rPlotFile.addEzsetCommand ("color 0");
+        rPlotFile.addEzsetCommand ("solid");
+       rPlotFile.addEzsetCommand ("box");
       rPlotFile.addEzsetCommand ("grid");
-      rPlotFile.setCurveSize (2, nx);
+      rPlotFile.setCurveSize (4, ny);
       rPlotFile.addColumn (0, pX);
-      rPlotFile.addColumn (1, pY);
+      rPlotFile.addColumn (1, pYReal);
+      rPlotFile.addColumn (2, pYImag);
+      rPlotFile.addColumn (3, pYMag);
+      for (int iL = 0; iL < rIF.nLabels(); iL++)
+        rPlotFile.addDescription (rIF.labelGet(iL).getLabelString().c_str());
+      os << " FFT Plot of " << GetDocument()->GetFirstView()->GetFrame()->GetTitle().c_str();
+      *theApp->getLog() << os.str().c_str() << "\n";
+      rPlotFile.addDescription (os.str().c_str());
     }
     delete pX;
-    delete pY;
+    delete pYReal;
+    delete pYImag;
+    delete pYMag;
+    delete pdTemp;
+    delete [] pcIn;
+
     if (theApp->getSetModifyNewDocs())
       pPlotDoc->Modify(true);
     pPlotDoc->UpdateAllViews();
@@ -1403,26 +1455,26 @@ ImageFileView::OnCompareCol (wxCommandEvent& event)
         rPlotFile.addEzsetCommand ("dash 5");
         rPlotFile.addEzsetCommand ("box");
         rPlotFile.addEzsetCommand ("grid");
-        rPlotFile.setCurveSize (3, nx);
+        rPlotFile.setCurveSize (3, ny);
         rPlotFile.addColumn (0, pX);
         rPlotFile.addColumn (1, pY1);
         rPlotFile.addColumn (2, pY2);
 
         unsigned int iL;
         for (iL = 0; iL < rIF.nLabels(); iL++) {
-          std::string s = GetDocument()->GetFirstView()->GetFrame()->GetTitle();
+          std::string s = GetDocument()->GetFirstView()->GetFrame()->GetTitle().c_str();
           s += ": ";
           s += rIF.labelGet(iL).getLabelString();
           rPlotFile.addDescription (s.c_str());
         }
-        for (iL = 0; iL < rIF.nLabels(); iL++) {
-          std::string s = pCompareDoc->GetFirstView()->GetFrame()->GetTitle();
+        for (iL = 0; iL < rCompareIF.nLabels(); iL++) {
+          std::string s = pCompareDoc->GetFirstView()->GetFrame()->GetTitle().c_str();
           s += ": ";
           s += rCompareIF.labelGet(iL).getLabelString();
           rPlotFile.addDescription (s.c_str());
         }
-        os << " Between " << GetDocument()->GetFirstView()->GetFrame()->GetTitle() << " and "
-          << pCompareDoc->GetFirstView()->GetFrame()->GetTitle();
+        os << " Between " << GetDocument()->GetFirstView()->GetFrame()->GetTitle().c_str() << " and "
+          << pCompareDoc->GetFirstView()->GetFrame()->GetTitle().c_str();
         *theApp->getLog() << os.str().c_str() << "\n";
         rPlotFile.addDescription (os.str().c_str());
       }
@@ -1494,25 +1546,25 @@ ImageFileView::OnCompareRow (wxCommandEvent& event)
         rPlotFile.addEzsetCommand ("dash 5");
         rPlotFile.addEzsetCommand ("box");
         rPlotFile.addEzsetCommand ("grid");
-        rPlotFile.setCurveSize (3, ny);
+        rPlotFile.setCurveSize (3, nx);
         rPlotFile.addColumn (0, pX);
         rPlotFile.addColumn (1, pY1);
         rPlotFile.addColumn (2, pY2);
         unsigned int iL;
         for (iL = 0; iL < rIF.nLabels(); iL++) {
-          std::string s = GetDocument()->GetFirstView()->GetFrame()->GetTitle();
+          std::string s = GetDocument()->GetFirstView()->GetFrame()->GetTitle().c_str();
           s += ": ";
           s += rIF.labelGet(iL).getLabelString();
           rPlotFile.addDescription (s.c_str());
         }
-        for (iL = 0; iL < rIF.nLabels(); iL++) {
-          std::string s = pCompareDoc->GetFirstView()->GetFrame()->GetTitle();
+        for (iL = 0; iL < rCompareIF.nLabels(); iL++) {
+          std::string s = pCompareDoc->GetFirstView()->GetFrame()->GetTitle().c_str();
           s += ": ";
           s += rCompareIF.labelGet(iL).getLabelString();
           rPlotFile.addDescription (s.c_str());
         }
-        os << " Between " << GetDocument()->GetFirstView()->GetFrame()->GetTitle() << " and "
-          << pCompareDoc->GetFirstView()->GetFrame()->GetTitle();
+        os << " Between " << GetDocument()->GetFirstView()->GetFrame()->GetTitle().c_str() << " and "
+          << pCompareDoc->GetFirstView()->GetFrame()->GetTitle().c_str();
         *theApp->getLog() << os.str().c_str() << "\n";
         rPlotFile.addDescription (os.str().c_str());
       }
@@ -1570,16 +1622,16 @@ ImageFileView::OnPlotHistogram (wxCommandEvent& event)
     rPlotFile.addEzsetCommand ("ylabel Count");
     rPlotFile.addEzsetCommand ("box");
     rPlotFile.addEzsetCommand ("grid");
-    rPlotFile.setCurveSize (2, nx);
+    rPlotFile.setCurveSize (2, NUMBER_HISTOGRAM_BINS);
     rPlotFile.addColumn (0, pX);
     rPlotFile.addColumn (1, pY);
     for (unsigned int iL = 0; iL < rIF.nLabels(); iL++) {
-      std::string s = GetDocument()->GetFirstView()->GetFrame()->GetTitle();
+      std::string s = GetDocument()->GetFirstView()->GetFrame()->GetTitle().c_str();
       s += ": ";
       s += rIF.labelGet(iL).getLabelString();
       rPlotFile.addDescription (s.c_str());
     }
-    os << " Plot of " << GetDocument()->GetFirstView()->GetFrame()->GetTitle();
+    os << " Plot of " << GetDocument()->GetFirstView()->GetFrame()->GetTitle().c_str();
     *theApp->getLog() << os.str().c_str() << "\n";
     rPlotFile.addDescription (os.str().c_str());
     delete pX;
@@ -1957,7 +2009,7 @@ ProjectionFileView::ProjectionFileView(void)
   m_dDefaultFilterParam = 1.;
 #if HAVE_FFTW
   m_iDefaultFilterMethod = ProcessSignal::FILTER_METHOD_RFFTW;
-  m_iDefaultFilterGeneration = ProcessSignal::FILTER_GENERATION_INVERSE_FOURIER;
+  m_iDefaultFilterGeneration = ProcessSignal::FILTER_GENERATION_DIRECT;
 #else
   m_iDefaultFilterMethod = ProcessSignal::FILTER_METHOD_CONVOLUTION;
   m_iDefaultFilterGeneration = ProcessSignal::FILTER_GENERATION_DIRECT;
@@ -2095,11 +2147,6 @@ ProjectionFileView::OnReconstructFBP (wxCommandEvent& event)
       const Projections& rProj = GetDocument()->getProjections();
       imageFile.setArraySize (m_iDefaultNX, m_iDefaultNY);
       
-      if (m_iDefaultFilterMethod != ProcessSignal::FILTER_METHOD_CONVOLUTION && m_iDefaultFilterGeneration == ProcessSignal::FILTER_GENERATION_DIRECT && rProj.geometry() != Scanner::GEOMETRY_PARALLEL) {
-        wxMessageBox ("Sorry!\nCurrently, frequency-based filtering with direct filter generation is not support for geometries other than parallel.\nAborting command.", "Not Supported", wxOK | wxICON_WARNING, m_frame);
-        return;
-      }
-
       Reconstructor* pReconstruct = new Reconstructor (rProj, imageFile, optFilterName.c_str(), m_dDefaultFilterParam, optFilterMethodName.c_str(), m_iDefaultZeropad, optFilterGenerationName.c_str(), optInterpName.c_str(), m_iDefaultInterpParam, optBackprojectName.c_str(), m_iDefaultTrace);
 
       Timer timerRecon;
Computed Tomography simulator and reconstructor