diff --git a/simulation/g4simulation/g4cemc/RawTower.h b/simulation/g4simulation/g4cemc/RawTower.h
index 5a20c0caee..213e01f17c 100644
--- a/simulation/g4simulation/g4cemc/RawTower.h
+++ b/simulation/g4simulation/g4cemc/RawTower.h
@@ -20,6 +20,8 @@ class RawTower : public PHObject {
   virtual int get_bineta() const { PHOOL_VIRTUAL_WARN("get_ieta()"); return -1; }
   virtual int get_binphi() const { PHOOL_VIRTUAL_WARN("get_iphi()");; return -1; }
   virtual float get_energy() const { PHOOL_VIRTUAL_WARN("get_energy()"); return 0.0; }
+  virtual void set_light_yield(float l)  { PHOOL_VIRTUAL_WARN("set_light_yield()"); return ; }
+  virtual float get_light_yield() const { PHOOL_VIRTUAL_WARN("get_light_yield()"); return 0.0; }
 
   virtual bool is_adjacent(RawTower& tower) { PHOOL_VIRTUAL_WARNING; return false; }
 
diff --git a/simulation/g4simulation/g4cemc/RawTowerBuilder.cc b/simulation/g4simulation/g4cemc/RawTowerBuilder.cc
index d592077d16..0c766afae9 100644
--- a/simulation/g4simulation/g4cemc/RawTowerBuilder.cc
+++ b/simulation/g4simulation/g4cemc/RawTowerBuilder.cc
@@ -102,6 +102,7 @@ RawTowerBuilder::process_event(PHCompositeNode *topNode)
 
       if (verbosity > 2)
         {
+          std::cout << PHWHERE << " print out the cell:" << std::endl;
           cell->identify();
         }
 
@@ -113,6 +114,7 @@ RawTowerBuilder::process_event(PHCompositeNode *topNode)
 	   _towers->AddTower(cell->get_binz(), cell->get_binphi(), tower);
 	 }
        tower->add_ecell(cell->get_cell_id(), cell->get_edep());
+       tower->set_light_yield( tower->get_light_yield() +  cell->get_light_yield()  );
        rawtowergeom =  findNode::getClass<RawTowerGeom>(topNode, TowerGeomNodeName.c_str());
        if (verbosity > 2)
 	 {
@@ -224,6 +226,16 @@ RawTowerBuilder::CreateNodes(PHCompositeNode *topNode)
             {
               _netabins = cellgeo->get_zbins();// bin eta in the same number of z bins
             }
+          else if (_cell_binning == PHG4CylinderCellDefs::spacalbinning)
+            {
+              // use eta definiton for each row of towers
+              _netabins = cellgeo->get_etabins();
+            }
+          else
+            {
+              cout <<"RawTowerBuilder::CreateNodes::"<<Name()
+                  <<" - Fatal Error - unsupported cell binning method "<<_cell_binning<<endl;
+            }
 	  inner_radius = cellgeo->get_radius();
 	  first_layer = cellgeo->get_layer();
 	  ifirst = 0;
@@ -322,7 +334,9 @@ RawTowerBuilder::CreateNodes(PHCompositeNode *topNode)
     }
 
   if (_cell_binning == PHG4CylinderCellDefs::etaphibinning
-      or _cell_binning == PHG4CylinderCellDefs::etaslatbinning)
+      or _cell_binning == PHG4CylinderCellDefs::etaslatbinning
+      or _cell_binning == PHG4CylinderCellDefs::spacalbinning
+      )
     {
 //  rawtowergeom->set_etamin(_etamin);
 //  rawtowergeom->set_etastep(_etastep);
diff --git a/simulation/g4simulation/g4cemc/RawTowerGeomv2.cc b/simulation/g4simulation/g4cemc/RawTowerGeomv2.cc
index cc39d24abb..3253c6da96 100644
--- a/simulation/g4simulation/g4cemc/RawTowerGeomv2.cc
+++ b/simulation/g4simulation/g4cemc/RawTowerGeomv2.cc
@@ -89,17 +89,48 @@ RawTowerGeomv2::get_etabin(const double eta) const
 
   int ibin = -1;
   int i = 0;
+//  for (bound_map_t::const_iterator iter = eta_bound_map.begin();
+//      iter != eta_bound_map.end(); ++iter)
+//    {
+//      if (eta >= iter->first && eta < iter->second)
+//        {
+//          ibin = i;
+//          break;
+//        }
+//
+//      i++;
+//    }
+
+  // switch to search for the closest bin
+  // since in a realistic calorimeter, there could be gaps
+  double min_deta = 10;
+
   for (bound_map_t::const_iterator iter = eta_bound_map.begin();
       iter != eta_bound_map.end(); ++iter)
     {
+      const double mean_eta = 0.5*(iter->first + iter->second);
+
       if (eta >= iter->first && eta < iter->second)
         {
+          // found the bin that the hit belong
+          min_deta = 0;
           ibin = i;
           break;
         }
+      else
+        {
+          const double deta = fabs(mean_eta - eta);
+          if (deta < min_deta)
+            {
+              min_deta = deta;
+              ibin = i;
+            }// keep searching
+        }
 
       i++;
     }
+
+
   if (ibin < 0 || ibin >= nphibins)
     {
       cout
diff --git a/simulation/g4simulation/g4cemc/RawTowerv1.cc b/simulation/g4simulation/g4cemc/RawTowerv1.cc
index e861f7b50f..eb007684a1 100644
--- a/simulation/g4simulation/g4cemc/RawTowerv1.cc
+++ b/simulation/g4simulation/g4cemc/RawTowerv1.cc
@@ -15,7 +15,8 @@ RawTowerv1::RawTowerv1() :
 
 RawTowerv1::RawTowerv1(const int ieta, const int iphi) :
   bineta(ieta),
-  binphi(iphi)
+  binphi(iphi),
+  light_yield(0)
 {}
 
 RawTowerv1::~RawTowerv1()
diff --git a/simulation/g4simulation/g4cemc/RawTowerv1.h b/simulation/g4simulation/g4cemc/RawTowerv1.h
index 1fc0ab7efd..935df3a83c 100644
--- a/simulation/g4simulation/g4cemc/RawTowerv1.h
+++ b/simulation/g4simulation/g4cemc/RawTowerv1.h
@@ -21,6 +21,9 @@ class RawTowerv1 : public RawTower {
   int get_binphi() const { return binphi; }
   float get_energy() const;
 
+  void set_light_yield(float l)  { light_yield = l; }
+  float get_light_yield() const { return light_yield; };
+
   bool is_adjacent(RawTower& tower);
 
   std::pair< std::map<unsigned int,float>::const_iterator, std::map<unsigned int,float>::const_iterator > get_g4cells()
@@ -30,10 +33,11 @@ class RawTowerv1 : public RawTower {
  protected:
   int bineta;
   int binphi;
+  float light_yield;
 
   std::map<unsigned int, float> ecells;
 
-  ClassDef(RawTowerv1,1)
+  ClassDef(RawTowerv1,2)
 };
 
 #endif /* RAWTOWERV1_H_ */
diff --git a/simulation/g4simulation/g4detectors/Makefile.am b/simulation/g4simulation/g4detectors/Makefile.am
index 922951b09f..69625555af 100644
--- a/simulation/g4simulation/g4detectors/Makefile.am
+++ b/simulation/g4simulation/g4detectors/Makefile.am
@@ -33,6 +33,7 @@ pkginclude_HEADERS = \
   PHG4BlockGeomContainer.h \
   PHG4CylinderCell.h \
   PHG4CylinderCellv1.h \
+  PHG4CylinderCell_Spacalv1.h \
   PHG4CylinderCellv2.h \
   PHG4CylinderCellContainer.h \
   PHG4CylinderGeom.h \
@@ -42,6 +43,7 @@ pkginclude_HEADERS = \
   PHG4CylinderGeomContainer.h \
   PHG4CylinderCellDefs.h \
   PHG4CylinderCellGeom.h \
+  PHG4CylinderCellGeom_Spacalv1.h \
   PHG4CylinderCellGeomContainer.h
 
 libg4detectors_io_la_SOURCES = \
@@ -77,12 +79,16 @@ libg4detectors_io_la_SOURCES = \
   PHG4CylinderCell_Dict.cc \
   PHG4CylinderCellv1.cc \
   PHG4CylinderCellv1_Dict.cc \
+  PHG4CylinderCell_Spacalv1.cc \
+  PHG4CylinderCell_Spacalv1_Dict.cc \
   PHG4CylinderCellv2.cc \
   PHG4CylinderCellv2_Dict.cc \
   PHG4CylinderCellContainer.cc \
   PHG4CylinderCellContainer_Dict.cc \
   PHG4CylinderCellGeom.cc \
   PHG4CylinderCellGeom_Dict.cc \
+  PHG4CylinderCellGeom_Spacalv1.cc \
+  PHG4CylinderCellGeom_Spacalv1_Dict.cc \
   PHG4CylinderCellGeomContainer.cc \
   PHG4CylinderCellGeomContainer_Dict.cc
 
@@ -184,6 +190,8 @@ libg4detectors_la_SOURCES = \
   PHG4SpacalDetector.cc \
   PHG4ProjSpacalDetector.cc \
   PHG4FullProjSpacalDetector.cc \
+  PHG4FullProjSpacalCellReco.cc \
+  PHG4FullProjSpacalCellReco_Dict.cc \
   PHG4SpacalSteppingAction.cc \
   PHG4SpacalSubsystem.cc \
   PHG4SpacalSubsystem_Dict.cc
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCell.h b/simulation/g4simulation/g4detectors/PHG4CylinderCell.h
index bdf955ae68..69694225f1 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderCell.h
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCell.h
@@ -28,19 +28,25 @@ class PHG4CylinderCell : public PHObject
   }
   
   virtual void add_edep(const PHG4HitDefs::keytype g4hitid, const float edep) {return;}
+  virtual void add_edep(const PHG4HitDefs::keytype g4hitid, const float edep, const float light_yield) {return;}
+
   virtual void set_cell_id(const PHG4CylinderCellDefs::keytype id) {return;}
   virtual void set_layer(const unsigned int i) {return;}
- 
+
   virtual double get_edep() const {return NAN;}
   virtual unsigned int get_layer() const {return 0xFFFFFFFF;}
   virtual PHG4CylinderCellDefs::keytype get_cell_id() const {return 0xFFFFFFFF;}
   virtual int get_binz() const {return -1;}
   virtual int get_binphi() const {return -1;}
   virtual int get_bineta() const {return -1;}
+  virtual float  get_light_yield() const  {    return NAN;  }
+  virtual int    get_fiber_ID() const {return -1;}
 
   virtual void set_phibin(const int i) {return;}
   virtual void set_zbin(const int i) {return;}
   virtual void set_etabin(const int i) {return;}
+  virtual void  set_light_yield(float lightYield)  {   return;  }
+  virtual void  set_fiber_ID(int fiberId) {   return;  }
 
   virtual void set_sensor_index(const std::string &si) {return;}
   virtual std::string get_sensor_index() const {return "";}
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCellDefs.h b/simulation/g4simulation/g4detectors/PHG4CylinderCellDefs.h
index d3bbe6de21..2d740ab002 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderCellDefs.h
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCellDefs.h
@@ -6,7 +6,7 @@ namespace PHG4CylinderCellDefs
   typedef unsigned int keytype;
   static int keybits = 8;
   static int cell_idbits = 32 - keybits;
-  enum {undefined = 0, sizebinning = 1, etaphibinning = 2, etaslatbinning = 3};
+  enum {undefined = 0, sizebinning = 1, etaphibinning = 2, etaslatbinning = 3, spacalbinning = 4};
 }
 
 #endif
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom.cc b/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom.cc
index 98d840ae46..9f1b68ba17 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom.cc
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom.cc
@@ -151,7 +151,7 @@ PHG4CylinderCellGeom::set_etabins(const int i)
 void
 PHG4CylinderCellGeom::identify(std::ostream& os) const
 {
-  os << "layer: " << layer
+  os << "PHG4CylinderCellGeom::identify - layer: " << layer
      << ", radius: " << radius
      << ", thickness: " << thickness;
   switch (binning)
@@ -168,8 +168,11 @@ PHG4CylinderCellGeom::identify(std::ostream& os) const
       break;
     case PHG4CylinderCellDefs::etaslatbinning:
       os << ", etabins: " << nzbins
-	 << ", etamin: " << zmin
-	 << ", etastepsize: " << zstep;
+   << ", etamin: " << zmin
+   << ", etastepsize: " << zstep;
+      break;
+    case PHG4CylinderCellDefs::spacalbinning:
+      os << ", etabins: " << nzbins<<" for Spacal";
       break;
     default:
       os << "no valid binning method: " << binning << endl;
@@ -317,6 +320,9 @@ PHG4CylinderCellGeom::methodname(const int i) const
     case PHG4CylinderCellDefs::etaslatbinning:
       return "Eta/numslat bins";
       break;
+    case PHG4CylinderCellDefs::spacalbinning:
+      return "SPACAL Tower bins";
+      break;
     default:
       break;
     }
@@ -340,7 +346,8 @@ void
 PHG4CylinderCellGeom::check_binning_method_eta(const std::string & src) const
 {
   if (binning != PHG4CylinderCellDefs::etaphibinning && 
-      binning != PHG4CylinderCellDefs::etaslatbinning)
+      binning != PHG4CylinderCellDefs::etaslatbinning&&
+      binning != PHG4CylinderCellDefs::spacalbinning)
     {
       if (src.size())
         cout << src<<" : ";
@@ -348,6 +355,7 @@ PHG4CylinderCellGeom::check_binning_method_eta(const std::string & src) const
       cout << "different binning method used " << methodname(binning)
            << ", not : " << methodname(PHG4CylinderCellDefs::etaphibinning)
 	   << " or " << methodname(PHG4CylinderCellDefs::etaslatbinning)
+     << " or " << methodname(PHG4CylinderCellDefs::spacalbinning)
            << endl;
       exit(1);
     }
@@ -359,14 +367,17 @@ PHG4CylinderCellGeom::check_binning_method_phi(const std::string & src) const
 {
   if (binning != PHG4CylinderCellDefs::etaphibinning && 
       binning != PHG4CylinderCellDefs::sizebinning &&
-      binning != PHG4CylinderCellDefs::etaslatbinning)
+      binning != PHG4CylinderCellDefs::etaslatbinning &&
+      binning != PHG4CylinderCellDefs::spacalbinning)
     {
       if (src.size())
         cout << src<<" : ";
 
       cout << "different binning method used " << methodname(binning)
            << ", not : " << methodname(PHG4CylinderCellDefs::etaphibinning)
-	   << " or " << methodname(PHG4CylinderCellDefs::sizebinning)
+           << " or " << methodname(PHG4CylinderCellDefs::sizebinning)
+           << " or " << methodname(PHG4CylinderCellDefs::etaslatbinning)
+           << " or " << methodname(PHG4CylinderCellDefs::spacalbinning)
            << endl;
       exit(1);
     }
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom.h b/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom.h
index 724dd00aa9..62b9647da3 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom.h
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom.h
@@ -28,16 +28,16 @@ class PHG4CylinderCellGeom: public PHObject
   double get_etastep() const;
   double get_etamin() const;
 
-  std::pair<double, double> get_zbounds(const int ibin) const;
-  std::pair<double, double> get_phibounds(const int ibin) const;
-  std::pair<double, double> get_etabounds(const int ibin) const;
-  double get_etacenter(const int ibin) const;
-  double get_zcenter(const int ibin) const;
-  double get_phicenter(const int ibin) const;
+  virtual std::pair<double, double> get_zbounds(const int ibin) const;
+  virtual std::pair<double, double> get_phibounds(const int ibin) const;
+  virtual std::pair<double, double> get_etabounds(const int ibin) const;
+  virtual double get_etacenter(const int ibin) const;
+  virtual double get_zcenter(const int ibin) const;
+  virtual double get_phicenter(const int ibin) const;
 
-  int get_etabin(const double eta) const;
-  int get_zbin(const double z) const;
-  int get_phibin(const double phi) const;
+  virtual int get_etabin(const double eta) const;
+  virtual int get_zbin(const double z) const;
+  virtual int get_phibin(const double phi) const;
 
    void set_layer(const int i) {layer = i;}
    void set_binning(const int i) {binning = i;}
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom_Spacalv1.cc b/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom_Spacalv1.cc
new file mode 100644
index 0000000000..f009fe42a6
--- /dev/null
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom_Spacalv1.cc
@@ -0,0 +1,192 @@
+// $Id: $                                                                                             
+
+/*!
+ * \file PHG4CylinderCellGeomSpacalv1.cc
+ * \brief 
+ * \author Jin Huang <jhuang@bnl.gov>
+ * \version $Revision:   $
+ * \date $Date: $
+ */
+
+#include "PHG4CylinderCellGeom_Spacalv1.h"
+#include "PHG4CylinderCellDefs.h"
+
+#include <stdexcept>
+#include <sstream>
+#include <cassert>
+#include <boost/foreach.hpp>
+#include <cmath>
+#include <cstdlib>
+#include <iostream>
+#include <sstream>
+ClassImp(PHG4CylinderCellGeom_Spacalv1);
+
+using namespace std;
+
+PHG4CylinderCellGeom_Spacalv1::PHG4CylinderCellGeom_Spacalv1()
+{
+
+}
+
+PHG4CylinderCellGeom_Spacalv1::~PHG4CylinderCellGeom_Spacalv1()
+{
+  // TODO Auto-generated destructor stub
+}
+
+void
+PHG4CylinderCellGeom_Spacalv1::identify(std::ostream& os) const
+{
+  PHG4CylinderCellGeom::identify(os);
+
+  cout << "PHG4CylinderCellGeom_Spacalv1::identify - Tower mapping:" << endl;
+  BOOST_FOREACH(const tower_z_ID_eta_bin_map_t::value_type& tower_z_ID_eta_bin,
+      get_tower_z_ID_eta_bin_map())
+    {
+      cout << "\t" << "Tower Z ID[" << tower_z_ID_eta_bin.first
+          << "] \t-> Eta Bin " << tower_z_ID_eta_bin.second << endl;
+    }
+
+  cout << "PHG4CylinderCellGeom_Spacalv1::identify - Bin -> z range:" << endl;
+  BOOST_FOREACH(const bound_map_t::value_type& b, z_bound_map)
+    {
+      cout << "\t" << "bin[" << b.first << "] \t-> z = " << b.second.first
+          << " - " << b.second.second << endl;
+    }
+
+  cout << "PHG4CylinderCellGeom_Spacalv1::identify - Bin -> eta range:" << endl;
+  BOOST_FOREACH(const bound_map_t::value_type& b, eta_bound_map)
+    {
+      cout << "\t" << "bin[" << b.first << "] \t-> eta = " << b.second.first
+          << " - " << b.second.second << endl;
+    }
+  return;
+}
+
+void
+PHG4CylinderCellGeom_Spacalv1::map_consistency_check() const
+{
+  if ((size_t) nzbins != z_bound_map.size())
+    {
+      cout << "PHG4CylinderCellGeom_Spacalv1::map_consistency_check - "
+          << "z_bound_map.size() of " << z_bound_map.size()
+          << " in inconsistent with nzbins of " << nzbins << endl;
+      exit(1);
+    }
+  if ((size_t) nzbins != eta_bound_map.size())
+    {
+      cout << "PHG4CylinderCellGeom_Spacalv1::map_consistency_check - "
+          << "eta_bound_map.size() of " << eta_bound_map.size()
+          << " in inconsistent with nzbins of " << nzbins << endl;
+      exit(1);
+    }
+  if ((size_t) nzbins != tower_z_ID_eta_bin_map.size())
+    {
+      cout << "PHG4CylinderCellGeom_Spacalv1::map_consistency_check - "
+          << "tower_z_ID_eta_bin_map.size() of " << tower_z_ID_eta_bin_map.size()
+          << " in inconsistent with nzbins of " << nzbins << endl;
+      exit(1);
+    }
+}
+
+void
+PHG4CylinderCellGeom_Spacalv1::set_zbounds(const int ibin,
+    const std::pair<double, double> & bounds)
+{
+  assert(ibin>=0);
+  z_bound_map[ibin] = bounds;
+}
+
+void
+PHG4CylinderCellGeom_Spacalv1::set_etabounds(const int ibin,
+    const std::pair<double, double> & bounds)
+{
+  assert(ibin>=0);
+  eta_bound_map[ibin] = bounds;
+}
+
+pair<double, double>
+PHG4CylinderCellGeom_Spacalv1::get_zbounds(const int ibin) const
+{
+  map_consistency_check();
+  check_binning_method(PHG4CylinderCellDefs::spacalbinning);
+  bound_map_t ::const_iterator iter =
+  z_bound_map.find(ibin);
+
+  if (iter == z_bound_map.end())
+    {
+      cout
+          << "PHG4CylinderCellGeom_Spacalv1::get_zbounds - Fatal Error - Asking for invalid bin in z: "
+          << ibin<<". Print of content:" << endl;
+      identify();
+      exit(1);
+    }
+  return iter->second;
+}
+
+pair<double, double>
+PHG4CylinderCellGeom_Spacalv1::get_etabounds(const int ibin) const
+{
+  map_consistency_check();
+  check_binning_method(PHG4CylinderCellDefs::spacalbinning);
+
+  bound_map_t ::const_iterator iter =
+  eta_bound_map.find(ibin);
+
+  if (iter == eta_bound_map.end())
+    {
+      cout
+          << "PHG4CylinderCellGeom_Spacalv1::get_etabounds - Fatal Error - Asking for invalid bin in z: "
+          << ibin<<". Print of content:" << endl;
+      identify();
+      exit(1);
+    }
+  return iter->second;
+}
+
+int
+PHG4CylinderCellGeom_Spacalv1::get_zbin(const double z) const
+{
+  cout << "PHG4CylinderCellGeom_Spacalv1::get_zbin is invalid" << endl;
+  exit(1);
+  return -1;
+}
+
+int
+PHG4CylinderCellGeom_Spacalv1::get_etabin(const double eta) const
+{
+  cout << "PHG4CylinderCellGeom_Spacalv1::get_etabin is invalid" << endl;
+  exit(1);
+  return -1;
+}
+
+double
+PHG4CylinderCellGeom_Spacalv1::get_zcenter(const int ibin) const
+{
+  pair<double, double> bound = get_zbounds(ibin);
+  return 0.5 * (bound.first + bound.second);
+}
+
+double
+PHG4CylinderCellGeom_Spacalv1::get_etacenter(const int ibin) const
+{
+  pair<double, double> bound = get_etabounds(ibin);
+  return 0.5 * (bound.first + bound.second);
+}
+
+int PHG4CylinderCellGeom_Spacalv1::get_etabin(const int tower_z_ID) const
+{
+  map_consistency_check();
+
+  tower_z_ID_eta_bin_map_t::const_iterator iter = tower_z_ID_eta_bin_map.find(tower_z_ID);
+
+  if (iter == tower_z_ID_eta_bin_map.end())
+    {
+      ostringstream o;
+
+      o <<"PHG4CylinderCellGeom_Spacalv1::get_etabin - Fatal Error - can not find tower_z_ID of "<<tower_z_ID<<".";
+
+      throw range_error(o.str());
+    }
+
+  return iter->second;
+}
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom_Spacalv1.h b/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom_Spacalv1.h
new file mode 100644
index 0000000000..671a7a09b1
--- /dev/null
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom_Spacalv1.h
@@ -0,0 +1,112 @@
+// $Id: $                                                                                             
+
+/*!
+ * \file PHG4CylinderCellGeomSpacalv1.h
+ * \brief 
+ * \author Jin Huang <jhuang@bnl.gov>
+ * \version $Revision:   $
+ * \date $Date: $
+ */
+
+#ifndef PHG4CYLINDERCELLGEOMSPACALV1_H_
+#define PHG4CYLINDERCELLGEOMSPACALV1_H_
+
+#include <PHG4CylinderCellGeom.h>
+#include <map>
+
+/*!
+ * \brief PHG4CylinderCellGeom_Spacalv1
+ */
+class PHG4CylinderCellGeom_Spacalv1 : public PHG4CylinderCellGeom
+{
+public:
+  PHG4CylinderCellGeom_Spacalv1();
+  virtual
+  ~PHG4CylinderCellGeom_Spacalv1();
+  void
+  identify(std::ostream& os = std::cout) const;
+
+  virtual std::pair<double, double>
+  get_zbounds(const int ibin) const;
+  virtual std::pair<double, double>
+  get_etabounds(const int ibin) const;
+
+  virtual double
+  get_etacenter(const int ibin) const;
+  virtual double
+  get_zcenter(const int ibin) const;
+
+  virtual int
+  get_etabin(const double eta) const;
+  virtual int
+  get_zbin(const double z) const;
+
+  void
+  set_zbounds(const int ibin, const std::pair<double, double> & bounds);
+  void
+  set_etabounds(const int ibin, const std::pair<double, double> & bounds);
+
+  typedef std::pair<double, double> bound_t;
+  typedef std::map<int, bound_t> bound_map_t;
+
+  const bound_map_t &
+  get_eta_bound_map() const
+  {
+    map_consistency_check();
+    return eta_bound_map;
+  }
+
+  void
+  set_eta_bound_map(const bound_map_t & etaBoundMap)
+  {
+    eta_bound_map = etaBoundMap;
+  }
+
+  const bound_map_t &
+  get_z_bound_map() const
+  {
+    map_consistency_check();
+    return z_bound_map;
+  }
+
+  void
+  set_z_bound_map(const bound_map_t & boundMap)
+  {
+    z_bound_map = boundMap;
+  }
+
+  //! map tower_z_ID -> eta_bin number
+  typedef std::map<int, int> tower_z_ID_eta_bin_map_t;
+
+  //! map tower_z_ID -> eta_bin number
+  const tower_z_ID_eta_bin_map_t &
+  get_tower_z_ID_eta_bin_map() const
+  {
+    return tower_z_ID_eta_bin_map;
+  }
+
+  virtual int
+  get_etabin(const int tower_z_ID) const;
+
+  //! map tower_z_ID -> eta_bin number
+  void
+  set_tower_z_ID_eta_bin_map(const tower_z_ID_eta_bin_map_t & m)
+  {
+    tower_z_ID_eta_bin_map = m;
+  }
+
+protected:
+
+  void
+  map_consistency_check() const;
+
+  bound_map_t z_bound_map;
+  bound_map_t eta_bound_map;
+
+  tower_z_ID_eta_bin_map_t tower_z_ID_eta_bin_map;
+
+ClassDef(PHG4CylinderCellGeom_Spacalv1,1)
+
+};
+
+#endif /* PHG4CYLINDERCELLGEOMSPACALV1_H_ */
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom_Spacalv1LinkDef.h b/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom_Spacalv1LinkDef.h
new file mode 100644
index 0000000000..a23d410c97
--- /dev/null
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCellGeom_Spacalv1LinkDef.h
@@ -0,0 +1,5 @@
+#ifdef __CINT__
+
+#pragma link C++ class PHG4CylinderCellGeom_Spacalv1+;
+
+#endif /* __CINT__ */
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCellReco.cc b/simulation/g4simulation/g4detectors/PHG4CylinderCellReco.cc
index 541accab2b..80b1ef8276 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderCellReco.cc
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCellReco.cc
@@ -454,7 +454,7 @@ PHG4CylinderCellReco::process_event(PHCompositeNode *topNode)
                       it->second->set_layer(*layer);
                       it->second->set_phibin(iphibin);
                       it->second->set_etabin(ietabin);		      
-		      it->second->add_edep(hiter->first, hiter->second->get_edep()*vdedx[i1]);
+		      it->second->add_edep(hiter->first, hiter->second->get_edep()*vdedx[i1], hiter->second->get_light_yield());
 		    }
 
 		  // just a sanity check - we don't want to mess up by having Nan's or Infs in our energy deposition
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCell_Spacalv1.cc b/simulation/g4simulation/g4detectors/PHG4CylinderCell_Spacalv1.cc
new file mode 100644
index 0000000000..efd454ef59
--- /dev/null
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCell_Spacalv1.cc
@@ -0,0 +1,32 @@
+// $Id: $                                                                                             
+
+/*!
+ * \file PHG4CylinderCellSpacalv1.cc
+ * \brief 
+ * \author Jin Huang <jhuang@bnl.gov>
+ * \version $Revision:   $
+ * \date $Date: $
+ */
+
+#include "PHG4CylinderCell_Spacalv1.h"
+#include "PHG4CylinderCellDefs.h"
+
+#include <cmath>
+#include <cstdlib>
+
+
+ClassImp(PHG4CylinderCell_Spacalv1);
+
+using namespace std;
+
+PHG4CylinderCell_Spacalv1::PHG4CylinderCell_Spacalv1():
+    fiber_ID(-1)
+{
+  // TODO Auto-generated constructor stub
+
+}
+
+PHG4CylinderCell_Spacalv1::~PHG4CylinderCell_Spacalv1()
+{
+  // TODO Auto-generated destructor stub
+}
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCell_Spacalv1.h b/simulation/g4simulation/g4detectors/PHG4CylinderCell_Spacalv1.h
new file mode 100644
index 0000000000..d20e5f664e
--- /dev/null
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCell_Spacalv1.h
@@ -0,0 +1,50 @@
+// $Id: $                                                                                             
+
+/*!
+ * \file PHG4CylinderCellSpacalv1.h
+ * \brief 
+ * \author Jin Huang <jhuang@bnl.gov>
+ * \version $Revision:   $
+ * \date $Date: $
+ */
+
+#ifndef PHG4CYLINDERCELLSPACALV1_H_
+#define PHG4CYLINDERCELLSPACALV1_H_
+
+#include <PHG4CylinderCellv1.h>
+/*!
+ * \brief PHG4CylinderCell_Spacalv1
+ */
+class PHG4CylinderCell_Spacalv1 : public PHG4CylinderCellv1
+{
+public:
+  PHG4CylinderCell_Spacalv1();
+  virtual
+  ~PHG4CylinderCell_Spacalv1();
+
+  //! Group hit into each fiber, so we allow fiber-fiber variation studies off-Geant production.
+  int
+  get_fiber_ID() const
+  {
+    return fiber_ID;
+  }
+
+  //! Group hit into each fiber, so we allow fiber-fiber variation studies off-Geant production.
+  void
+  set_fiber_ID(int fiberId)
+  {
+    fiber_ID = fiberId;
+  }
+
+
+protected:
+
+  //! Group hit into each fiber, so we allow fiber-fiber variation studies off-Geant production.
+  int fiber_ID;
+
+
+ClassDef(PHG4CylinderCell_Spacalv1,1)
+
+};
+
+#endif /* PHG4CYLINDERCELLSPACALV1_H_ */
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCell_Spacalv1LinkDef.h b/simulation/g4simulation/g4detectors/PHG4CylinderCell_Spacalv1LinkDef.h
new file mode 100644
index 0000000000..de8bea3ea0
--- /dev/null
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCell_Spacalv1LinkDef.h
@@ -0,0 +1,5 @@
+#ifdef __CINT__
+
+#pragma link C++ class PHG4CylinderCell_Spacalv1+;
+
+#endif /* __CINT__ */
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCellv1.cc b/simulation/g4simulation/g4detectors/PHG4CylinderCellv1.cc
index ad80189e50..76c74d0061 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderCellv1.cc
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCellv1.cc
@@ -10,7 +10,7 @@ PHG4CylinderCellv1::PHG4CylinderCellv1():
   cellid(0xFFFFFFFF),
   binz(-1),
   binphi(-1),
-  edeps()
+  edeps(), light_yield(0)
 {}
 
 void
@@ -26,6 +26,13 @@ PHG4CylinderCellv1::add_edep(const PHG4HitDefs::keytype g4hitid, const float ede
     }
 }
 
+void
+PHG4CylinderCellv1::add_edep(const PHG4HitDefs::keytype g4hitid, const float edep, const float ly)
+{
+  add_edep(g4hitid, edep);
+  light_yield += ly;
+}
+
 double PHG4CylinderCellv1::get_edep() const {
   
   double esum = 0.0;
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderCellv1.h b/simulation/g4simulation/g4detectors/PHG4CylinderCellv1.h
index 1ccf3a4e88..b0a69cee9f 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderCellv1.h
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderCellv1.h
@@ -18,7 +18,8 @@ class PHG4CylinderCellv1 : public PHG4CylinderCell
 
    EdepConstRange get_g4hits()
     {return make_pair(edeps.begin(), edeps.end());}
-  void add_edep(const PHG4HitDefs::keytype g4hitid, const float edep);
+   void add_edep(const PHG4HitDefs::keytype g4hitid, const float edep);
+   void add_edep(const PHG4HitDefs::keytype g4hitid, const float edep, const float light_yield);
   void set_cell_id(const PHG4CylinderCellDefs::keytype id) {cellid = id;}
   void set_layer(const unsigned int i) {layer = i;}
   double get_edep() const;
@@ -27,10 +28,13 @@ class PHG4CylinderCellv1 : public PHG4CylinderCell
   int get_binz() const {return binz;}
   int get_binphi() const {return binphi;}
   int get_bineta() const {return get_binz();}
+  float  get_light_yield() const  {    return light_yield;  }
+
 
   void set_zbin(const int i) {binz = i;}
   void set_etabin(const int i) {set_zbin(i);}
   void set_phibin(const int i) {binphi = i;}
+  void  set_light_yield(float lightYield)  {    light_yield = lightYield;  }
 
  protected:
 
@@ -39,8 +43,9 @@ class PHG4CylinderCellv1 : public PHG4CylinderCell
   int binz;
   int binphi;
   EdepMap edeps;
+  float light_yield;
    
-  ClassDef(PHG4CylinderCellv1,1)
+  ClassDef(PHG4CylinderCellv1,2)
 };
 
 #endif
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv1.cc b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv1.cc
index 4083510dd5..a6ff770fdf 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv1.cc
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv1.cc
@@ -43,26 +43,31 @@ PHG4CylinderGeom_Spacalv1::Print(Option_t *) const
 {
   identify(cout);
 
-  cout <<"Configuration is #"<<get_config()<<":"<<endl;
-  switch(get_config())
-  {
+  cout << "Configuration is #" << get_config() << ":" << endl;
+  switch (get_config())
+    {
   case kNonProjective:
-    cout <<"fiber always placed radially"<<endl;
+    cout << "fiber always placed radially" << endl;
     break;
   case kProjective_PolarTaper:
-    cout <<"Block constructed with taper in polar direction, non-taper in azimuthal direction. "
-    <<"The final layout is approximately projective in both azimuthal and polar directions."<<endl;
+    cout
+        << "Block constructed with taper in polar direction, non-taper in azimuthal direction. "
+        << "The final layout is approximately projective in both azimuthal and polar directions."
+        << endl;
     break;
   case kFullProjective_2DTaper:
-    cout <<"Fully projective spacal with 2D tapered modules"<<endl;
+    cout << "Fully projective spacal with 2D tapered modules" << endl;
     break;
   case kFullProjective_2DTaper_SameLengthFiberPerTower:
-    cout <<"Fully projective spacal with 2D tapered modules. To speed up construction, same-length fiber is used cross one tower"<<endl;
+    cout
+        << "Fully projective spacal with 2D tapered modules. To speed up construction, same-length fiber is used cross one tower"
+        << endl;
     break;
   default:
-    cout <<"PHG4CylinderGeom_Spacalv1::Print - ERROR - unknown configuration #"<<get_config()<<endl;
+    cout << "PHG4CylinderGeom_Spacalv1::Print - ERROR - unknown configuration #"
+        << get_config() << endl;
     break;
-  }
+    }
 
   cout << "\t" << "get_max_radius() = " << get_max_radius() << endl;
   cout << "\t" << "get_half_radius() = " << get_half_radius() << endl;
@@ -70,11 +75,11 @@ PHG4CylinderGeom_Spacalv1::Print(Option_t *) const
   cout << "\t" << "get_*pos() = " << get_xpos() << ", " << get_ypos() << ", "
       << get_zpos() << endl;
 
-  cout << "\t" << "get_azimuthal_n_sec() = " << get_azimuthal_n_sec() << endl;
+  cout << "\t" << "get_azimuthal_n_sec() = " << get_azimuthal_n_sec()
+      <<", "<<sector_map.size()<<"/"<< get_azimuthal_n_sec()<<" azimuthal sectors would be filled with SPACAL."<< endl;
   cout << "\t" << "get_azimuthal_distance() = " << get_azimuthal_distance()
       << endl;
-  cout << "\t" << "get_z_distance() = " << get_z_distance()
-      << endl;
+  cout << "\t" << "get_z_distance() = " << get_z_distance() << endl;
   cout << "\t" << "get_fiber_outer_r() = " << get_fiber_outer_r() << endl;
   cout << "\t" << "get_fiber_clading_thickness() = "
       << get_fiber_clading_thickness() << endl;
@@ -116,17 +121,20 @@ PHG4CylinderGeom_Spacalv1::SetDefault()
   zpos = 0;
 
   fiber_clading_thickness = 0.003 / 2;
-  fiber_core_diameter = 0.047 - fiber_clading_thickness*2;
+  fiber_core_diameter = 0.047 - fiber_clading_thickness * 2;
   fiber_distance = 0.1;
 
   virualize_fiber = false;
   construction_verbose = 0;
+
+  init_default_sector_map();
 }
 
 int
 PHG4CylinderGeom_Spacalv1::get_azimuthal_n_sec() const
 {
-  return std::floor(get_half_radius() * twopi / (get_fiber_distance() * sqrt(3.)));
+  return std::floor(
+      get_half_radius() * twopi / (get_fiber_distance() * sqrt(3.)));
 }
 
 double
@@ -138,5 +146,22 @@ PHG4CylinderGeom_Spacalv1::get_azimuthal_distance() const
 double
 PHG4CylinderGeom_Spacalv1::get_z_distance() const
 {
-  return get_fiber_distance()/2.;
+  return get_fiber_distance() / 2.;
+}
+
+//! load a default map that populate all the sectors
+void
+PHG4CylinderGeom_Spacalv1::init_default_sector_map()
+{
+  sector_map.clear();
+
+  for (int sec = 0; sec < get_azimuthal_n_sec(); ++sec)
+    {
+      const double rot = twopi / (double) (get_azimuthal_n_sec())
+          * ((double) (sec));
+
+      sector_map[sec] = rot;
+    }
 }
+
+
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv1.h b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv1.h
index d098179686..e9fc8941dd 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv1.h
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv1.h
@@ -13,6 +13,7 @@
 #include "PHG4CylinderGeomv2.h"
 #include <string>
 #include <cmath>
+#include <map>
 
 class PHG4CylinderGeom_Spacalv1 : public PHG4CylinderGeomv2
 {
@@ -114,6 +115,27 @@ class PHG4CylinderGeom_Spacalv1 : public PHG4CylinderGeomv2
   double
   get_z_distance() const;
 
+  //! sector map sector_ID -> azimuthal rotation.
+  typedef std::map<int,double> sector_map_t;
+
+  //! sector map sector_ID -> azimuthal rotation.
+  const sector_map_t &
+  get_sector_map() const
+  {
+    return sector_map;
+  }
+
+  //! sector map sector_ID -> azimuthal rotation.
+  sector_map_t &
+  get_sector_map()
+  {
+    return sector_map;
+  }
+
+  //! load a default map that populate all the sectors
+  void
+  init_default_sector_map();
+
   ///@}
 
   /** @name Fiber geometry
@@ -307,7 +329,11 @@ class PHG4CylinderGeom_Spacalv1 : public PHG4CylinderGeomv2
   bool virualize_fiber;
   int construction_verbose;
 
-ClassDef(PHG4CylinderGeom_Spacalv1,1)
+  //! sector map sector_ID -> azimuthal rotation.
+  sector_map_t sector_map;
+
+
+ClassDef(PHG4CylinderGeom_Spacalv1,2)
 
 };
 
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv2.cc b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv2.cc
index ae61b137fd..84eb2530b5 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv2.cc
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv2.cc
@@ -74,6 +74,8 @@ PHG4CylinderGeom_Spacalv2::SetDefault()
   polar_taper_ratio = 1 + 1.1 / 42.;
   assembly_spacing = 0.0001; // order ~1um clearance around all structures
 
+  init_default_sector_map();
+
 }
 
 double
@@ -120,6 +122,7 @@ PHG4CylinderGeom_Spacalv2::set_azimuthal_n_sec(int azimuthalNSec)
     }
 
   azimuthal_n_sec = azimuthalNSec;
+  init_default_sector_map();
 }
 
 bool
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3.cc b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3.cc
index 2a70571996..d8a9cdbeda 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3.cc
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3.cc
@@ -14,11 +14,13 @@
 #include <Geant4/G4PhysicalConstants.hh>
 
 #include <cmath>
+#include <cassert>
 #include <limits>       // std::numeric_limits
 #include <iostream>
 
 ClassImp(PHG4CylinderGeom_Spacalv3)
 ClassImp(PHG4CylinderGeom_Spacalv3::geom_tower)
+ClassImp(PHG4CylinderGeom_Spacalv3::scint_id_coder)
 
 using namespace std;
 
@@ -54,6 +56,8 @@ PHG4CylinderGeom_Spacalv3::Print(Option_t *opt) const
   cout << "\t" << "get_sidewall_thickness() = " << get_sidewall_thickness()
       << endl;
   cout << "\t" << "get_sidewall_mat() = " << get_sidewall_mat() << endl;
+  cout << "\t" << "get_max_phi_bin_in_sec() = " << get_max_phi_bin_in_sec()
+      << endl;
   cout << "Containing " << sector_tower_map.size()
       << " unique towers per sector." << endl;
 
@@ -78,7 +82,7 @@ PHG4CylinderGeom_Spacalv3::SetDefault()
   sidewall_thickness = 0.075000;
   sidewall_outer_torr = 0.030000;
   sidewall_mat = "SS310";
-
+  max_phi_bin_in_sec = 8;
 }
 
 PHG4CylinderGeom_Spacalv3::geom_tower::geom_tower() :
@@ -111,6 +115,52 @@ PHG4CylinderGeom_Spacalv3::geom_tower::identify(std::ostream& os) const
       << centralZ << " cm" << endl;
 }
 
+PHG4CylinderGeom_Spacalv3::scint_id_coder::scint_id_coder(int scint_id) :
+    scint_ID(scint_id)
+{
+  sector_ID = (scint_ID >> (kfiber_bit + ktower_bit))
+      & ((1 << ksector_bit) - 1);
+  tower_ID = (scint_ID >> kfiber_bit) & ((1 << ktower_bit) - 1);
+  fiber_ID = (scint_ID) & ((1 << kfiber_bit) - 1);
+}
+
+PHG4CylinderGeom_Spacalv3::scint_id_coder::scint_id_coder(int sector_id,
+    int tower_id, int fiber_id) :
+    sector_ID(sector_id), tower_ID(tower_id), fiber_ID(fiber_id)
+{
+  assert(fiber_ID<(1<<kfiber_bit) and fiber_ID>=0);
+  assert(tower_ID<(1<<ktower_bit) and tower_ID>=0);
+  assert(sector_ID<(1<<ksector_bit) and sector_ID>=0);
+
+  scint_ID = (((sector_ID << ktower_bit) | tower_ID) << kfiber_bit) | fiber_ID;
+}
+
+std::pair<int, int>
+PHG4CylinderGeom_Spacalv3::get_tower_z_phi_ID(const int tower_ID,
+    const int sector_ID) const
+{
+  // tower_ID to eta/z within a sector
+  const int z_bin = floor(tower_ID / 10);
+
+  // colume ID is from -x to +x at the top of the detector, which is reverse of the phi bin direction.
+  const int phi_bin_in_sec = max_phi_bin_in_sec - (tower_ID % 10);
+
+  if (!(phi_bin_in_sec < max_phi_bin_in_sec and phi_bin_in_sec >= 0))
+    {
+      cout
+          << "PHG4CylinderGeom_Spacalv3::get_tower_z_phi_ID - Fatal Error - invalid in put with "
+          << "tower_ID = " << tower_ID << ", sector_ID = " << sector_ID
+          << ". Dump object:" << endl;
+      Print();
+    }
+
+  assert(phi_bin_in_sec < max_phi_bin_in_sec and phi_bin_in_sec>=0);
+
+  const int phi_bin = sector_ID * max_phi_bin_in_sec + phi_bin_in_sec;
+
+  return make_pair<int, int>(z_bin, phi_bin);
+}
+
 void
 PHG4CylinderGeom_Spacalv3::load_demo_sector_tower_map1()
 {
@@ -123,6 +173,10 @@ PHG4CylinderGeom_Spacalv3::load_demo_sector_tower_map1()
   zmin = 149.470000;
   zmax = -zmin;
   azimuthal_n_sec = 32;
+  max_phi_bin_in_sec = 8;
+  sector_map.clear();
+  sector_map[0] = 0; // only install one sector
+
   azimuthal_tilt = 0;
   azimuthal_seg_visible = false;
   polar_taper_ratio = 1.;
@@ -195,6 +249,9 @@ PHG4CylinderGeom_Spacalv3::load_demo_sector_tower_map2()
   zmin = 149.470000;
   zmax = -zmin;
   azimuthal_n_sec = 32;
+  max_phi_bin_in_sec = 8;
+  sector_map.clear();
+  sector_map[0] = 0; // only install one sector
   azimuthal_tilt = 0;
   azimuthal_seg_visible = false;
   polar_taper_ratio = 1.;
@@ -409,6 +466,8 @@ PHG4CylinderGeom_Spacalv3::load_demo_sector_tower_map3()
   zmin = 149.470000;
   zmax = -zmin;
   azimuthal_n_sec = 32;
+  max_phi_bin_in_sec = 8;
+  init_default_sector_map();
   azimuthal_tilt = 0;
   azimuthal_seg_visible = false;
   polar_taper_ratio = 1.;
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3.h b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3.h
index bac1203e54..bb243b7acc 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3.h
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3.h
@@ -13,6 +13,8 @@
 #include "PHG4CylinderGeom_Spacalv2.h"
 #include <string>
 #include <map>
+#include <utility>      // std::pair, std::make_pair
+
 
 class PHG4CylinderGeom_Spacalv3 : public PHG4CylinderGeom_Spacalv2
 {
@@ -68,6 +70,17 @@ class PHG4CylinderGeom_Spacalv3 : public PHG4CylinderGeom_Spacalv2
     sidewall_mat = absorberMat;
   }
 
+  int get_max_phi_bin_in_sec() const
+    {
+      return max_phi_bin_in_sec;
+    }
+
+  void set_max_phi_bin_in_sec(int maxPhiBinInSec)
+    {
+      max_phi_bin_in_sec = maxPhiBinInSec;
+    }
+
+  
 
   class geom_tower
   {
@@ -124,15 +137,55 @@ class PHG4CylinderGeom_Spacalv3 : public PHG4CylinderGeom_Spacalv2
 //  {
 //    geom_super_tower_map = geomSuperTowerMap;
 //  }
+
+//! compact ID of each fiber in 32bit PHG4Hit::set_scint_id(). Buffer the result for repeated use.
+  class scint_id_coder
+  {
+
+  public:
+
+    scint_id_coder(int scint_id);
+    scint_id_coder(int sector_id, int tower_id, int fiber_id);
+    virtual
+    ~scint_id_coder()
+    {
+    }
+
+    virtual void
+      identify(std::ostream& os = std::cout) const
+    {
+      os <<"scint_id_coder with "
+          <<"scint_ID("<<scint_ID<<") = "
+          <<"sector_ID("<<sector_ID<<"), "
+          <<"tower_ID("<<tower_ID<<"), "
+          <<"fiber_ID("<<fiber_ID<<")"<<std::endl;
+    }
+
+    int scint_ID;
+    int sector_ID;
+    int tower_ID;
+    int fiber_ID;
+
+    static const int kfiber_bit = 12;
+    static const int ktower_bit = 12;
+    static const int ksector_bit = 8;
+
+    ClassDef(PHG4CylinderGeom_Spacalv3::scint_id_coder,1)
+  };
+
+  //! convert tower_ID + sector ID to eta and z bins as in other cylindrical calorimeters
+  //! @return: a std::pair of zbin and phibin number
+  virtual std::pair<int,int> get_tower_z_phi_ID(const int tower_ID, const int sector_ID) const;
+
 protected:
   double sidewall_thickness;
   double sidewall_outer_torr;
   std::string sidewall_mat;
-
+  int max_phi_bin_in_sec;
 
   tower_map_t sector_tower_map;
 
-ClassDef(PHG4CylinderGeom_Spacalv3,2)
+ClassDef(PHG4CylinderGeom_Spacalv3,3)
 
 };
 
diff --git a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3LinkDef.h b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3LinkDef.h
index 03ffd41e30..5995591f4c 100644
--- a/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3LinkDef.h
+++ b/simulation/g4simulation/g4detectors/PHG4CylinderGeom_Spacalv3LinkDef.h
@@ -2,5 +2,6 @@
 
 #pragma link C++ class PHG4CylinderGeom_Spacalv3+;
 #pragma link C++ class PHG4CylinderGeom_Spacalv3::geom_tower+;
+#pragma link C++ class PHG4CylinderGeom_Spacalv3::scint_id_coder+;
 
 #endif /* __CINT__ */
diff --git a/simulation/g4simulation/g4detectors/PHG4FullProjSpacalCellReco.cc b/simulation/g4simulation/g4detectors/PHG4FullProjSpacalCellReco.cc
new file mode 100644
index 0000000000..10c6da0079
--- /dev/null
+++ b/simulation/g4simulation/g4detectors/PHG4FullProjSpacalCellReco.cc
@@ -0,0 +1,444 @@
+#include "PHG4FullProjSpacalCellReco.h"
+
+#include "PHG4CylinderGeomContainer.h"
+#include "PHG4CylinderGeom_Spacalv3.h"
+#include "PHG4CylinderCell_Spacalv1.h"
+#include "PHG4CylinderCellGeomContainer.h"
+#include "PHG4CylinderCellGeom.h"
+#include "PHG4CylinderCellGeom_Spacalv1.h"
+#include "PHG4CylinderCellContainer.h"
+#include "PHG4CylinderCellDefs.h"
+
+#include <g4main/PHG4Hit.h>
+#include <g4main/PHG4HitContainer.h>
+#include <fun4all/Fun4AllReturnCodes.h>
+#include <fun4all/Fun4AllServer.h>
+#include <phool/PHNodeIterator.h>
+#include <phool/PHCompositeNode.h>
+#include <phool/PHIODataNode.h>
+#include <fun4all/getClass.h>
+
+#include <cmath>
+#include <cstdlib>
+#include <iostream>
+#include <sstream>
+#include <cassert>
+#include <boost/foreach.hpp>
+#include <exception>
+
+using namespace std;
+
+PHG4FullProjSpacalCellReco::PHG4FullProjSpacalCellReco(const string &name) :
+    SubsysReco(name), _timer(PHTimeServer::get()->insert_new(name.c_str())), chkenergyconservation(
+        0)
+{
+}
+
+int
+PHG4FullProjSpacalCellReco::InitRun(PHCompositeNode *topNode)
+{
+  PHNodeIterator iter(topNode);
+
+  // Looking for the DST node
+  PHCompositeNode *dstNode;
+  dstNode = dynamic_cast<PHCompositeNode*>(iter.findFirst("PHCompositeNode",
+      "DST"));
+  if (!dstNode)
+    {
+      std::cout << PHWHERE << "DST Node missing, doing nothing." << std::endl;
+      exit(1);
+    }
+  hitnodename = "G4HIT_" + detector;
+  PHG4HitContainer *g4hit = findNode::getClass<PHG4HitContainer>(topNode,
+      hitnodename.c_str());
+  if (!g4hit)
+    {
+      cout << "Could not locate g4 hit node " << hitnodename << endl;
+      exit(1);
+    }
+  cellnodename = "G4CELL_" + detector;
+  PHG4CylinderCellContainer *cells = findNode::getClass<
+      PHG4CylinderCellContainer>(topNode, cellnodename);
+  if (!cells)
+    {
+      cells = new PHG4CylinderCellContainer();
+      PHIODataNode<PHObject> *newNode = new PHIODataNode<PHObject>(cells,
+          cellnodename.c_str(), "PHObject");
+      dstNode->addNode(newNode);
+    }
+
+  geonodename = "CYLINDERGEOM_" + detector;
+  PHG4CylinderGeomContainer *geo =
+      findNode::getClass<PHG4CylinderGeomContainer>(topNode,
+          geonodename.c_str());
+  if (!geo)
+    {
+      cout << "Could not locate geometry node " << geonodename << endl;
+      exit(1);
+    }
+  if (verbosity > 0)
+    {
+      geo->identify();
+    }
+  seggeonodename = "CYLINDERCELLGEOM_" + detector;
+  PHG4CylinderCellGeomContainer *seggeo = findNode::getClass<
+      PHG4CylinderCellGeomContainer>(topNode, seggeonodename.c_str());
+  if (!seggeo)
+    {
+      seggeo = new PHG4CylinderCellGeomContainer();
+      PHCompositeNode *runNode = dynamic_cast<PHCompositeNode*>(iter.findFirst(
+          "PHCompositeNode", "RUN"));
+      PHIODataNode<PHObject> *newNode = new PHIODataNode<PHObject>(seggeo,
+          seggeonodename.c_str(), "PHObject");
+      runNode->addNode(newNode);
+    }
+
+  map<int, PHG4CylinderGeom *>::const_iterator miter;
+  pair<map<int, PHG4CylinderGeom *>::const_iterator,
+      map<int, PHG4CylinderGeom *>::const_iterator> begin_end =
+      geo->get_begin_end();
+  map<int, std::pair<double, int> >::iterator sizeiter;
+  for (miter = begin_end.first; miter != begin_end.second; ++miter)
+    {
+      const PHG4CylinderGeom *layergeom_raw = miter->second;
+      assert(layergeom_raw);
+
+      // a special implimentation of PHG4CylinderGeom is required here.
+      const PHG4CylinderGeom_Spacalv3 *layergeom =
+          dynamic_cast<const PHG4CylinderGeom_Spacalv3 *>(layergeom_raw);
+
+      if (!layergeom)
+        {
+          cout << "PHG4FullProjSpacalCellReco::InitRun - Fatal Error -"
+              << " require to work with a version of SPACAL geometry of PHG4CylinderGeom_Spacalv3 or higher. "
+              << "However the incoming geometry has version "
+              << layergeom_raw->ClassName() << endl;
+          exit(1);
+        }
+
+//      int layer = layergeom->get_layer();
+
+// create geo object and fill with variables common to all binning methods
+      PHG4CylinderCellGeom_Spacalv1 *layerseggeo =
+          new PHG4CylinderCellGeom_Spacalv1();
+      layerseggeo->set_layer(layergeom->get_layer());
+      layerseggeo->set_radius(layergeom->get_radius());
+      layerseggeo->set_thickness(layergeom->get_thickness());
+
+      if (verbosity > 1)
+        {
+          layergeom->identify();
+        }
+
+      layerseggeo->set_binning(PHG4CylinderCellDefs::spacalbinning);
+
+      // construct a map to convert tower_ID into the older eta bins.
+
+      const PHG4CylinderGeom_Spacalv3::tower_map_t & tower_map =
+          layergeom->get_sector_tower_map();
+      const PHG4CylinderGeom_Spacalv3::sector_map_t & sector_map =
+          layergeom->get_sector_map();
+      const int nphibin = layergeom->get_azimuthal_n_sec()
+          * layergeom->get_max_phi_bin_in_sec();
+      const double deltaphi = 2. * M_PI / nphibin;
+
+      typedef map<double, int> map_z_tower_z_ID_t;
+      map_z_tower_z_ID_t map_z_tower_z_ID;
+      double phi_min = NAN;
+
+      BOOST_FOREACH(const PHG4CylinderGeom_Spacalv3::tower_map_t::value_type& tower_pair, tower_map)
+        {
+
+          const int & tower_ID = tower_pair.first;
+          const PHG4CylinderGeom_Spacalv3::geom_tower & tower =
+              tower_pair.second;
+
+          // inspect index in sector 0
+          std::pair<int, int> tower_z_phi_ID = layergeom->get_tower_z_phi_ID(
+              tower_ID, 0);
+
+          const int & tower_ID_z = tower_z_phi_ID.first;
+          const int & tower_ID_phi = tower_z_phi_ID.second;
+
+          if (tower_ID_phi == 0)
+            {
+              //assign phi min according phi bin 0
+              phi_min = atan2(tower.centralY, tower.centralX + 0.25*(tower.pDx1+tower.pDx2+tower.pDx3+tower.pDx4))
+                  + sector_map.begin()->second;
+            }
+
+          if (tower_ID_phi == layergeom->get_max_phi_bin_in_sec() / 2)
+            {
+              //assign eta min according phi bin 0
+              map_z_tower_z_ID[tower.centralZ] = tower_ID_z;
+            }
+          // ...
+        } //       BOOST_FOREACH(const PHG4CylinderGeom_Spacalv3::tower_map_t::value_type& tower_pair, tower_map)
+
+      assert(not isnan(phi_min));
+      layerseggeo->set_phimin(phi_min);
+      layerseggeo->set_phistep(deltaphi);
+      layerseggeo->set_phibins(nphibin);
+
+      PHG4CylinderCellGeom_Spacalv1::tower_z_ID_eta_bin_map_t tower_z_ID_eta_bin_map;
+      int eta_bin = 0;
+      BOOST_FOREACH( map_z_tower_z_ID_t::value_type& z_tower_z_ID, map_z_tower_z_ID)
+        {
+          tower_z_ID_eta_bin_map[z_tower_z_ID.second] = eta_bin;
+          eta_bin++;
+        }
+      layerseggeo->set_tower_z_ID_eta_bin_map(tower_z_ID_eta_bin_map);
+      layerseggeo->set_etabins(eta_bin);
+      layerseggeo->set_etamin(NAN);
+      layerseggeo->set_etastep(NAN);
+
+      //build eta bin maps
+      BOOST_FOREACH(const PHG4CylinderGeom_Spacalv3::tower_map_t::value_type& tower_pair, tower_map)
+        {
+
+          const int & tower_ID = tower_pair.first;
+          const PHG4CylinderGeom_Spacalv3::geom_tower & tower =
+              tower_pair.second;
+
+          // inspect index in sector 0
+          std::pair<int, int> tower_z_phi_ID = layergeom->get_tower_z_phi_ID(
+              tower_ID, 0);
+          const int & tower_ID_z = tower_z_phi_ID.first;
+          const int & tower_ID_phi = tower_z_phi_ID.second;
+          const int & etabin = tower_z_ID_eta_bin_map[tower_ID_z];
+
+          if (tower_ID_phi == layergeom->get_max_phi_bin_in_sec() / 2)
+            {
+              // half z-range
+              const double dz = fabs(
+                  0.5 * (tower.pDy1 + tower.pDy2) * sin(tower.pRotationAngleX));
+
+              const double eta_central = -log(
+                  tan(0.5 * atan2(tower.centralY, tower.centralZ)));
+              // half eta-range
+              const double deta =
+                  fabs(
+                      eta_central
+                          - (-log(
+                              tan(
+                                  0.5
+                                      * atan2(tower.centralY,
+                                          tower.centralZ + dz)))));
+
+              layerseggeo->set_etabounds(etabin,
+                  make_pair<double, double>(eta_central - deta,
+                      eta_central + deta));
+              layerseggeo->set_zbounds(etabin,
+                  make_pair<double, double>(tower.centralZ - dz,
+                      tower.centralZ + dz));
+            }
+          // ...
+        } //       BOOST_FOREACH(const PHG4CylinderGeom_Spacalv3::tower_map_t::value_type& tower_pair, tower_map)
+
+      // add geo object filled by different binning methods
+      seggeo->AddLayerCellGeom(layerseggeo);
+      if (verbosity > 1)
+        {
+          cout <<"PHG4FullProjSpacalCellReco::InitRun::"<<Name()<<" - Done layer"<< (layergeom->get_layer()) <<". Print out the cell geometry:"<<endl;
+          layerseggeo->identify();
+        }
+    }
+  return Fun4AllReturnCodes::EVENT_OK;
+}
+
+int
+PHG4FullProjSpacalCellReco::process_event(PHCompositeNode *topNode)
+{
+  _timer.get()->restart();
+  PHG4HitContainer *g4hit = findNode::getClass<PHG4HitContainer>(topNode,
+      hitnodename.c_str());
+  if (!g4hit)
+    {
+      cout
+          << "PHG4FullProjSpacalCellReco::process_event - Fatal Error - Could not locate g4 hit node "
+          << hitnodename << endl;
+      exit(1);
+    }
+  PHG4CylinderCellContainer *cells = findNode::getClass<
+      PHG4CylinderCellContainer>(topNode, cellnodename);
+  if (!cells)
+    {
+      cout
+          << "PHG4FullProjSpacalCellReco::process_event - Fatal Error - could not locate cell node "
+          << cellnodename << endl;
+      exit(1);
+    }
+
+  PHG4CylinderCellGeomContainer *seggeo = findNode::getClass<
+      PHG4CylinderCellGeomContainer>(topNode, seggeonodename.c_str());
+  if (!seggeo)
+    {
+      cout
+          << "PHG4FullProjSpacalCellReco::process_event - Fatal Error - could not locate cell geometry node "
+          << seggeonodename << endl;
+      exit(1);
+    }
+
+  PHG4CylinderGeomContainer *layergeo = findNode::getClass<
+      PHG4CylinderGeomContainer>(topNode, geonodename.c_str());
+  if (!layergeo)
+    {
+      cout
+          << "PHG4FullProjSpacalCellReco::process_event - Fatal Error - Could not locate sim geometry node "
+          << geonodename << endl;
+      exit(1);
+    }
+
+  PHG4HitContainer::LayerIter layer;
+  pair<PHG4HitContainer::LayerIter, PHG4HitContainer::LayerIter> layer_begin_end =
+      g4hit->getLayers();
+  for (layer = layer_begin_end.first; layer != layer_begin_end.second; ++layer)
+    {
+      // layer loop
+      PHG4HitContainer::ConstIterator hiter;
+      PHG4HitContainer::ConstRange hit_begin_end = g4hit->getHits(*layer);
+
+      PHG4CylinderCellGeom *geo_raw = seggeo->GetLayerCellGeom(*layer);
+      PHG4CylinderCellGeom_Spacalv1 * geo =
+          dynamic_cast<PHG4CylinderCellGeom_Spacalv1 *>(geo_raw);
+      assert(geo);
+
+      const PHG4CylinderGeom *layergeom_raw = layergeo->GetLayerGeom(*layer);
+      assert(layergeom_raw);
+      // a special implimentation of PHG4CylinderGeom is required here.
+      const PHG4CylinderGeom_Spacalv3 *layergeom =
+          dynamic_cast<const PHG4CylinderGeom_Spacalv3 *>(layergeom_raw);
+      assert(layergeom);
+
+      for (hiter = hit_begin_end.first; hiter != hit_begin_end.second; ++hiter)
+        {
+          // hit loop
+          int scint_id = hiter->second->get_scint_id();
+
+          unsigned int key = static_cast<unsigned int>(scint_id);
+          if (celllist.find(key) == celllist.end())
+            {
+              // decode scint_id
+              PHG4CylinderGeom_Spacalv3::scint_id_coder decoder(scint_id);
+
+              // convert to z_ID, phi_ID
+              std::pair<int, int> tower_z_phi_ID =
+                  layergeom->get_tower_z_phi_ID(decoder.tower_ID,
+                      decoder.sector_ID);
+              const int & tower_ID_z = tower_z_phi_ID.first;
+              const int & tower_ID_phi = tower_z_phi_ID.second;
+
+              // convert tower_ID_z to to eta bin number
+              int etabin = -1;
+              try
+                {
+                  etabin = geo->get_etabin(tower_ID_z);
+                }
+              catch (exception & e)
+                {
+                  cout << "Print cell geometry:" << endl;
+                  geo->identify();
+                  cout << "Print scint_id_coder:" << endl;
+                  decoder.identify();
+                  cout << "Print the hit:" << endl;
+                  hiter->second->print();
+                  cout << "PHG4FullProjSpacalCellReco::process_event::" << Name()
+                      << " - Fatal Error - " << e.what() << endl;
+                  exit(1);
+                }
+
+              celllist[key] = new PHG4CylinderCell_Spacalv1();
+              celllist[key]->set_layer(*layer);
+              celllist[key]->set_phibin(tower_ID_phi);
+              celllist[key]->set_etabin(etabin);
+              celllist[key]->set_fiber_ID(decoder.fiber_ID);
+            }
+
+          celllist[key]->add_edep(hiter->first, hiter->second->get_edep(),
+              hiter->second->get_light_yield());
+
+        } // end loop over g4hits
+      int numcells = 0;
+      for (map<unsigned int, PHG4CylinderCell *>::const_iterator mapiter =
+          celllist.begin(); mapiter != celllist.end(); ++mapiter)
+        {
+          cells->AddCylinderCell(*layer, mapiter->second);
+          numcells++;
+          if (verbosity > 1)
+            {
+              cout << "PHG4FullProjSpacalCellReco::process_event::" << Name()
+                  << " - " << "Adding cell in bin eta "
+                  << (mapiter->second->get_bineta())
+                  <<" phi "<< (mapiter->second->get_binphi())
+                  <<" fiber "<< (mapiter->second->get_fiber_ID())
+                  << ", energy dep: "
+                                    << mapiter->second->get_edep()<< ", light yield: "
+                                    << mapiter->second->get_light_yield() << endl;
+            }
+        }
+      celllist.clear();
+      if (verbosity > 0)
+        {
+          cout << "PHG4FullProjSpacalCellReco::process_event::" << Name() << " - "
+              << " found " << numcells << " fibers with energy deposition"
+              << endl;
+        }
+    }
+
+  if (chkenergyconservation or verbosity > 4)
+    {
+      CheckEnergy(topNode);
+    }
+  _timer.get()->stop();
+  return Fun4AllReturnCodes::EVENT_OK;
+}
+
+int
+PHG4FullProjSpacalCellReco::End(PHCompositeNode *topNode)
+{
+  return Fun4AllReturnCodes::EVENT_OK;
+}
+
+int
+PHG4FullProjSpacalCellReco::CheckEnergy(PHCompositeNode *topNode)
+{
+  PHG4HitContainer *g4hit = findNode::getClass<PHG4HitContainer>(topNode,
+      hitnodename.c_str());
+  PHG4CylinderCellContainer *cells = findNode::getClass<
+      PHG4CylinderCellContainer>(topNode, cellnodename);
+  double sum_energy_g4hit = 0.;
+  double sum_energy_cells = 0.;
+  PHG4HitContainer::ConstRange hit_begin_end = g4hit->getHits();
+  PHG4HitContainer::ConstIterator hiter;
+  for (hiter = hit_begin_end.first; hiter != hit_begin_end.second; ++hiter)
+    {
+      sum_energy_g4hit += hiter->second->get_edep();
+    }
+  PHG4CylinderCellContainer::ConstRange cell_begin_end =
+      cells->getCylinderCells();
+  PHG4CylinderCellContainer::ConstIterator citer;
+  for (citer = cell_begin_end.first; citer != cell_begin_end.second; ++citer)
+    {
+      sum_energy_cells += citer->second->get_edep();
+
+    }
+  // the fractional eloss for particles traversing eta bins leads to minute rounding errors
+  if (fabs(sum_energy_cells - sum_energy_g4hit) / sum_energy_g4hit > 1e-6)
+    {
+      cout << "PHG4FullProjSpacalCellReco::CheckEnergy - energy mismatch between cells: " << sum_energy_cells
+          << " and hits: " << sum_energy_g4hit
+          << " diff sum(cells) - sum(hits): "
+          << sum_energy_cells - sum_energy_g4hit << endl;
+      return -1;
+    }
+  else
+    {
+      if (verbosity > 0)
+        {
+          cout <<"PHG4FullProjSpacalCellReco::CheckEnergy::"<< Name() << " - total energy for this event: "
+              << sum_energy_g4hit << " GeV. Passed CheckEnergy" << endl;
+        }
+    }
+  return 0;
+}
+
diff --git a/simulation/g4simulation/g4detectors/PHG4FullProjSpacalCellReco.h b/simulation/g4simulation/g4detectors/PHG4FullProjSpacalCellReco.h
new file mode 100644
index 0000000000..c18d379dd8
--- /dev/null
+++ b/simulation/g4simulation/g4detectors/PHG4FullProjSpacalCellReco.h
@@ -0,0 +1,50 @@
+#ifndef PHG4FullProjSpacalCellReco_H
+#define PHG4FullProjSpacalCellReco_H
+
+#include <fun4all/SubsysReco.h>
+#include <phool/PHTimeServer.h>
+#include <string>
+#include <map>
+#include <vector>
+
+class PHCompositeNode;
+class PHG4CylinderCell;
+
+class PHG4FullProjSpacalCellReco : public SubsysReco
+{
+ public:
+
+  PHG4FullProjSpacalCellReco(const std::string &name = "HCALCELLRECO");
+
+  virtual ~PHG4FullProjSpacalCellReco(){}
+  
+  //! module initialization
+  int InitRun(PHCompositeNode *topNode);
+  
+    //! event processing
+  int process_event(PHCompositeNode *topNode);
+  
+  //! end of process
+  int End(PHCompositeNode *topNode);
+  
+  void Detector(const std::string &d) {detector = d;}
+
+  void checkenergy(const int i=1) {chkenergyconservation = i;}
+
+ protected:
+
+  int CheckEnergy(PHCompositeNode *topNode);
+
+
+  std::string detector;
+  std::string hitnodename;
+  std::string cellnodename;
+  std::string geonodename;
+  std::string seggeonodename;
+
+  PHTimeServer::timer _timer;
+  int chkenergyconservation;
+  std::map<unsigned int, PHG4CylinderCell *> celllist;
+};
+
+#endif
diff --git a/simulation/g4simulation/g4detectors/PHG4FullProjSpacalCellRecoLinkDef.h b/simulation/g4simulation/g4detectors/PHG4FullProjSpacalCellRecoLinkDef.h
new file mode 100644
index 0000000000..44478507b9
--- /dev/null
+++ b/simulation/g4simulation/g4detectors/PHG4FullProjSpacalCellRecoLinkDef.h
@@ -0,0 +1,5 @@
+#ifdef __CINT__
+
+#pragma link C++ class PHG4FullProjSpacalCellReco-!;
+
+#endif /* __CINT__ */
diff --git a/simulation/g4simulation/g4detectors/PHG4FullProjSpacalDetector.cc b/simulation/g4simulation/g4detectors/PHG4FullProjSpacalDetector.cc
index 7f23684ca0..ddbc6a3478 100644
--- a/simulation/g4simulation/g4detectors/PHG4FullProjSpacalDetector.cc
+++ b/simulation/g4simulation/g4detectors/PHG4FullProjSpacalDetector.cc
@@ -134,8 +134,14 @@ PHG4FullProjSpacalDetector::Construct_AzimuthalSeg()
       _geom->is_azimuthal_seg_visible() and (not _geom->is_virualize_fiber()));
   wall_VisAtt->SetForceSolid(true);
 
+  if (_geom->get_sidewall_thickness()>0)
     {
       // end walls
+      if (_geom->get_construction_verbose() >= 1)
+        {
+          cout << "PHG4FullProjSpacalDetector::Construct_AzimuthalSeg::" << GetName()
+              << " - construct end walls." << endl;
+        }
       G4Tubs* wall_solid = new G4Tubs(G4String(GetName() + string("_EndWall")),
           _geom->get_radius() * cm + _geom->get_sidewall_outer_torr() * cm,
           _geom->get_max_radius() * cm - _geom->get_sidewall_outer_torr() * cm,
@@ -178,8 +184,14 @@ PHG4FullProjSpacalDetector::Construct_AzimuthalSeg()
         }
     }
 
+  if (_geom->get_sidewall_thickness()>0)
     {
       // side walls
+      if (_geom->get_construction_verbose() >= 1)
+        {
+          cout << "PHG4FullProjSpacalDetector::Construct_AzimuthalSeg::" << GetName()
+              << " - construct side walls." << endl;
+        }
       G4Box* wall_solid = new G4Box(G4String(GetName() + string("_SideWall")),
           _geom->get_sidewall_thickness() * cm / 2.0,
           _geom->get_thickness() * cm / 2.
@@ -283,6 +295,7 @@ PHG4FullProjSpacalDetector::Construct_Fibers_SameLengthFiberPerTower(
 
   G4Vector3D v_zshift = G4Vector3D(tan(g_tower.pTheta) * cos(g_tower.pPhi),
       tan(g_tower.pTheta) * sin(g_tower.pPhi), 1) * g_tower.pDz;
+  int fiber_ID = 0;
   for (int ix = 0; ix < g_tower.NFiberX; ix++)
 //  int ix = 0;
     {
@@ -305,7 +318,6 @@ PHG4FullProjSpacalDetector::Construct_Fibers_SameLengthFiberPerTower(
           if ((ix + iy) % 2 == 1)
             continue; // make a triangle pattern
 
-          const int fiber_ID = ix * 1000 + iy;
 
           const double weighted_iy = static_cast<double>(iy)
               / (g_tower.NFiberY - 1.);
@@ -344,6 +356,8 @@ PHG4FullProjSpacalDetector::Construct_Fibers_SameLengthFiberPerTower(
           const G4double fiber_length = vector_fiber.mag();
 
           min_fiber_length = min(fiber_length, min_fiber_length);
+
+          ++fiber_ID;
         }
     }
 
@@ -431,9 +445,9 @@ PHG4FullProjSpacalDetector::Construct_Fibers(
 {
   assert(_geom);
 
-  int fiber_count = 0;
   G4Vector3D v_zshift = G4Vector3D(tan(g_tower.pTheta) * cos(g_tower.pPhi),
       tan(g_tower.pTheta) * sin(g_tower.pPhi), 1) * g_tower.pDz;
+  int fiber_ID = 0;
   for (int ix = 0; ix < g_tower.NFiberX; ix++)
     {
       const double weighted_ix = static_cast<double>(ix)
@@ -454,7 +468,6 @@ PHG4FullProjSpacalDetector::Construct_Fibers(
           if ((ix + iy) % 2 == 1)
             continue; // make a triangle pattern
 
-          const int fiber_ID = ix * 1000 + iy;
 
           const double weighted_iy = static_cast<double>(iy)
               / (g_tower.NFiberY - 1.);
@@ -530,16 +543,16 @@ PHG4FullProjSpacalDetector::Construct_Fibers(
               fiber_ID, overlapcheck_fiber);
           fiber_vol[fiber_physi] = fiber_ID;
 
-          fiber_count++;
+          ++fiber_ID;
         }
     }
 
   if (_geom->get_construction_verbose() >= 3)
     cout << "PHG4FullProjSpacalDetector::Construct_Fibers::" << GetName()
-        << " - constructed tower ID " << g_tower.id << " with " << fiber_count
+        << " - constructed tower ID " << g_tower.id << " with " << fiber_ID
         << " fibers" << endl;
 
-  return fiber_count;
+  return fiber_ID;
 }
 
 //! a block along z axis built with G4Trd that is slightly tapered in x dimension
diff --git a/simulation/g4simulation/g4detectors/PHG4HcalCellReco.cc b/simulation/g4simulation/g4detectors/PHG4HcalCellReco.cc
index 00b924a57c..2bff08854d 100644
--- a/simulation/g4simulation/g4detectors/PHG4HcalCellReco.cc
+++ b/simulation/g4simulation/g4detectors/PHG4HcalCellReco.cc
@@ -196,7 +196,7 @@ PHG4HcalCellReco::process_event(PHCompositeNode *topNode)
 	      celllist[key]->set_etabin(slatno);
 	    }
 
-	  celllist[key]->add_edep(hiter->first, hiter->second->get_edep());
+	  celllist[key]->add_edep(hiter->first, hiter->second->get_edep(),hiter->second->get_light_yield());
 	} // end loop over g4hits
       int numcells = 0;
       for (map<unsigned int, PHG4CylinderCell *>::const_iterator mapiter = celllist.begin();mapiter != celllist.end() ; ++mapiter)
diff --git a/simulation/g4simulation/g4detectors/PHG4SpacalDetector.cc b/simulation/g4simulation/g4detectors/PHG4SpacalDetector.cc
index 878f7235f0..66eec802a6 100644
--- a/simulation/g4simulation/g4detectors/PHG4SpacalDetector.cc
+++ b/simulation/g4simulation/g4detectors/PHG4SpacalDetector.cc
@@ -37,6 +37,7 @@
 #include <cassert>
 #include <cmath>
 #include <sstream>
+#include <boost/foreach.hpp>
 
 using namespace std;
 
@@ -88,12 +89,12 @@ PHG4SpacalDetector::IsInCylinderActive(const G4VPhysicalVolume * volume)
       if (fiber_vol.find(volume) != fiber_vol.end())
         return FIBER_CLADING;
 
-      if (calo_vol.find(volume) != calo_vol.end())
-        return ABSORBER;
-
       if (block_vol.find(volume) != block_vol.end())
         return ABSORBER;
 
+      if (calo_vol.find(volume) != calo_vol.end())
+        return SUPPORT;
+
     }
   return INACTIVE;
 }
@@ -143,38 +144,14 @@ PHG4SpacalDetector::Construct(G4LogicalVolume* logicWorld)
       logicWorld, false, 0, overlapcheck);
 
 
+  // install sectors
   std::pair<G4LogicalVolume *,G4Transform3D> psec = Construct_AzimuthalSeg();
   G4LogicalVolume *sec_logic = psec.first;
   const G4Transform3D & sec_trans = psec.second;
-
-//  int n_sec_construct =
-//      (_geom->is_virualize_fiber()) ? 1 : (_geom->is_azimuthal_seg_visible()?2:_geom->get_azimuthal_n_sec());
-  int n_sec_construct =
-      (_geom->is_virualize_fiber()) ? 1 : (_geom->get_azimuthal_n_sec());
-
-  if (_geom->is_virualize_fiber() or _geom->is_azimuthal_seg_visible())
-    {
-
-      cout
-          << "PHG4SpacalDetector::Construct - WARNING - "
-          <<"only construct "<<n_sec_construct<<" sectors for visualization purpose!!!"
-          << endl;
-
-    }
-  else
-    {
-      cout
-          << "PHG4SpacalDetector::Construct - INFO - "
-          <<"start construction of "<<n_sec_construct<<" sectors."
-          << endl;
-
-    }
-
-  for (int sec = 0; sec < n_sec_construct; ++sec)
+  BOOST_FOREACH(const SpacalGeom_t::sector_map_t::value_type& val, _geom->get_sector_map())
     {
-
-      const double rot = twopi / (double) (_geom->get_azimuthal_n_sec())
-          * ((double) (sec) - n_sec_construct/2);
+      const int sec = val.first;
+      const double rot = val.second;
 
       G4Transform3D sec_place = G4RotateZ3D(rot) * sec_trans;
 
@@ -187,7 +164,7 @@ PHG4SpacalDetector::Construct(G4LogicalVolume* logicWorld)
       calo_vol[calo_phys] = sec;
 
     }
-  _geom->set_nscint(_geom->get_nscint() * n_sec_construct);
+  _geom->set_nscint(_geom->get_nscint() *  _geom->get_sector_map().size());
 
   if (active)
     {
diff --git a/simulation/g4simulation/g4detectors/PHG4SpacalDetector.h b/simulation/g4simulation/g4detectors/PHG4SpacalDetector.h
index bdc60b7904..03ea8f20ca 100644
--- a/simulation/g4simulation/g4detectors/PHG4SpacalDetector.h
+++ b/simulation/g4simulation/g4detectors/PHG4SpacalDetector.h
@@ -133,6 +133,7 @@ class PHG4SpacalDetector : public PHG4Detector
     FIBER_CORE = 1,
     FIBER_CLADING = 0,
     ABSORBER = -1,
+    SUPPORT = -2,
     INACTIVE = -100
   };
 
diff --git a/simulation/g4simulation/g4detectors/PHG4SpacalSteppingAction.cc b/simulation/g4simulation/g4detectors/PHG4SpacalSteppingAction.cc
index 208111c60e..49198eb90c 100644
--- a/simulation/g4simulation/g4detectors/PHG4SpacalSteppingAction.cc
+++ b/simulation/g4simulation/g4detectors/PHG4SpacalSteppingAction.cc
@@ -8,6 +8,7 @@
 
 #include "PHG4SpacalSteppingAction.h"
 #include "PHG4SpacalDetector.h"
+#include "PHG4CylinderGeom_Spacalv3.h"
 
 #include <g4main/PHG4HitContainer.h>
 #include <g4main/PHG4Hitv1.h>
@@ -67,12 +68,55 @@ PHG4SpacalSteppingAction::UserSteppingAction(const G4Step* aStep, bool)
       G4StepPoint * prePoint = aStep->GetPreStepPoint();
       G4StepPoint * postPoint = aStep->GetPostStepPoint();
       int scint_id = -1;
-      if (isactive == PHG4SpacalDetector::FIBER_CORE)
-        scint_id = prePoint->GetTouchable()->GetReplicaNumber(2);
-      else if (isactive == PHG4SpacalDetector::FIBER_CLADING)
-        scint_id = prePoint->GetTouchable()->GetReplicaNumber(1);
-      else if (isactive == PHG4SpacalDetector::ABSORBER)
-        scint_id = prePoint->GetTouchable()->GetReplicaNumber(0);
+
+      if (//
+          detector_->get_geom()->get_config() == PHG4SpacalDetector::SpacalGeom_t::kFullProjective_2DTaper //
+          or //
+      detector_->get_geom()->get_config() == PHG4SpacalDetector::SpacalGeom_t::kFullProjective_2DTaper_SameLengthFiberPerTower//
+      )
+        {
+          //SPACAL ID that is associated with towers
+          int sector_ID =0;
+          int tower_ID = 0;
+          int fiber_ID = 0;
+
+          if (isactive == PHG4SpacalDetector::FIBER_CORE)
+            {
+
+              fiber_ID = prePoint->GetTouchable()->GetReplicaNumber(1);
+              tower_ID = prePoint->GetTouchable()->GetReplicaNumber(2);
+              sector_ID  = prePoint->GetTouchable()->GetReplicaNumber(3);
+
+            }
+
+          else if (isactive == PHG4SpacalDetector::FIBER_CLADING)
+            {
+              fiber_ID = prePoint->GetTouchable()->GetReplicaNumber(0);
+              tower_ID = prePoint->GetTouchable()->GetReplicaNumber(1);
+              sector_ID  = prePoint->GetTouchable()->GetReplicaNumber(2);
+            }
+
+          else if (isactive == PHG4SpacalDetector::ABSORBER)
+            {
+              tower_ID = prePoint->GetTouchable()->GetReplicaNumber(0);
+              sector_ID  = prePoint->GetTouchable()->GetReplicaNumber(1);
+            }
+
+          // compact the tower/sector/fiber ID into 32 bit scint_id, so we could save some space for SPACAL hits
+          scint_id = PHG4CylinderGeom_Spacalv3::scint_id_coder(sector_ID, tower_ID, fiber_ID).scint_ID;
+
+        }
+      else
+        {
+          // other configuraitons
+          if (isactive == PHG4SpacalDetector::FIBER_CORE)
+            scint_id = prePoint->GetTouchable()->GetReplicaNumber(2);
+          else if (isactive == PHG4SpacalDetector::FIBER_CLADING)
+            scint_id = prePoint->GetTouchable()->GetReplicaNumber(1);
+          else
+            scint_id = prePoint->GetTouchable()->GetReplicaNumber(0);
+        }
+
       //       cout << "track id " << aTrack->GetTrackID() << endl;
       //        cout << "time prepoint: " << prePoint->GetGlobalTime() << endl;
       //        cout << "time postpoint: " << postPoint->GetGlobalTime() << endl;
diff --git a/simulation/g4simulation/g4main/PHG4Hitv1.cc b/simulation/g4simulation/g4main/PHG4Hitv1.cc
index 5f892259ea..e57f9353cb 100644
--- a/simulation/g4simulation/g4main/PHG4Hitv1.cc
+++ b/simulation/g4simulation/g4main/PHG4Hitv1.cc
@@ -40,7 +40,7 @@ PHG4Hitv1::print() const {
     {
       PROPERTY prop_id = static_cast<PROPERTY>(i->first);
       pair<const string, PROPERTY_TYPE> property_info = get_property_info(prop_id);
-      cout << "\t" << i->first << ":\t" << property_info.first << " = \t";
+      cout << "\t" << prop_id << ":\t" << property_info.first << " = \t";
       switch(property_info.second)
 	{
 	case type_int: