r13: *** empty log message ***
[ctsim.git] / src / ctrec.c
index 56297d2af3dd494c7ba419ce48e0fc85e9b70cda..3c9138724fecc79cb2879a4ba56af41a281b64b3 100644 (file)
@@ -2,10 +2,14 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
-**  $Id: ctrec.c,v 1.1 2000/04/28 13:02:44 kevin Exp $
+**  $Id: ctrec.c,v 1.2 2000/04/29 23:24:56 kevin Exp $
 **  $Log: ctrec.c,v $
-**  Revision 1.1  2000/04/28 13:02:44  kevin
-**  Initial revision
+**  Revision 1.2  2000/04/29 23:24:56  kevin
+**  *** empty log message ***
+**
+**  Revision 1.1.1.1  2000/04/28 13:02:44  kevin
+**  Initial CVS import for first public release
+**
 **
 **
 **  This program is free software; you can redistribute it and/or modify
@@ -21,6 +25,7 @@
 **  along with this program; if not, write to the Free Software
 **  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 ******************************************************************************/
+
 /* FILE
  *   ctrec.c           Reconstruct an image from raysums
  *
@@ -104,25 +109,25 @@ usage (const char *program)
 
 
 #ifdef MPI_CT
-void mpi_scatter_rs (RAYSUM *rs_global, RAYSUM *rs_local, const int debug);
+static void mpi_scatter_rs (RAYSUM *rs_global, RAYSUM *rs_local, const int debug);
 #endif
-void print_raysum_info(const RAYSUM *rs);
+static void print_raysum_info(const RAYSUM *rs);
 
 int 
 main (const int argc, char *const argv[])
 {
-  IMAGE *im_global;
-  RAYSUM *rs_global;
-  char *rs_name, *im_filename;
+  IMAGE *im_global = NULL;
+  RAYSUM *rs_global = NULL;
+  char *rs_name, *im_filename = NULL;
   char remark[MAXREMARK];
   char filt_name[80];
   int nx, ny;
-  double time_start, time_end;
+  double time_start = 0, time_end = 0;
   char *endptr;
   int opt_verbose = 0;
   int opt_debug = 0;
   int opt_trace = TRACE_NONE;
-  double opt_filter_param;
+  double opt_filter_param = 1;
   int opt_filter = W_A_BANDLIMIT;
   int opt_interp = I_LINEAR;
   int opt_interp_param = 1;
@@ -239,10 +244,6 @@ main (const int argc, char *const argv[])
 
 #ifdef MPI_CT
   if (mpi_ct.my_rank == 0) {
-    if (file_exists(rs_name) == FALSE) {
-      fprintf (stderr, "File %s does not exist\n", rs_name);
-      exit(1);
-    } 
     rs_global = raysum_open (rs_name);
     raysum_read (rs_global);
     print_raysum_info(rs_global);
@@ -259,7 +260,6 @@ main (const int argc, char *const argv[])
     mpi_rotinc = rs_global->rot_inc;
   }
 
-
   mpi_t1 = MPI_Wtime();
   MPI_Bcast(&opt_verbose, 1, MPI_INT, 0, mpi_ct.comm);
   MPI_Bcast(&opt_debug, 1, MPI_INT, 0, mpi_ct.comm);
@@ -286,17 +286,24 @@ main (const int argc, char *const argv[])
 
   mpi_ct_calc_work_units(mpi_nview);
 
-  rs_local = raysum_create (NULL, mpi_ct.local_work_units[mpi_ct.my_rank], mpi_ndet, TRUE);
-  //  rs_local->ndet = mpi_ndet;
-  //  rs_local->nview = mpi_nview;
+  if (opt_debug) {
+    fprintf(stdout, "Calc'd local work units process %d: nviews=%d, local_work_units=%d, start_work_units=%d\n",
+           mpi_ct.my_rank, mpi_nview, mpi_ct.local_work_units[mpi_ct.my_rank], mpi_ct.start_work_unit[mpi_ct.my_rank]);
+    MPI_Barrier(mpi_ct.comm);
+  }
+
+  rs_local = raysum_create (NULL, mpi_ct.local_work_units[mpi_ct.my_rank], mpi_ndet);
+  if (opt_debug) {
+    fprintf(stdout, "Created rs_local %lx for process %d: local views=%4d, local dets=%4d\n", (unsigned long int) rs_local, mpi_ct.my_rank, mpi_ct.local_work_units[mpi_ct.my_rank], mpi_ndet);
+    MPI_Barrier(mpi_ct.comm);
+  }
+
+  rs_local->ndet = mpi_ndet;
+  rs_local->nview = mpi_nview;
   rs_local->det_inc = mpi_detinc;
   rs_local->piclen = mpi_piclen;
   rs_local->rot_inc = mpi_rotinc;
 
-  if (opt_debug)
-    fprintf(stderr, "Bcast'd rs vars, my nview=%d, local_work_units=%d, start_work_units=%d (process %d)\n",
-           mpi_nview, mpi_ct.local_work_units[mpi_ct.my_rank], mpi_ct.start_work_unit[mpi_ct.my_rank], mpi_ct.my_rank);
-
   if (opt_verbose)
     mpi_t1 = MPI_Wtime();
   mpi_scatter_rs(rs_global, rs_local, opt_debug);
@@ -387,11 +394,11 @@ main (const int argc, char *const argv[])
 
 
 #ifdef MPI_CT
-void mpi_scatter_rs (RAYSUM *rs_global, RAYSUM *rs_local, const int opt_debug)
+static void mpi_scatter_rs (RAYSUM *rs_global, RAYSUM *rs_local, const int opt_debug)
 {
   int iproc;
   int end_work_unit;
-  int iw;
+  int iw = 0;
 
   if (mpi_ct.my_rank == 0) {
     for (iproc = 0; iproc < mpi_ct.nproc; iproc++) {
@@ -403,27 +410,48 @@ void mpi_scatter_rs (RAYSUM *rs_global, RAYSUM *rs_local, const int opt_debug)
        fprintf(stdout, "Sending rs data to process %d\n", iproc);
 
       for (iw = mpi_ct.start_work_unit[iproc]; iw <= end_work_unit; iw++) {
-       MPI_Send(&rs_global->view[iw]->view_angle, 1, MPI_DOUBLE, iproc, 0, mpi_ct.comm);
+        if (opt_debug)
+         fprintf(stdout, "Sending from process %2d to process %2d: view_angle=%8f, ndet=%5d\n", mpi_ct.my_rank, iproc, rs_global->view[iw]->view_angle, rs_global->view[iw]->ndet);
        MPI_Send(&rs_global->view[iw]->ndet, 1, MPI_INT, iproc, 0, mpi_ct.comm);
+       MPI_Send(&rs_global->view[iw]->view_angle, 1, MPI_DOUBLE, iproc, 0, mpi_ct.comm);
        MPI_Send(rs_global->view[iw]->detval, rs_global->ndet, MPI_FLOAT, iproc, 0, mpi_ct.comm);
       }
     }
   }
 
+  if (opt_debug) {
+    fprintf(stdout, "Receiving rs data in process %d\n", mpi_ct.my_rank);
+    MPI_Barrier(mpi_ct.comm);
+  }
+
   end_work_unit = mpi_ct.local_work_units[mpi_ct.my_rank] - 1;
   for (iw = 0; iw <= end_work_unit; iw++) {
     MPI_Status status;
 
-    MPI_Recv(&rs_local->view[iw]->view_angle, 1, MPI_DOUBLE, 0, 0, mpi_ct.comm, &status);
-    MPI_Recv(&rs_local->view[iw]->ndet, 1, MPI_INT, 0, 0, mpi_ct.comm, &status);
-    MPI_Recv(rs_local->view[iw]->detval, rs_global->ndet, MPI_FLOAT, 0, 0, mpi_ct.comm, &status);
+    if (opt_debug) {
+      fprintf(stdout,"Receiving into rs_local from mpi_scatter_rs, process %2d: rs_local=%lx, ", mpi_ct.my_rank, (unsigned long int) rs_local);
+      fprintf(stdout,"iw=%4d, nrot=%4d, ndet=%4d, view=%lx, detval=%lx\n", iw, rs_local->nview, rs_local->ndet, (unsigned long int) rs_local->view[iw], (unsigned long int) rs_local->view[iw]->detval);
+    }
+
+    // MPI_Recv(&rs_local->view[iw]->ndet, 1, MPI_INT, 0, 0, mpi_ct.comm, &status);
+    // MPI_Recv(&rs_local->view[iw]->view_angle, 1, MPI_DOUBLE, 0, 0, mpi_ct.comm, &status);
+    //    MPI_Recv(rs_local->view[iw]->detval, rs_global->ndet, MPI_FLOAT, 0, 0, mpi_ct.comm, &status);
+
+    if (opt_debug) {
+      fprintf(stdout, "Received view_angle=%8f, ndet=%5d\n", rs_local->view[iw]->view_angle, rs_local->view[iw]->ndet);
+    }
   }
   rs_local->nview = mpi_ct.local_work_units[mpi_ct.my_rank];
+  if (opt_debug) {
+    MPI_Barrier(MPI_COMM_WORLD);
+    fprintf(stdout, "Done with mpi_scatter_rs in process %2d\n", mpi_ct.my_rank);
+    exit(0);
+  }
 }
 
 #endif
 
-void print_raysum_info(const RAYSUM *rs)
+static void print_raysum_info(const RAYSUM *rs)
 {
   printf ("Number of detectors: %d\n", rs->ndet);
   printf ("    Number of views: %d\n", rs->nview);