r328: *** empty log message ***
[ctsim.git] / tools / phm2if.cpp
index 9b8fd22f89263b3be6405ff8f51463d5a1c8b835..e0ca9c81238e97137226410af2f2ed05bff77ab5 100644 (file)
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
-**  $Id: phm2if.cpp,v 1.1 2000/07/13 07:01:35 kevin Exp $
+**  $Id: phm2if.cpp,v 1.17 2001/01/02 05:34:57 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
@@ -50,91 +50,89 @@ static struct option my_options[] =
   {0, 0, 0, 0}
 };
 
+static const char* g_szIdStr = "$Id: phm2if.cpp,v 1.17 2001/01/02 05:34:57 kevin Exp $";
+
 void 
 phm2if_usage (const char *program)
 {
-  cout << "phm2if_usage: " << fileBasename(program) << " outfile nx ny [--phantom phantom-name] [--phmfile filename] [--filter filter-name] [OPTIONS]" << endl;
-  cout << "Generate phantom image from a predefined --phantom or a --phmfile" << endl;
-  cout << endl;
-  cout << "     outfile         Name of output file for image" << endl;
-  cout << "     nx              Number of pixels X-axis" << endl;
-  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 << "        unitpulse    Unit pulse phantom" << endl;
-  cout << "     --phmfile       Generate Phantom from phantom file" << endl;
-  cout << "     --filter        Generate Phantom from a filter function" << endl;
-  cout << "       abs_bandlimit Abs * Bandlimiting" << endl;
-  cout << "       abs_sinc      Abs * Sinc" << endl;
-  cout << "       abs_cos       Abs * Cosine" << endl;
-  cout << "       abs_hamming   Abs * Hamming" << endl;
-  cout << "       shepp         Shepp-Logan" << endl;
-  cout << "       bandlimit     Bandlimiting" << endl;
-  cout << "       sinc          Sinc" << endl;
-  cout << "       cos           Cosine" << endl;
-  cout << "       triangle      Triangle" << endl;
-  cout << "       hamming       Hamming" << endl;
-  cout << "     --filter-param  Alpha level for Hamming filter" << endl;
-  cout << "     --filter-domain Set domain of filter" << endl;
-  cout << "         spatial     Spatial domain (default)" << endl;
-  cout << "         freq        Frequency domain" << endl;
-  cout << "     --filter-bw     Filter bandwidth (default = 1)" << endl;
-  cout << "     --desc          Description of raysum" << endl;
-  cout << "     --nsample       Number of samples per axis per pixel (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" << endl;
-  cout << "        rays         Trace rays" << endl;
-  cout << "        plot         Trace plot" << endl;
-  cout << "        clipping     Trace clipping" << endl;
-  cout << "     --debug         Debug mode" << endl;
-  cout << "     --verbose       Verbose mode" << endl;
-  cout << "     --version       Print version" << endl;
-  cout << "     --help          Print this help message" << endl;
+  std::cout << "phm2if_usage: " << fileBasename(program) << " outfile nx ny [--phantom phantom-name] [--phmfile filename] [--filter filter-name] [OPTIONS]\n";
+  std::cout << "Generate phantom image from a predefined --phantom or a --phmfile\n";
+  std::cout << std::endl;
+  std::cout << "     outfile         Name of output file for image\n";
+  std::cout << "     nx              Number of pixels X-axis\n";
+  std::cout << "     ny              Number of pixels Y-axis\n";
+  std::cout << "     --phantom       Phantom to use for projection\n";
+  std::cout << "        herman       Herman head phantom\n";
+  std::cout << "        herman-b     Herman head phantom (Bordered)\n";
+  std::cout << "        shepp-logan  Shepp-Logan head phantom\n";
+  std::cout << "        shepp-logan-b Shepp-Logan head phantom (Bordered)\n";
+  std::cout << "        unitpulse    Unit pulse phantom\n";
+  std::cout << "     --phmfile       Generate Phantom from phantom file\n";
+  std::cout << "     --filter        Generate Phantom from a filter function\n";
+  std::cout << "       abs_bandlimit Abs * Bandlimiting\n";
+  std::cout << "       abs_sinc      Abs * Sinc\n";
+  std::cout << "       abs_cos       Abs * Cosine\n";
+  std::cout << "       abs_hamming   Abs * Hamming\n";
+  std::cout << "       shepp         Shepp-Logan\n";
+  std::cout << "       bandlimit     Bandlimiting\n";
+  std::cout << "       sinc          Sinc\n";
+  std::cout << "       cos           Cosine\n";
+  std::cout << "       triangle      Triangle\n";
+  std::cout << "       hamming       Hamming\n";
+  std::cout << "     --filter-param  Alpha level for Hamming filter\n";
+  std::cout << "     --filter-domain Set domain of filter\n";
+  std::cout << "         spatial     Spatial domain (default)\n";
+  std::cout << "         freq        Frequency domain\n";
+  std::cout << "     --filter-bw     Filter bandwidth (default = 1)\n";
+  std::cout << "     --desc          Description of raysum\n";
+  std::cout << "     --nsample       Number of samples per axis per pixel (default = 1)\n";
+  std::cout << "     --trace         Trace level to use\n";
+  std::cout << "        none         No tracing (default)\n";
+  std::cout << "        console      Trace text level\n";
+  std::cout << "        phantom      Trace phantom\n";
+  std::cout << "        proj         Trace projections\n";
+  std::cout << "        plot         Trace plot\n";
+  std::cout << "        clipping     Trace clipping\n";
+  std::cout << "     --debug         Debug mode\n";
+  std::cout << "     --verbose       Verbose mode\n";
+  std::cout << "     --version       Print version\n";
+  std::cout << "     --help          Print this help message\n";
 }
 
 #ifdef HAVE_MPI
-void mpi_gather_image (MPIWorld& mpiWorld, ImageFile* imGlobal, ImageFile* imLocal, const int opt_debug);
+void mpi_gather_image (MPIWorld& mpiWorld, ImageFile* pImGlobal, ImageFile* pImLocal, const int optDebug);
 #endif
 
 int 
 phm2if_main (int argc, char* argv[])
 {
-  ImageFile* imGlobal = NULL;
+  ImageFile* pImGlobal = NULL;
   Phantom phm;
-  int opt_nx = 0, opt_ny = 0;
+  std::string optPhmName;
+  std::string optFilterName;
+  std::string optDomainName (SignalFilter::convertDomainIDToName (SignalFilter::DOMAIN_SPATIAL));
+  std::string optOutFilename;
+  std::string optDesc;
+  std::string optPhmFilename;
+  int opt_nx = 0;
+  int opt_ny = 0;
   int opt_nsample = 1;
-  string optPhmName = Phantom::PHM_HERMAN_STR;
-  string optFilterName;
-  string optDomainName = SignalFilter::DOMAIN_SPATIAL_STR;
-  char *opt_outfile = NULL;
-  int opt_debug = 0;
-  string opt_desc;
-  string opt_phmFileName;
-  char *endstr, *endptr;
-  double opt_filter_param = 1;
-  double opt_filter_bw = 1.;
-  int opt_trace = TRACE_NONE;
-  bool opt_verbose = false;
-#ifdef HAVE_MPI
-  ImageFile* imLocal = NULL;
-  MPIWorld mpiWorld (argc, argv);
-#endif
-
+  double optFilterParam = 1.;
+  double optFilterBW = 1.;
+  int optTrace = Trace::TRACE_NONE;
+  bool optVerbose = false;
+  bool optDebug = false;
+  char *endptr = NULL;
+  char *endstr;
   Timer timerProgram;
-  
+
 #ifdef HAVE_MPI
+  ImageFile* pImLocal = NULL;
+  MPIWorld mpiWorld (argc, argv);
   if (mpiWorld.getRank() == 0) {
 #endif
     while (1) {
       int c = getopt_long(argc, argv, "", my_options, NULL);
-      char *endptr = NULL;
-      char *endstr;
-      
       if (c == -1)
        break;
       
@@ -143,16 +141,16 @@ phm2if_main (int argc, char* argv[])
        optPhmName = optarg;
        break;
       case O_PHMFILE:
-       opt_phmFileName = optarg;
+       optPhmFilename = optarg;
        break;
       case O_VERBOSE:
-       opt_verbose = true;
+       optVerbose = true;
        break;
       case O_DEBUG:
-       opt_debug = 1;
+       optDebug = true;
        break;
       case O_TRACE:
-       if ((opt_trace = convertTraceNameToID(optarg)) == TRACE_INVALID) {
+       if ((optTrace = Trace::convertTraceNameToID(optarg)) == Trace::TRACE_INVALID) {
          phm2if_usage(argv[0]);
          return (1);
        }
@@ -164,10 +162,10 @@ phm2if_main (int argc, char* argv[])
        optDomainName = optarg;
        break;
       case O_DESC:
-       opt_desc =  optarg;
+       optDesc =  optarg;
        break;
       case O_FILTER_BW:
-       opt_filter_bw = strtod(optarg, &endptr);
+       optFilterBW = strtod(optarg, &endptr);
        endstr = optarg + strlen(optarg);
        if (endptr != endstr) {
          sys_error(ERR_SEVERE,"Error setting --filter-bw to %s\n", optarg);
@@ -176,7 +174,7 @@ phm2if_main (int argc, char* argv[])
        }
        break;
       case O_FILTER_PARAM:
-       opt_filter_param = strtod(optarg, &endptr);
+       optFilterParam = strtod(optarg, &endptr);
        endstr = optarg + strlen(optarg);
        if (endptr != endstr) {
          sys_error(ERR_SEVERE,"Error setting --filter-param to %s\n", optarg);
@@ -195,9 +193,9 @@ phm2if_main (int argc, char* argv[])
        break;
         case O_VERSION:
 #ifdef VERSION
-         cout <<  "Version " << VERSION << endl;
+         std::cout << "Version " << VERSION << std::endl << g_szIdStr << std::endl;
 #else
-          cerr << "Unknown version number" << endl;
+          std::cerr << "Unknown version number\n";
 #endif
       case O_HELP:
       case '?':
@@ -209,8 +207,8 @@ phm2if_main (int argc, char* argv[])
       }
     }
     
-    if (optPhmName == "" && optFilterName == "" && opt_phmFileName == "") {
-      cerr << "No phantom defined" << endl << endl;
+    if (optPhmName == "" && optFilterName == "" && optPhmFilename == "") {
+      std::cerr << "No phantom defined\n" << std::endl;
       phm2if_usage(argv[0]);
       return (1);
     }
@@ -219,7 +217,7 @@ phm2if_main (int argc, char* argv[])
       phm2if_usage(argv[0]);
       return (1);
     }
-    opt_outfile = argv[optind];
+    optOutFilename = argv[optind];
     opt_nx = strtol(argv[optind+1], &endptr, 10);
     endstr = argv[optind+1] + strlen(argv[optind+1]);
     if (endptr != endstr) {
@@ -235,39 +233,39 @@ phm2if_main (int argc, char* argv[])
       return (1);
     }
     
-    ostringstream oss;
+    std::ostringstream oss;
     oss << "phm2if: nx=" << opt_nx << ", ny=" << opt_ny << ", nsample=" << opt_nsample << ", ";
-    if (opt_phmFileName != "")
-      oss << "phantom=" << opt_phmFileName;
+    if (optPhmFilename != "")
+      oss << "phantomFile=" << optPhmFilename;
     else if (optPhmName != "")
       oss << "phantom=" << optPhmName;
     else if (optFilterName != "") {
       oss << "filter=" << optFilterName << " - " << optDomainName;
     }
-    if (opt_desc != "")
-      oss << ": " << opt_desc;
-    opt_desc = oss.str();
+    if (optDesc != "")
+      oss << ": " << optDesc;
+    optDesc = oss.str();
     
     if (optPhmName != "") {
       phm.createFromPhantom (optPhmName.c_str());
       if (phm.fail()) {
-       cout << phm.failMessage() << endl << endl;
+       std::cout << phm.failMessage() << std::endl << std::endl;
        phm2if_usage(argv[0]);
        return (1);
       }
     }
 
-    if (opt_phmFileName != "") {
-      phm.createFromFile(opt_phmFileName.c_str());
+    if (optPhmFilename != "") {
+      phm.createFromFile(optPhmFilename.c_str());
 #ifdef HAVE_MPI
       if (mpiWorld.getRank() == 0) 
-       cerr << "Can't use phantom from file in MPI mode" << endl;
+       std::cerr << "Can't use phantom from file in MPI mode\n";
       return (1);
 #endif
     }
 
-    if (opt_verbose)
-      cout << "Rasterize Phantom to Image" << endl << endl;
+    if (optVerbose)
+      std::cout << "Rasterize Phantom to Image\n" << std::endl;
 #ifdef HAVE_MPI
   }
 #endif
@@ -275,19 +273,19 @@ phm2if_main (int argc, char* argv[])
 #ifdef HAVE_MPI
   TimerCollectiveMPI timerBcast (mpiWorld.getComm());
   mpiWorld.BcastString (optPhmName);
-  mpiWorld.getComm().Bcast (&opt_verbose, 1, MPI::INT, 0);
-  mpiWorld.getComm().Bcast (&opt_debug, 1, MPI::INT, 0);
-  mpiWorld.getComm().Bcast (&opt_trace, 1, MPI::INT, 0);
+  mpiWorld.getComm().Bcast (&optVerbose, 1, MPI::INT, 0);
+  mpiWorld.getComm().Bcast (&optDebug, 1, MPI::INT, 0);
+  mpiWorld.getComm().Bcast (&optTrace, 1, MPI::INT, 0);
   mpiWorld.getComm().Bcast (&opt_nx, 1, MPI::INT, 0);
   mpiWorld.getComm().Bcast (&opt_ny, 1, MPI::INT, 0);
   mpiWorld.getComm().Bcast (&opt_nsample, 1, MPI::INT, 0);
-  mpiWorld.getComm().Bcast (&opt_filter_param, 1, MPI::DOUBLE, 0);
-  mpiWorld.getComm().Bcast (&opt_filter_bw, 1, MPI::DOUBLE, 0);
+  mpiWorld.getComm().Bcast (&optFilterParam, 1, MPI::DOUBLE, 0);
+  mpiWorld.getComm().Bcast (&optFilterBW, 1, MPI::DOUBLE, 0);
 
   mpiWorld.BcastString (optFilterName);
   mpiWorld.BcastString (optDomainName);
 
-  if (opt_verbose)
+  if (optVerbose)
     timerBcast.timerEndAndReport ("Time to broadcast variables");
 
   mpiWorld.setTotalWorkUnits (opt_nx);
@@ -296,17 +294,17 @@ phm2if_main (int argc, char* argv[])
       phm.createFromPhantom (optPhmName.c_str());
 
   if (mpiWorld.getRank() == 0) {
-    imGlobal = new ImageFile (opt_nx, opt_ny);
+    pImGlobal = new ImageFile (opt_nx, opt_ny);
   }
-  imLocal = new ImageFile (opt_nx, opt_ny);
+  pImLocal = new ImageFile (opt_nx, opt_ny);
 #else
-  imGlobal = new ImageFile (opt_nx, opt_ny);
+  pImGlobal = new ImageFile (opt_nx, opt_ny);
 #endif
 
-  ImageFileArray v;
+  ImageFileArray v = NULL;
 #ifdef HAVE_MPI
   if (mpiWorld.getRank() == 0)
-    v = imGlobal->getArray ();
+    v = pImGlobal->getArray ();
 
   if (phm.getComposition() == P_UNIT_PULSE) {
     if (mpiWorld.getRank() == 0) {
@@ -314,31 +312,31 @@ phm2if_main (int argc, char* argv[])
     }
   } else if (optFilterName != "") {
     if (mpiWorld.getRank() == 0) {
-      imGlobal->filterResponse (optDomainName.c_str(), opt_filter_bw, optFilterName.c_str(), opt_filter_param);
+      pImGlobal->filterResponse (optDomainName.c_str(), optFilterBW, optFilterName.c_str(), optFilterParam);
     }
   } else {
     TimerCollectiveMPI timerRasterize (mpiWorld.getComm());
-    phm.convertToImagefile (*imLocal, opt_nsample, opt_trace, mpiWorld.getMyStartWorkUnit(), mpiWorld.getMyLocalWorkUnits());
-    if (opt_verbose)
+    phm.convertToImagefile (*pImLocal, opt_nsample, optTrace, mpiWorld.getMyStartWorkUnit(), mpiWorld.getMyLocalWorkUnits(), false);
+    if (optVerbose)
       timerRasterize.timerEndAndReport ("Time to rasterize phantom");
 
     TimerCollectiveMPI timerGather (mpiWorld.getComm());
-    mpi_gather_image (mpiWorld, imGlobal, imLocal, opt_debug);
-    if (opt_verbose)
+    mpi_gather_image (mpiWorld, pImGlobal, pImLocal, optDebug);
+    if (optVerbose)
       timerGather.timerEndAndReport ("Time to gather image");
   }
 #else
-  v = imGlobal->getArray ();
+  v = pImGlobal->getArray ();
   if (phm.getComposition() == P_UNIT_PULSE) {
     v[opt_nx/2][opt_ny/2] = 1.;
   } else if (optFilterName != "") {
-    imGlobal->filterResponse (optDomainName.c_str(), opt_filter_bw, optFilterName.c_str(), opt_filter_param);
+    pImGlobal->filterResponse (optDomainName.c_str(), optFilterBW, optFilterName.c_str(), optFilterParam);
   } else {
 #if HAVE_SGP
-      if (opt_trace >= TRACE_PHM)
+      if (optTrace >= Trace::TRACE_PHANTOM)
        phm.show();
 #endif
-      phm.convertToImagefile (*imGlobal, opt_nsample, opt_trace);
+      phm.convertToImagefile (*pImGlobal, opt_nsample, optTrace);
   }
 #endif
 
@@ -347,46 +345,52 @@ phm2if_main (int argc, char* argv[])
 #endif
   {
     double calctime = timerProgram.timerEnd ();
-    imGlobal->labelAdd (Array2dFileLabel::L_HISTORY, opt_desc.c_str(), calctime);
-    imGlobal->fileWrite (opt_outfile);
-    if (opt_verbose)
-      cout << "Time to rasterized phantom: " << calctime << " seconds" << endl;
+    pImGlobal->labelAdd (Array2dFileLabel::L_HISTORY, optDesc.c_str(), calctime);
+    pImGlobal->fileWrite (optOutFilename.c_str());
+    if (optVerbose)
+      std::cout << "Time to rasterized phantom: " << calctime << " seconds\n";
 
-    if (opt_trace >= TRACE_PHM) {
+    if (optTrace >= Trace::TRACE_PHANTOM) {
       double dmin, dmax;
       int nscale;
       
-      printf ("Enter display size scale (nominal = 1): ");
-      scanf ("%d", &nscale);
-      printf ("Enter minimum and maximum densities (min, max): ");
-      scanf ("%lf %lf", &dmin, &dmax);
-      imGlobal->displayScaling (nscale, dmin, dmax);
+      std::cout << "Enter display size scale (nominal = 1): ";
+      std::cin >> nscale;
+      std::cout ("Enter minimum and maximum densities (min, max): ";
+      std::cin >> dmin;\r
+      std::cin >> dmax;
+      pImGlobal->displayScaling (nscale, dmin, dmax);
     }
   }
 
+  delete pImGlobal;
+#ifdef HAVE_MPI
+  delete pImLocal;
+#endif
+
   return (0);
 }
 
 
 
 #ifdef HAVE_MPI
-void mpi_gather_image (MPIWorld& mpiWorld, ImageFile* imGlobal, ImageFile* imLocal, const int opt_debug)
+void mpi_gather_image (MPIWorld& mpiWorld, ImageFile* pImGlobal, ImageFile* pImLocal, const int optDebug)
 {
-  ImageFileArray vLocal = imLocal->getArray();
+  ImageFileArray vLocal = pImLocal->getArray();
   ImageFileArray vGlobal = NULL;
-  int nyLocal = imLocal->ny();
+  int nyLocal = pImLocal->ny();
 
   if (mpiWorld.getRank() == 0)
-    vGlobal = imGlobal->getArray();
+    vGlobal = pImGlobal->getArray();
   
   for (int iw = 0; iw < mpiWorld.getMyLocalWorkUnits(); iw++)
-    mpiWorld.getComm().Send(vLocal[iw], nyLocal, imLocal->getMPIDataType(), 0, 0);
+    mpiWorld.getComm().Send(vLocal[iw], nyLocal, pImLocal->getMPIDataType(), 0, 0);
 
   if (mpiWorld.getRank() == 0) {
     for (int iProc = 0; iProc < mpiWorld.getNumProcessors(); iProc++) {
       for (int iw = mpiWorld.getStartWorkUnit(iProc); iw <= mpiWorld.getEndWorkUnit(iProc); iw++) {
        MPI::Status status;
-       mpiWorld.getComm().Recv(vGlobal[iw], nyLocal, imLocal->getMPIDataType(), iProc, 0, status);
+       mpiWorld.getComm().Recv(vGlobal[iw], nyLocal, pImLocal->getMPIDataType(), iProc, 0, status);
       }
     }
   }
@@ -403,9 +407,9 @@ main (int argc, char* argv[])
   try {
     retval = phm2if_main(argc, argv);
   } catch (exception e) {
-    cerr << "Exception: " << e.what() << endl;
+    std::cerr << "Exception: " << e.what() << std::endl;
   } catch (...) {
-    cerr << "Unknown exception" << endl;
+    std::cerr << "Unknown exception\n";
   }
 
   return (retval);