diff -u -r1.1 -r1.2
--- SimpleResiduals.java	4 May 2012 15:45:44 -0000	1.1
+++ SimpleResiduals.java	8 May 2012 17:20:29 -0000	1.2
@@ -289,7 +289,7 @@

         Hep3Vector ecal_cl;
         Hep3Vector origin_pos;
         List<Hep3Vector> ecal_cls;

-        FastTrack fastTrack;

+        FastTrack fastTrack = new FastTrack(debug);

         List<SiTrackerHitStrip1D> stripList;
         boolean isaxial;
         String name;

@@ -313,8 +313,7 @@

             origin_pos = getFastTrackOrigin(event);
             
             //Create the fast track 

-            fastTrack = new FastTrack(origin_pos,ecal_cl,debug);
-

+            fastTrack.setTrack(origin_pos, ecal_cl);

             
             if ( debug ) System.out.println(fastTrack.toString());
            

diff -u -r1.1 -r1.2
--- FastTrack.java	4 May 2012 15:45:44 -0000	1.1
+++ FastTrack.java	8 May 2012 17:20:29 -0000	1.2
@@ -4,6 +4,7 @@

  */
 package org.lcsim.hps.users.phansson;
 

+import hep.physics.vec.BasicHep3Vector;

 import hep.physics.vec.Hep3Vector;
 import java.util.List;
 import org.lcsim.detector.identifier.IIdentifier;

@@ -22,22 +23,79 @@

     private Hep3Vector _endVec;
     private boolean _debug;
     private double _slope;

+    private double _const;
+    double _beamGapCorr;

-    public FastTrack(Hep3Vector startPos, Hep3Vector endPos, boolean debug) {

+    public FastTrack(boolean debug) {
+        _debug = debug; 
+        _beamGapCorr = 0.0;
+    }
+    
+    
+    
+    public FastTrack(Hep3Vector startPos, Hep3Vector endPos, double corr, boolean debug) {

         //create the track using a straight line in the non-bending plane
         //The non-bending plane is y-z:
         //    z is beam direction
         //    x is parallel to axial strips
         //    y is measurement direction/orthogonal to strip direction 
         _debug = debug;

+        _beamGapCorr = corr;    
+        setTrack(startPos,endPos);
+        
+    }
+    
+     public FastTrack(Hep3Vector startPos, double[] endPos, double corr, boolean debug) {
+        //create the track using a straight line in the non-bending plane
+        //The non-bending plane is y-z:
+        //    z is beam direction
+        //    x is parallel to axial strips
+        //    y is measurement direction/orthogonal to strip direction 
+        _debug = debug;
+        _beamGapCorr = corr;    
+        setTrack(startPos,endPos);
+        
+    }
+
+     
+    public void setTrack(Hep3Vector startPos, double[] endPos) {
+        //create the track using a straight line in the non-bending plane
+        //The non-bending plane is y-z:
+        //    z is beam direction
+        //    x is parallel to axial strips
+        //    y is measurement direction/orthogonal to strip direction 
+        Hep3Vector endPos_vec = new BasicHep3Vector(endPos[0],endPos[1],endPos[2]);
+        setTrack(startPos,endPos_vec);
+        
+    }
+     
+    
+    public double getCorrY(double y) {
+        double ynew = y;
+        if (y>0) ynew += _beamGapCorr;
+        else ynew -= _beamGapCorr;
+        return ynew;
+    }
+      
+    public void setTrack(Hep3Vector startPos, Hep3Vector endPos) {
+        //create the track using a straight line in the non-bending plane
+        //The non-bending plane is y-z:
+        //    z is beam direction
+        //    x is parallel to axial strips
+        //    y is measurement direction/orthogonal to strip direction 

         _startVec = startPos;

-        _endVec = endPos;

+        //Correct for a possible beamgap difference in the geometry of the ecal
+        double y = getCorrY(endPos.y());
+        Hep3Vector endPos_corr = new BasicHep3Vector(endPos.x(),y,endPos.z());
+        _endVec = endPos_corr;

         _slope = getSlope();

+        _const = getConst();

         
         
     }

-

+     
+     

     private double getSlope() {
         double z1 = _startVec.z();  
         double z2 = _endVec.z();  

@@ -45,45 +103,39 @@

         double y2 = _endVec.y();  
         return (y2-y1)/(z2-z1);
     }

+
+    private double getConst() {
+        return -1.0*(_slope*_startVec.z());  
+    }
+

     
     public double eval(double x) {
         // non-bend plane is y-z

-        return _slope*x;

+        return _slope*x+_const;

     }
     
     public double getFastTrackResidual(SiTrackerHitStrip1D stripCluster) {

-        //Get the measured position on the sensor which is in the non-bend plane        
-        double yHit;
-        double zPlane;
-        double yFT;
-        double res;
-        Hep3Vector posVec;
-
-        posVec = stripCluster.getPositionAsVector();
-        zPlane = posVec.z();
-        yHit = posVec.y();
-        yFT = this.eval(zPlane);
-        res = yHit - yFT;
-        
-         if(_debug) {
-            System.out.println("1D strip cluster pos as vec: " + posVec.toString());
-            System.out.println("zPlane " + zPlane + " xFT " + yFT + " yHit " + yHit + " --> res " + res);
-            
-        }
-        
-        
-       
-        
-        return res;

+        Hep3Vector posVec = stripCluster.getPositionAsVector();
+        double res = posVec.y() - this.eval(posVec.z());

+        if(_debug) {
+            System.out.println("Hit pos " + posVec.toString());
+            System.out.println("FastTrack position " + this.eval(posVec.z()) + " and z " + posVec.z() + " --> res " + res);
+        }
+        return res;        

     }
             

+    public void setEcalBeamgapCorr(double offset) {
+        this._beamGapCorr = offset;    
+    }

     
     public String toString() {
         String s = "FastTrack info: ";
         s += " startPos: " + _startVec.toString() + "\n";
         s += " endPos: " + _endVec.toString() + "\n";

+        s += " ECalOffsetY: " + _beamGapCorr + "\n";

         s += "slope: " + _slope;

+        s += "const: " + _const;

         return s;
     }
     

diff -u -r1.2 -r1.3
--- FastTrackResidualDriver.java	4 May 2012 18:23:35 -0000	1.2
+++ FastTrackResidualDriver.java	8 May 2012 17:20:29 -0000	1.3
@@ -25,7 +25,7 @@

 
 /**
  *

- * @author phansson

+ * @author phansson+

  */
 public class FastTrackResidualDriver extends Driver {
     

@@ -43,6 +43,8 @@

     double conversionZ;
     private static int crystalCols;
     private static int crystalRows;

+    double ecalBeamgapCorr;
+    int ecalClusterSel;

 
     private String outputPlotFileName = "";
     

@@ -52,9 +54,13 @@

     private List< List<IHistogram1D> > resy_org_layallhit = new ArrayList<List<IHistogram1D> >();
     private List< List<IHistogram1D> > resy_org_lay1hit = new ArrayList<List<IHistogram1D> >();
     private List< IHistogram1D > nhits_tracker = new ArrayList<IHistogram1D>();

+    private List< List<IHistogram1D> > nhits_tracker_layer = new ArrayList<List<IHistogram1D> >();

     private List< IHistogram1D > ncl_ecal = new ArrayList<IHistogram1D>();

+    private List< IHistogram1D > selcl_ecal_e = new ArrayList<IHistogram1D>(); 
+    private List< IHistogram1D > cl_ecal_e = new ArrayList<IHistogram1D>();

     private List< IHistogram2D > ncl_ecal_map = new ArrayList<IHistogram2D>(); 
     private List< IHistogram2D > nselcl_ecal_map = new ArrayList<IHistogram2D>(); 

+    private List< IHistogram2D> resy_2d_org_layallhit = new ArrayList<IHistogram2D>();

     
     
     public void startOfData() {

@@ -75,26 +81,49 @@

         crystalCols = 46;
         crystalRows = 5;
         

+        ecalBeamgapCorr = 20.0;
+        

         // Position of the conversion
         conversionZ = -1000.0; 
         

+        ecalClusterSel = 1;
+        

         IHistogramFactory hf = aida.histogramFactory();
         String side;
         for (int iSide=0;iSide<2;++iSide) {
             if(iSide==0) side="up";
             else side="down";

+            
+            ncl_ecal.add(hf.createHistogram1D("FT_ncl_ecal_" + side , 20, 0, 20));
+            selcl_ecal_e.add(hf.createHistogram1D("FT_selcl_ecal_e" + side , 100, 0, 20000));
+            cl_ecal_e.add(hf.createHistogram1D("FT_cl_ecal_e" + side , 100, 0, 20000));
+           
+
+            nhits_tracker.add(hf.createHistogram1D("FT_nhits_tracker_" + side , 15, 0, 15));
+
+            List<IHistogram1D> nhl_list  = new ArrayList<IHistogram1D>();
+            nhits_tracker_layer.add(nhl_list);
+            for (int i=0;i<5;++i) {
+                nhl_list.add(hf.createHistogram1D("FT_nhits_tracker_" + side + "_layer" + i , 15, 0, 15));
+            }
+
+                        
+            double res_axis[] = {-40,40};
+
+            
+            resy_2d_org_layallhit.add(hf.createHistogram2D("FT_resy_2d_org_layallhit_" + side, 5, 0.5, 5.5, 60, res_axis[0], res_axis[1]));
+            
+                   

             List<IHistogram1D> list  = new ArrayList<IHistogram1D>();
             resy_org.add(list);
             List<IHistogram1D> listLayAllHit  = new ArrayList<IHistogram1D>();
             resy_org_layallhit.add(listLayAllHit);
             List<IHistogram1D> listLay1Hit  = new ArrayList<IHistogram1D>();
             resy_org_lay1hit.add(listLay1Hit);            

-            
-            nhits_tracker.add(hf.createHistogram1D("FT_nhits_tracker_" + side , 40, 0, 40));
-            ncl_ecal.add(hf.createHistogram1D("FT_ncl_ecal_" + side , 60, 0, 60));
-            
-            
-             //Setup the ecal 2D plot

+
+                
+            //Setup the ecal 2D plot
+                        

             if (side == "up") {
                 ncl_ecal_map.add(hf.createHistogram2D("FT_ecal_hitmap_" + side, 92, -46.0, 46.0, 5, 0.0, 5.0));
                 nselcl_ecal_map.add(hf.createHistogram2D("FT_sel_ecal_hitmap_" + side, 92, -46.0, 46.0, 5, 0.0, 5.0));

@@ -105,13 +134,13 @@

 
             for (int iLayer=1;iLayer<6;++iLayer) {
                 

-                IHistogram1D h = hf.createHistogram1D("FT_res_" + side + "_l"+iLayer, 50, -25, 25);

+                IHistogram1D h = hf.createHistogram1D("FT_res_" + side + "_l"+iLayer, 50, res_axis[0], res_axis[1]);

                 list.add(h);
 

-                IHistogram1D hAll = hf.createHistogram1D("FT_res_LayAllHit_" + side + "_l"+iLayer, 50, -25, 25);

+                IHistogram1D hAll = hf.createHistogram1D("FT_res_LayAllHit_" + side + "_l"+iLayer, 50,res_axis[0], res_axis[1]);

                 listLayAllHit.add(hAll);
                 

-                IHistogram1D hLay1 = hf.createHistogram1D("FT_res_Lay1Hit_" + side + "_l"+iLayer, 50, -25, 25);

+                IHistogram1D hLay1 = hf.createHistogram1D("FT_res_Lay1Hit_" + side + "_l"+iLayer, 50,res_axis[0], res_axis[1]);

                 listLay1Hit.add(hLay1);
                 
             }

@@ -130,6 +159,11 @@

     public void setConversionZ(double z) {
         this.conversionZ = z;
     }

+    
+    public void setecalBeamgapCorr(double val) {
+        this.ecalBeamgapCorr = val;
+        System.out.println("offset set: "  + ecalBeamgapCorr);
+    }

 
     public void setDebug(boolean flag) {
         this.debug = flag;

@@ -143,6 +177,11 @@

         this.EcalZPosition = val;
     }
     

+    public void setEcalClusterSel(int id) {
+        this.ecalClusterSel = id;
+    }
+    
+    

     public void process(EventHeader event) {
         ++nevents;
         if( debug ) {

@@ -177,17 +216,20 @@

         int nhits;
         int nhitsInTracker;
         int nhitsInLayer1;

-        Hep3Vector ecal_cl;

         Hep3Vector origin_pos;

-        List<Hep3Vector> ecal_cls;
-        FastTrack fastTrack;
-        List<SiTrackerHitStrip1D> stripList;

+        List<Integer> ecal_cls;
+        FastTrack fastTrack = new FastTrack(debug);
+        fastTrack.setEcalBeamgapCorr(ecalBeamgapCorr); //Correct for geometry
+
+        //List<SiTrackerHitStrip1D> stripList;

         boolean isaxial;
         String name;
         SiSensor siSensor;
         int layer;
         String si_side;
         double res;

+        int selclids[];
+        int sel_ecal_idx;

         int layerIndex = -1;
         for (int iSide=0;iSide<2;++iSide) {
             

@@ -196,31 +238,44 @@

             
             ecal_cls = getEcalClustersForFastTracking(ecal_all_clusters, sides[iSide]);
             

-            if (debug) System.out.println("Found  " + ecal_cls.size() +" Ecal clusters on the " + sides[iSide]);

+            if (debug) {
+                    System.out.println("Found  " + ecal_cls.size() +" Ecal clusters on the " + sides[iSide] + ": " + ecal_cls.toString());
+            }
+            

             
             //if ( 1==1 ) return;
             
             ncl_ecal.get(iSide).fill(ecal_cls.size());
             
             if( ecal_cls.size() ==0 ) {

-                System.out.println("No clusters on this side...");

+                if(debug) System.out.println("No clusters on this side...");

                 continue;
             }
             
             
             //Fill map of Ecal hits

-            for( Hep3Vector cl : ecal_cls) {
-                int clpos[] = getCrystalPair(cl);

+            for( int icl=0; icl<ecal_cls.size(); ++icl) {
+                int clid = ecal_cls.get(icl); 
+                double p_cl[] = ecal_all_clusters.get(clid).getPosition();
+                int clpos[] = getCrystalPair(p_cl);

                 ncl_ecal_map.get(iSide).fill(clpos[0], clpos[1]);

+                cl_ecal_e.get(iSide).fill(ecal_all_clusters.get(icl).getEnergy());

             }
             

-            ecal_cl = selectCluster(ecal_cls);  
-
-            if (debug) System.out.println("Selected   " + ecal_cl.toString() +" cluster to use as pointer of fast track ");

+            sel_ecal_idx = selectCluster(ecal_cls,ecal_all_clusters); 
+            
+            if (debug) System.out.println("Selected clid " + sel_ecal_idx + " is cluster " + ecal_all_clusters.get(sel_ecal_idx).toString() + " and will be used as pointer of fast track ");
+            
+            if (sel_ecal_idx < 0) {
+                if (debug) System.out.println("No selected cluster!");
+                continue;
+            }

-             int selclpos[] = getCrystalPair(ecal_cl);
-             nselcl_ecal_map.get(iSide).fill(selclpos[0], selclpos[1]);

+            selclids = getCrystalPair(ecal_all_clusters.get(sel_ecal_idx).getPosition());
+            nselcl_ecal_map.get(iSide).fill(selclids[0], selclids[1]);
+            
+            selcl_ecal_e.get(iSide).fill(ecal_all_clusters.get(sel_ecal_idx).getEnergy());

              
             //Get "target" position i.e. the origin of the radiation
             origin_pos = getFastTrackOrigin(event);

@@ -228,9 +283,9 @@

             if (debug) System.out.println("Conversion started at " + origin_pos.toString());
             
             //Create the fast track 

-            fastTrack = new FastTrack(origin_pos,ecal_cl,debug);
-

+            fastTrack.setTrack(origin_pos,ecal_all_clusters.get(sel_ecal_idx).getPosition());

+

             if ( debug ) System.out.println(fastTrack.toString());
            
             

@@ -241,17 +296,25 @@

             if( debug ) System.out.println("There are " + nhitsInTracker + " hits on this side of the tracker");
             
             //most upstream layer nr is different for top and bottom

-            int firstLayer = 1;
-            if (sides[iSide] == "down") firstLayer=2;
-            nhitsInLayer1 = getNAxialHitsInLayers(trackerHits,sides[iSide],firstLayer);

+//            int firstLayer = 1;
+//            if (sides[iSide] == "down") firstLayer=2;
+//            nhitsInLayer1 = getNAxialHitsInLayers(trackerHits,sides[iSide],firstLayer);

-            if( debug ) System.out.println("There are " + nhitsInLayer1 + " hits on the first layer for this side (layer=" + firstLayer +")");

+            int nhitsInLayers[] = getNAxialHitsPerLayer(trackerHits,sides[iSide]);

-            nhits_tracker.get(iSide).fill(nhitsInTracker);
-            

+            if( debug ) {
+                System.out.println("Number of hits per layer:");
+                for (int i=0;i<5;++i){
+                    System.out.println("Layer " + (i+1) + " has " + nhitsInLayers[i] + " axial hits");                    
+                }
+            }

-            stripList = new ArrayList<SiTrackerHitStrip1D>();

+            nhits_tracker.get(iSide).fill(nhitsInTracker);
+            for (int i=0;i<5;++i){
+                nhits_tracker_layer.get(iSide).get(i).fill(nhitsInLayers[i]);
+            }
+            //stripList = new ArrayList<SiTrackerHitStrip1D>();

             for ( SiTrackerHitStrip1D stripCluster : trackerHits ) {
                 
                 isaxial = isAxialHit(stripCluster);

@@ -286,6 +349,8 @@

                 
                 res = fastTrack.getFastTrackResidual(stripCluster);
             

+                if (debug) System.out.println("-->res " + res);
+                

 //                1 2 3 4 5 6 7 8 9 10
 //                0   1   2   3   4
 //                  0   1   2   3   4

@@ -298,9 +363,10 @@

                 
                 if (nhitsInTracker==5) {
                     resy_org_layallhit.get(iSide).get(layerIndex).fill(res);

+                    resy_2d_org_layallhit.get(iSide).fill((double)layerIndex+1, res, 1);

                 }
                 

-                if (nhitsInLayer1>0) {

+                if (nhitsInLayers[0]>0) {

                     resy_org_lay1hit.get(iSide).get(layerIndex).fill(res);
                 }
             

@@ -336,6 +402,9 @@

         return pos;
     }
         

+    
+    
+    

 //    private List<HPSEcalCluster> getAllEcalClusters(EventHeader event) {
 //        
 //        List<HPSEcalCluster> clusters = event.get(HPSEcalCluster.class, "EcalReadoutHits"); 

@@ -350,72 +419,101 @@

  
 
     private int[] getCrystalPair(Hep3Vector cluster) {

-        double x = cluster.x();
-        double y = cluster.y();

+        double pos[] = {cluster.x(),cluster.y(),0.};
+        return getCrystalPair(pos);        
+    }
+    
+     private int[] getCrystalPair(double[] cluster) {
+        double x = cluster[0];
+        double y = cluster[1];

         int position[] = new int[2];
         position[0] = (int)Math.floor(x/crystalX);

-        //position[0] = crystalX * ix;
-        position[1] = (int) Math.floor( (y - beamGap * Math.signum(y) ) / crystalY);
-        //position[1] = crystalY * iy + beamGap * Math.signum(iy);
-        //System.out.println("x " + x + " y " + y + " -->  ix " + position[0] + " iy " + position[1]);
-        return position; 

-        //IDDecoder dec = dec.getSubdetector("Ecal").getIDDecoder();
-        //dec.setID(hit.getCellID());
-        //int ix = dec.getValue("ix");
-        //int iy = dec.getValue("iy");

+        position[1] = (int) Math.floor( (y - beamGap * Math.signum(y) ) / crystalY);
+         return position; 

         
     }
 
     private double[] getClusterPosition(HPSEcalCluster cluster) {
         CalorimeterHit hit = cluster.getSeedHit();
         //return cluster.getPosition();

-

+        

         //IDDecoder dec = dec.getSubdetector("Ecal").getIDDecoder();

-        dec.setID(hit.getCellID());
-        int ix = dec.getValue("ix");
-        int iy = dec.getValue("iy");
-        double position[] = new double[2];
-        position[0] = crystalX * ix;
-        position[1] = crystalY * iy + beamGap * Math.signum(iy);
-        return position; 

+//        dec.setID(hit.getCellID());
+//        int ix = dec.getValue("ix");
+//        int iy = dec.getValue("iy");
+//        double position[] = new double[2];
+//        position[0] = crystalX * ix;
+//        position[1] = crystalY * iy + beamGap * Math.signum(iy);
+        if (debug) System.out.println("Getting cluster position");
+        double pos[] = cluster.getPosition();
+        if (debug) System.out.println("Found " + pos[0] + "," + pos[1] + "," + pos[2]);
+        return pos; 

     }
     

-    private List<Hep3Vector> getEcalClustersForFastTracking(List<HPSEcalCluster> clusters, String side) {
-        if(!side.equalsIgnoreCase("up") && !side.equalsIgnoreCase("down")) {

+    private List<Integer> getEcalClustersForFastTracking(List<HPSEcalCluster> clusters, String side) {
+        if(side!="up" && side!="down") {

             throw new RuntimeException("This ecal side" + side + " do not exist!!");
         }

-        List<Hep3Vector> cls = new ArrayList<Hep3Vector>();

+        List<Integer> cls = new ArrayList<Integer>();

         boolean save;
         double [] pos;

-        for ( HPSEcalCluster cl : clusters) {

+        
+        for (int i=0;i<clusters.size();++i) {

             save=false;
             

-            pos = getClusterPosition(cl);

+            pos = getClusterPosition(clusters.get(i));

-            if(pos[1]>=0 && side.equalsIgnoreCase("up")) {

+            if(pos[1]>=0 && side=="up") {

                 save=true;

-            } else if(pos[1]<0 && side.equalsIgnoreCase("down")) {

+            } 
+            if(pos[1]<0 && side=="down") {

                 save=true;
             } 
             if(save) {

-                cls.add(new BasicHep3Vector(pos[0],pos[1],EcalZPosition));
-                //cls.add(new BasicHep3Vector(pos[0],pos[1],pos[2]));

+                //cls.add(new BasicHep3Vector(pos[0],pos[1],EcalZPosition));
+                cls.add(i);
+                //new BasicHep3Vector(pos[0],pos[1],pos[2]));

             }
             
         }
         return cls;
     }
      

-    private Hep3Vector selectCluster(List<Hep3Vector> clusters) {

+    private int selectCluster(List<Integer> ids, List<HPSEcalCluster> clusters) {

         //need to decide which cluster to take

-        if (clusters.size()>0) {
-            Hep3Vector pos = clusters.get(0);
-            return pos;
-        } else {
-            throw new RuntimeException("No cluser positions to choose from!");

+
+        if (clusters.size()==0) {
+            throw new RuntimeException("No clusters to select from!!!");
+        }
+
+        if (ids.size()==0) {
+            throw new RuntimeException("No idx to clusters to select from!!!");

         }
         

+        int sel_id = -1;
+        
+        switch (ecalClusterSel) {
+            case 1: 
+                //Require at least 3000MeV cluster
+                double E = -1.0;
+                for (int i=0;i<ids.size();++i) {
+                    if (clusters.get(ids.get(i)).getEnergy()>3000.0 && clusters.get(ids.get(i)).getEnergy()>E ) {
+                        sel_id = ids.get(i);
+                        E = clusters.get(ids.get(i)).getEnergy();
+                    }
+                }
+                break;
+            default:
+                //Use the first in list
+                if (ids.size()>0 && ids.get(0)<clusters.size()) {
+                    sel_id = ids.get(0);
+                }
+                break;
+        }
+                
+        return sel_id;        
+           

     }
 
     private List<HPSEcalCluster> getAllEcalClustersForFastTracking(EventHeader event) {

@@ -490,6 +588,48 @@

                           
         return nhits;
     }   

+     
+      private int[] getNAxialHitsPerLayer(List<SiTrackerHitStrip1D> trackerHits, String side) {
+        int nhits=0;
+        String si_side;
+        SiSensor siSensor;
+        String name;
+        int l;
+        int i;
+        int n[] = {0,0,0,0,0};
+        for ( SiTrackerHitStrip1D stripCluster : trackerHits ) {
+            
+            if(isAxialHit(stripCluster)==false) continue;
+           
+            si_side = getSideFromSiCluster(stripCluster);
+
+            if( side == si_side) {
+                
+                siSensor = stripCluster.getSensor();
+                name = siSensor.getName();
+                if ( name.length() < 14) {
+                    System.err.println("This name is too short!!");
+                    throw new RuntimeException("SiSensor name " + name + " is invalid?");
+                }
+            
+                l = getLayerFromSensorName(name);
+                
+                //Turn this into the layer scheme I use
+                if(side=="down") {
+                    //even numbers are used
+                    i = (l/2) - 1;
+                    n[i] = n[i] + 1;
+                } else {
+                    i = ((l+1)/2) - 1;
+                    n[i] = n[i] + 1;                    
+                }
+            
+            } 
+
+        }
+                          
+        return n;
+    }   

     
      private String getSideFromSiCluster(SiTrackerHitStrip1D stripCluster) {       
         Hep3Vector posVec = stripCluster.getPositionAsVector();

@@ -531,7 +671,7 @@

         IPlotter plotter_org = af.createPlotterFactory().create("HPS SVT Fast Track Residuals (All tracks)");
         plotter_org.createRegions(5,2,0);
         IPlotter plotter_org_LayAllHit = af.createPlotterFactory().create("HPS SVT Fast Track Residuals (All layers has hits)");

-        plotter_org_LayAllHit.createRegions(5,2,0);

+        plotter_org_LayAllHit.createRegions(6,2,0);

         IPlotter plotter_org_Lay1Hit = af.createPlotterFactory().create("HPS SVT Fast Track Residuals (1st layer has hits)");
         plotter_org_Lay1Hit.createRegions(5,2,0);
         

@@ -546,30 +686,58 @@

                 plotter_org.region(idx).plot(resy_org.get(iSide).get(iLayer));
                 plotter_org_LayAllHit.region(idx).plot(resy_org_layallhit.get(iSide).get(iLayer));
                 plotter_org_Lay1Hit.region(idx).plot(resy_org_lay1hit.get(iSide).get(iLayer));

-              

+
+

             }

+            //IPlotterStyle style = plotter_org_LayAllHit.style();
+            //style.regionBoxStyle().
+                    
+            plotter_org_LayAllHit.region(10+iSide).style().setParameter("hist2DStyle", "colorMap");
+            plotter_org_LayAllHit.region(10+iSide).style().setParameter("colorMapScheme", "rainbow");
+            plotter_org_LayAllHit.region(10+iSide).plot(resy_2d_org_layallhit.get(iSide));

         }
         plotter_org.show();
         plotter_org_LayAllHit.show();
         plotter_org_Lay1Hit.show();
         

-        //Hit mupltiplicity

+        
+        
+        //Hit multiplicity

         IPlotter plotter_hitmult = af.createPlotterFactory().create();
         plotter_hitmult.createRegions(1,2,0);
         plotter_hitmult.region(0).plot(nhits_tracker.get(0));
         plotter_hitmult.region(1).plot(nhits_tracker.get(1));
         plotter_hitmult.show();
         

+        IPlotter plotter_hitmult_layer = af.createPlotterFactory().create();
+        plotter_hitmult_layer.createRegions(5,2,0);
+        for ( int iSide=0;iSide<2;++iSide) {
+            int nLayers = resy_org.get(iSide).size();
+            for (int iLayer=0;iLayer<nLayers;++iLayer) {
+                int idx = 2*iLayer+iSide;
+                
+                //0 1 2 3 4
+                //1 3 5 7 9
+                plotter_hitmult_layer.region(idx).plot(nhits_tracker_layer.get(iSide).get(iLayer));
+              
+            }
+        }
+        plotter_hitmult_layer.show();
+        

         
         //ECal plots
         IPlotter plotter_ecalhitmult = af.createPlotterFactory().create();

-        plotter_ecalhitmult.createRegions(3,2,0);

+        plotter_ecalhitmult.createRegions(5,2,0);

         plotter_ecalhitmult.region(0).plot(ncl_ecal.get(0));
         plotter_ecalhitmult.region(1).plot(ncl_ecal.get(1));
         plotter_ecalhitmult.region(2).plot(ncl_ecal_map.get(0));
         plotter_ecalhitmult.region(3).plot(ncl_ecal_map.get(1));
         plotter_ecalhitmult.region(4).plot(nselcl_ecal_map.get(0));
         plotter_ecalhitmult.region(5).plot(nselcl_ecal_map.get(1));

+        plotter_ecalhitmult.region(6).plot(cl_ecal_e.get(0));
+        plotter_ecalhitmult.region(7).plot(cl_ecal_e.get(1));
+        plotter_ecalhitmult.region(8).plot(selcl_ecal_e.get(0));
+        plotter_ecalhitmult.region(9).plot(selcl_ecal_e.get(1));

         plotter_ecalhitmult.show();