Commit in lcsim/src/org/lcsim/mc/CCDSim on MAIN
CCD.java+6-51.4 -> 1.5
CCDClusterCenter.java+4-41.3 -> 1.4
CCDSim.java+237-591.5 -> 1.6
RawTrackSegment.java+68-361.1 -> 1.2
+315-104
4 modified files
Nick Sinev. Commiting changes in CCDSim package

lcsim/src/org/lcsim/mc/CCDSim
CCD.java 1.4 -> 1.5
diff -u -r1.4 -r1.5
--- CCD.java	17 Feb 2007 02:00:18 -0000	1.4
+++ CCD.java	3 Mar 2007 02:45:58 -0000	1.5
@@ -8,7 +8,7 @@
  * including CCD parameters and all signal parameters, including digitized signals
  * in each CCD pixel
  * @author sinev U of Oregon; [log in to unmask] ; SLAC x2970 
- * @version $Id: CCD.java,v 1.4 2007/02/17 02:00:18 sinev Exp $
+ * @version $Id: CCD.java,v 1.5 2007/03/03 02:45:58 sinev Exp $
  */
 public class CCD
 {
@@ -97,7 +97,7 @@
     }
 /**
  * constructor of any type CCD using specified CCD dimensions and specified CCDSpec and CCDElectronicsSpec
- * @param ar1 <code>double</code> inner radius of the cylindrical CCD, or  hole radius for endcap plane CCD (in mm)
+ * @param ar1 <code>double</code> middle of the sens.lr radius of the cylindrical CCD, or  hole radius for endcap plane CCD (in mm)
  * @param ar2 <code>double</code> minimum value of Z for cylindrical CCD  or external radius for endcap plane CCD (in mm)
  * @param ar3 <code>double</code> maximum value of Z for cylindrical CCD or z position of endcap plane CCD (in mm)
  * @param ec <code>boolean</code>, if true endcap CCD is built, if false - cylindrical CCD is built
@@ -124,7 +124,7 @@
       Stype = 0;                 
       minPhi = 0.;
       maxPhi = 2.0 * Math.PI;
-      minRad = ar1;
+      minRad = ar1 - 0.5*pixToL*EpiDepth;
       maxRad = minRad + pixToL * EpiDepth;
       depDir = 1; 
       minZ = ar2;
@@ -134,11 +134,11 @@
      {
       isEndcap=true;
       Stype = 1;
-      depDir = 0;
+      depDir = 1;
       direction[0] = 0.;
       direction[1] = 0.;
       direction[2]=1.;
-      if(ar3 < 0.) direction[2]=-1.;                 
+      if(ar3 < 0.) { direction[2]=-1.; depDir = -1; }                 
       minPhi = 0.;
       maxPhi = 2.0 * Math.PI;
       minRad = ar1;
@@ -481,6 +481,7 @@
      * and active CCD media (epi layer) extends in the direction of increasing radius or Z 
      */
     public int getDepDir() { return depDir; }
+    public double[] getOrientation() { return direction; }  
     /**
      * set layer number for this CCD
      * @param ln layer number

lcsim/src/org/lcsim/mc/CCDSim
CCDClusterCenter.java 1.3 -> 1.4
diff -u -r1.3 -r1.4
--- CCDClusterCenter.java	17 Feb 2007 02:00:18 -0000	1.3
+++ CCDClusterCenter.java	3 Mar 2007 02:45:58 -0000	1.4
@@ -6,7 +6,7 @@
  * class implementing algorithm of cluster center finding in CCD active pixels cluster.
  * Center finding algorithm is choosen from 3 currently implemented by setting Method parameter.
  * @author sinev ; U of Oregon; SLAC x2970; [log in to unmask]
- * @version $Id: CCDClusterCenter.java,v 1.3 2007/02/17 02:00:18 sinev Exp $
+ * @version $Id: CCDClusterCenter.java,v 1.4 2007/03/03 02:45:58 sinev Exp $
  * @see CCDCluster
  */
 public class CCDClusterCenter
@@ -61,7 +61,7 @@
      ncc = ccd.getNColumns();
      double wrp = (double) ncc;
      double maxcd = ncc*maxcfr;
-     dep = 0.0005 * (ccd.getDepDepth() + ccd.getEpiDepth());
+     dep = 0.0005 * ccd.getEpiDepth();
      // note: factor 10E-3 is due to micron - mm convertion
      totSignal=0.; 
      acr = 0.;
@@ -191,8 +191,8 @@
 		{
 		 r = ccd.getSurfRad();
 		 arclen = cl.getCenterCol() * ccd.getPixelSizeX() * 0.001;
-//		 r = r + ccd.getDepDir() * cl.getCenterDepth();
 		 phi = ccd.getMinPhi() + arclen/r;
+		 r = r + ccd.getDepDir() * cl.getCenterDepth();
 		 pos[0] = r * Math.cos(phi);
 		 pos[1] = r * Math.sin(phi);
 		 pos[2] = ccd.getMinZ() + cl.getCenterRow() * ccd.getPixelSizeY() * 0.001;
@@ -201,7 +201,7 @@
 		{
 		 pos[0] = -ccd.getMaxRad() + cl.getCenterCol() * ccd.getPixelSizeX() * 0.001;
 		 pos[1] = -ccd.getMaxRad() + cl.getCenterRow() * ccd.getPixelSizeY() * 0.001;
-		 pos[2] = ccd.getMinZ();
+		 pos[2] = ccd.getMinZ() + ccd.getDepDir() * cl.getCenterDepth();
                  type = 1;
 		}
 		if(lclh != null)

lcsim/src/org/lcsim/mc/CCDSim
CCDSim.java 1.5 -> 1.6
diff -u -r1.5 -r1.6
--- CCDSim.java	17 Feb 2007 02:00:18 -0000	1.5
+++ CCDSim.java	3 Mar 2007 02:45:58 -0000	1.6
@@ -1,11 +1,12 @@
 package org.lcsim.mc.CCDSim;
 import org.lcsim.util.aida.AIDA;
-import java.util.Random;
+import java.util.*;
 import org.lcsim.math.distribution.*;
 import org.lcsim.event.*;
 import hep.physics.vec.Hep3Vector;
 import hep.physics.vec.VecOp;
 import hep.physics.vec.BasicHep3Vector;
+import hep.aida.ITree;
 
 /**
  * A class to simulate the response of a ccd detector .
@@ -17,12 +18,14 @@
  *                 November 2006 - Nick Sinev - added diffusion in depleted region,
  *                                 Lorentz Angle, output TrackerHit, not SimTrackerHit
  * @author  sinev  U of Oregon, SLAC. x2970 <BR>
- * @version $Id: CCDSim.java,v 1.5 2007/02/17 02:00:18 sinev Exp $
+ * @version $Id: CCDSim.java,v 1.6 2007/03/03 02:45:58 sinev Exp $
 *    
 */
 public class CCDSim
 {
    private boolean debug = false;
+   private boolean doHist = true;
+   private AIDA aida = AIDA.defaultInstance();
    private CCDSpec ccdspec = null;
    private LorentzAngle _LA = new LorentzAngle();
 /*.....................................
@@ -122,6 +125,11 @@
     private Random _ran = new Random();
     private CCDElectronicsSpec cespec;
     private RawTrackSegment parent = null;
+    double Ef = 0.;                               // electric field in v/cm in depleted region
+    double cspd = Math.sqrt(_D);                  // carriers speed (in cm/s) (if there is no el.field)
+
+        float[][] tpxls = new float[NPXMAXX][NPXMAXY];
+        float[][] pxls = new float[NPXMAXX][NPXMAXY];
     
     // start with everything static and hardcoded, will progress to more complicated models later (if necessary)
     /** 
@@ -522,7 +530,7 @@
  * Simulates high energy particle hit using RawTrackSegment as input
  *
  * @param ccd  <code>CCD</code> object in which simulated digitized hit will be created
- * @param th   <code>SimTrackerHit</code> - initial (not digitised) MC simulated hit
+ * @param th   <code>RawTrackSegment</code> - collection of SimTrackerHit of the same MC part. in same layer
  * @return 1 if success, -1 if failed (hit outside CCD or wrong type hit)
  */    
     public int simTrack(CCD ccd, RawTrackSegment ts)
@@ -561,9 +569,148 @@
             {
              System.out.println("Simulating hit in lr "+ccd.getLayerNumber()+" "+x+" "+y+" "+z+" dir: "+cx+" "+cy+" "+cz);
             } 
-           } 
-	   return simTrack(ccd,x,y,z,cx,cy,cz);
+           }
+           if((doHist) && (dir != null) && (!ccd.isEndcap()))
+           {
+            MCParticle mcp = ts.getMCParticle();
+            Hep3Vector mcpm = mcp.getMomentum();
+            double mcPt = Math.sqrt(mcpm.x()*mcpm.x()+mcpm.y()*mcpm.y()); 
+            double tlm = mcpm.z()/mcPt;
+            double tls = cz/Math.sqrt(cx*cx+cy*cy);
+            double deltl = tls-tlm;
+            if((tlm>0.) && (Math.abs(deltl)<0.2)) aida.cloud1D("CCDsim/ diff ts tl and MC tl, pos tl").fill(deltl); 
+            if((tlm<0.) && (Math.abs(deltl)<0.2)) aida.cloud1D("CCDsim/ diff ts tl and MC tl, neg tl").fill(deltl); 
+           }
+	   double cdr = 0.;
+           if(!ccd.isEndcap()) cdr = (x*cx + y*cy)/(Math.sqrt(x*x+y*y) * Math.sqrt(cx*cx+cy*cy+cz*cz));
+           int status =  0;
+           if(!ccd.isEndcap() && Math.abs(cdr) < 0.01) status = simBlob(ccd,ts); 
+           if((ts.isBlob() || (ts.getHits().size() < 3))&&(status!=1)) status = simBlob(ccd,ts);
+           if(!ts.isBlob() && (status != 1)) 
+             status = simTrack(ccd,x,y,z,cx,cy,cz);
+           if(status != 1)
+           {
+            if(!ts.isBlob()) System.out.println("Failed to simulate track segment consisting of hits with r,z: ");
+            if(ts.isBlob()) System.out.println("Failed to simulate blob of hits consisting of hits with r,z: ");
+            List<SimTrackerHit> hits = ts.getHits();
+            for(SimTrackerHit hit:hits)
+            {
+             System.out.println(Math.sqrt(hit.getPoint()[0]*hit.getPoint()[0]+hit.getPoint()[1]*hit.getPoint()[1])+"  "+hit.getPoint()[2]);
+            } 
+           }  
+	   return status;
 	}
+
+     public int simBlob(CCD ccd, RawTrackSegment ts)
+     {
+        int  i,j,k,ofsx,ofsy,coln,rown;
+        float ttheta,fracdep,eloss,efar,ddep,depth,fracst;
+        float exsv,exss,dl,xlc,ylc;
+        double D,dis,ang0,ang,minz,xc,yc,zc;
+        double pxsizx,pxsizy,minx,miny;
+
+//        System.out.println("Processing blob!"); 
+        ddep = ccd.getDepDir();        
+        pxsizx = 0.001 * ccd.getPixelSizeX();              // convert um into mm
+        pxsizy = 0.001 * ccd.getPixelSizeY();
+	minz = ccd.getMinZ();
+	minx = - ccd.getMaxRad();
+	miny = minx;
+        Ef = _bV/(DEPDEP*0.0001);                    // electric field in v/cm in depleted region
+        cspd = _mu*Ef;                             // carriers speed (in cm/s)
+        _LorA = _LA.getLorentzAngle(_Bfield,Ef,_T,_typ);
+        double tlza = Math.tan(_LorA);
+        double rref = 0.;
+        for(i=0; i<NPXMAXX; i++)
+            for(j=0; j<NPXMAXY; j++)
+                pxls[i][j] = 0.f;
+        List<SimTrackerHit> hits = ts.getHits();
+        double _laxd=0.;
+        rown = 0;
+        coln = 0;
+        ofsx = 0;
+        ofsy = 0;
+        boolean firstht = true; 
+        for(SimTrackerHit hit:hits)
+        {
+         xc = hit.getPoint()[0];
+         yc = hit.getPoint()[1];
+         zc = hit.getPoint()[2];
+         double hr = Math.sqrt(xc*xc+yc*yc);
+         if(!ccd.isEndcap())
+         {
+          if(ccd.getStype() != 0) return -1;
+          rref = ccd.getSurfRad();
+          depth = (float) (ddep * (hr - rref) * 1000.);
+          if((depth > 0.) && (depth < EPIDEP))
+          {
+           ang0 = Math.atan2(yc,xc);
+           if(firstht)
+           {  
+            coln = (int) Math.floor(ang0*rref/pxsizx);
+            rown = (int) Math.floor((zc-minz)/pxsizy);
+            ofsx = coln - _icenx;
+            ofsy = rown - _iceny;
+            firstht = false;
+           }
+           if(depth<=DEPDEP)  _laxd = depth * tlza;
+           if(depth > DEPDEP) _laxd = DEPDEP*tlza;
+           xlc = (float) (((ang0*rref - coln * pxsizx) * 1000.) + _laxd);
+           ylc = (float) ((zc - minz - rown * pxsizy) * 1000.);
+           eloss = (float) ((hit.getdEdx()*1.E9)/LOSSPERPAIR);
+           addSegment(xlc,ylc,depth,eloss);
+          }  
+         }
+         if(ccd.isEndcap())
+         {
+          if(ccd.getStype() == 0) return -1;
+          depth = (float) (ddep * (zc - minz) * 1000.);
+//          System.out.println("CCD is endcap, depth: "+depth+" ddep= "+ddep+" zc= "+zc+" r "+hr);
+          if((depth > 0.) && (depth < EPIDEP))
+          {
+           if(firstht)
+           {  
+            coln = (int) Math.floor((xc-minx)/pxsizx);
+            rown = (int) Math.floor((yc-miny)/pxsizy);
+            ofsx = coln - _icenx;
+            ofsy = rown - _iceny;
+            firstht=false;
+           }
+	   xlc = (float) ((xc - minx - coln*pxsizx) * 1000.); // microns!
+	   ylc = (float) ((yc - miny - rown*pxsizy) * 1000.); // microns!
+//           System.out.println("CCD is endcap, minz is: "+minz+" col: "+coln+" row: "+rown+" ofsx "+ofsx+" ofsy "+ofsy);
+           eloss = (float) ((hit.getdEdx()*1.E9)/LOSSPERPAIR);
+           addSegment(xlc,ylc,depth,eloss);
+          }  
+         }
+        }   
+        int napad = 0;
+        for(i=0; i<NPXMAXX; i++)
+        {
+            int ig = i+ofsx;
+            for(j=0; j<NPXMAXY; j++)
+            {
+                int jg =j+ofsy;
+                exsv = pxls[i][j];
+                exss = (float) (Math.sqrt((double)exsv));
+                pxls[i][j]= exsv + exss * (float) _ran.nextGaussian();
+                if((ig >= ccd.getNColumns()) && !ccd.isEndcap()) ig-=ccd.getNColumns();
+                if((ig < 0) && !ccd.isEndcap()) ig+=ccd.getNColumns();
+                if((pxls[i][j] > 1.)&&(ig >= 0)&&(ig < ccd.getNColumns())
+                 &&(jg >= 0)&&(jg < ccd.getNRows()))
+                {   
+                  ccd.addSignal(jg,ig,pxls[i][j],parent);
+                  napad++;
+                }
+            }
+        }
+        if(napad > 0) { ccd.incPartCount(); return 1;}
+        if(napad == 0) 
+        {
+         System.out.println("No signal was generated by blob started in CCD lr "+ccd.getLayerNumber()+" col "+coln+" row "+rown);
+        } 
+      return -1;
+     }  
     
 /**
  * Simulates high energy particle hit, using coordinates and directions
@@ -579,51 +726,51 @@
     public int simTrack(CCD ccd, double x, double y, double z, double ic, double jc, double kc)
     {
 
-        float[][] tpxls = new float[NPXMAXX][NPXMAXY];
-        float[][] pxls = new float[NPXMAXX][NPXMAXY];
-        int delta_flg,i,j,k,ofsx,ofsy,coln,rown,pthr,clthr,adc;
-        float lfact,ttheta,fracdep,eloss,edelta=0.f,efar,xc,yc,ddep,depth,fracst;
-        float sigman,sigmaw,frnar,exsv,exss,dl,xlc,ylc;
-        double r,rma,rmi,rref,D,dis,dL,dL1,dL2,x0,y0,z0,ino,jno,kno,ang0,ang,minz;
-        double ita,jta,kta,iyl,jyl,kyl,il,jl,kl,len,pxsizx,pxsizy,adcsc,noise,minx,miny;
+        int  i,j,k,ofsx,ofsy,coln,rown;
+        float lfact,ttheta,fracdep,eloss,efar,xc,yc,ddep,depth,fracst;
+        float  exsv,exss,dl,xlc,ylc;
+        double D,dis,dL1,dL2,ino,jno,kno,ang0,ang,minz;
+        double ita,jta,kta,iyl,jyl,kyl,len,pxsizx,pxsizy,minx,miny;
            
-                
         ddep = ccd.getDepDir();        
-        adcsc = cespec.getADCscale();
-        noise = cespec.getNoiseRMS();
-        pthr = cespec.getSinglePixThr();
-        clthr = cespec.getClustThr();
         pxsizx = 0.001 * ccd.getPixelSizeX();              // convert um into mm
         pxsizy = 0.001 * ccd.getPixelSizeY();
 	minz = ccd.getMinZ();
 	minx = - ccd.getMaxRad();
 	miny = minx;
-        double Ef = _bV/(DEPDEP*0.0001);                    // electric field in v/cm in depleted region
-        double cspd = _mu*Ef;                             // carriers speed (in cm/s)
+        Ef = _bV/(DEPDEP*0.0001);                    // electric field in v/cm in depleted region
+        cspd = _mu*Ef;                             // carriers speed (in cm/s)
         _LorA = _LA.getLorentzAngle(_Bfield,Ef,_T,_typ);
-	il = ic;
-	jl = jc;
-	kl = kc;
+        double x0=0.;
+        double y0=0.;
+        double z0=0.;
+	double il = ic;
+	double jl = jc;
+	double kl = kc;
 	xlc = 0.f;
 	ylc = 0.f;
 	coln = 0;
 	rown = 0;
         double tlza = Math.tan(_LorA);
+        double r=0.;
+        double cdr = 0.; 
+        double rref = 0.;
+        double dL=0.; 
         if(ccd.isEndcap()) tlza = 0.;
 //        System.out.println("Tangent of Lorentz Angle is: "+tlza+" for electric field "+Ef+" V/cm and B field "+_Bfield+" Tesla");  
         if(!ccd.isEndcap())
         {
-         if(ccd.getStype() != 0) return -1;
+         if(ccd.getStype() != 0) { System.out.println("Non-cyl surface in barrel!"); return -1; }
          rref = ccd.getSurfRad();  
-         rma = ccd.getMaxRad();
-         rmi = ccd.getMinRad();
          r = Math.sqrt(x*x + y*y);
          if(Math.abs(r-rref) > 0.1)
          {
            System.out.println("CcdSim: hit R "+r+" while reference CCD rad: "+rref);  
              return -1; // hit is too far away from CCD
          }
-	 double cdr = (x*ic + y*jc+z*kc)/(Math.sqrt(x*x+y*y+z*z) * Math.sqrt(ic*ic+jc*jc+kc*kc));
+//	 double cdr = (x*ic + y*jc+z*kc)/(Math.sqrt(x*x+y*y+z*z) * Math.sqrt(ic*ic+jc*jc+kc*kc));
+	 cdr = (x*ic + y*jc)/(Math.sqrt(x*x+y*y) * Math.sqrt(ic*ic+jc*jc+kc*kc));
+         if(Math.abs(cdr) < 0.01) { System.out.println("Attempt to simulate track segment in barrel with r-proj < 0.01"); return -1;}
 	 // this is cosine of angle between R and hit direction
 	 dL = (rref - r)/cdr; // length of the vector along hit direction to reach rref
          x0 = x + dL * ic;   //  Coordinates of extrapolation of 
@@ -639,7 +786,7 @@
          // directed in the direction of increasing phi 
          ang0 = Math.atan2(y0,x0);
          ang = ang0 + 0.5 * Math.PI;
-		 if(ang0 < 0.) ang0 = ang0 + 2. * Math.PI;
+	 if(ang0 < 0.) ang0 = ang0 + 2. * Math.PI;
          ita = Math.cos(ang);
          jta = Math.sin(ang);
          kta = 0.;
@@ -695,6 +842,7 @@
 		  jl = -jl;
 		 } 
 	}
+
         ofsx = coln - _icenx;
         ofsy = rown - _iceny;
         for(i=0; i<NPXMAXX; i++)
@@ -712,9 +860,67 @@
          depth = (float) (k * dl * kl);
          if(depth<=DEPDEP)
          {  
-//            _laxd = depth * tlza;  
+            _laxd = depth * tlza;  
             xc = (float) (xlc + k * dl * il + _laxd);
             yc = (float) (ylc + k * dl * jl);
+            addSegment(xc,yc,depth,eloss);
+         }    
+         if(depth > DEPDEP)
+         {
+               _laxd = DEPDEP*tlza;
+               xc = (float) (xlc + k * dl * il + _laxd);
+               yc = (float) (ylc + k * dl * jl);
+            addSegment(xc,yc,depth,eloss);
+         }
+        }
+        int napad = 0;
+        for(i=0; i<NPXMAXX; i++)
+        {
+            int ig = i+ofsx;
+            for(j=0; j<NPXMAXY; j++)
+            {
+                int jg =j+ofsy;
+                exsv = pxls[i][j] * lfact;
+                exss = (float) (Math.sqrt((double)exsv));
+                pxls[i][j]= exsv + exss * (float) _ran.nextGaussian();
+                if((ig >= ccd.getNColumns()) && !ccd.isEndcap()) ig-=ccd.getNColumns();
+                if((ig < 0) && !ccd.isEndcap()) ig+=ccd.getNColumns();
+                if((pxls[i][j] > 1.)&&(ig >= 0)&&(ig < ccd.getNColumns())
+                 &&(jg >= 0)&&(jg < ccd.getNRows()))
+                {   
+                  ccd.addSignal(jg,ig,pxls[i][j],parent);
+                  napad++;
+                }
+            }
+        }
+        if(napad > 0) { ccd.incPartCount(); return 1;}
+        if(napad == 0)
+        {  
+         if(!ccd.isEndcap())
+         { 
+          System.out.println("Failed hit on barr. lr "+ccd.getLayerNumber()+" col. offs. "+
+                ofsx+" row offs. "+ofsy);
+          System.out.println(" x0,y0,z0 "+x0+" "+y0+" "+z0+" original x,y,z "+x+" "+y+" "+z);
+          System.out.println(" local dir: "+il+" "+jl+" "+kl+" original dir: "+ic+" "+jc+" "+kc);
+          System.out.println(" or. hit r: "+r+" rref "+rref+" cdr: "+cdr+" dL "+dL);
+         }
+         if(ccd.isEndcap()) System.out.println("Failed hit on ec lr "+ccd.getLayerNumber()+" col. offs. "+
+                ofsx+" row offs. "+ofsy+", x0,y0,z0 "+x0+" "+y0+" "+z0+" local dir: "+il+" "+jl+" "+kl);
+          return -1;
+        }
+        System.out.println("How do we get here !?");
+        return 1;
+    }
+    
+
+    private void addSegment(float x, float y, float depth, float eloss)
+    {
+         float xc = x;
+         float yc = y;
+         int i=0;
+         int j=0; 
+         if(depth<=DEPDEP)
+         {  
             double trt = (depth*0.0001)/cspd;     // charge travel time (in seconds). Depth is converted into cm here
             double difsgm = 10000.*1.375 * Math.sqrt(_D*trt);    // diffusion sigma (10000. is to convert cm into microns) 
             // _icenx and _iceny are the indexies of the central pixel in
@@ -759,15 +965,13 @@
          }
          if(depth > DEPDEP)
          {
-               xc = (float) (xlc + k * dl * il + _laxd);
-               yc = (float) (ylc + k * dl * jl);
                //
                // we don't need to watch for track part leaving pixel area,
                // as getpix() takes care about it
                // 
-               sigman = SIGCO * (depth - DEPDEP);
-               sigmaw = SIGCO * (2*EPIDEP - depth - DEPDEP);
-               frnar = (float) (1.f-(depth-DEPDEP)/(EPIDEP-DEPDEP));
+               float sigman = SIGCO * (depth - DEPDEP);
+               float sigmaw = SIGCO * (2*EPIDEP - depth - DEPDEP);
+               float frnar = (float) (1.f-(depth-DEPDEP)/(EPIDEP-DEPDEP));
                getpix(tpxls,xc,yc,sigman);
                for(i=0; i<NPXMAXX; i++)
                     for(j=0; j<NPXMAXY; j++)
@@ -777,35 +981,9 @@
                     for(j=0; j<NPXMAXY; j++)
                         pxls[i][j]+= (1.-frnar) * eloss * tpxls[i][j];
          }
-        }
-        int napad = 0;
-        for(i=0; i<NPXMAXX; i++)
-        {
-            int ig = i+ofsx;
-            for(j=0; j<NPXMAXY; j++)
-            {
-                int jg =j+ofsy;
-                exsv = pxls[i][j] * lfact;
-                exss = (float) (Math.sqrt((double)exsv));
-                pxls[i][j]= exsv + exss * (float) _ran.nextGaussian();
-                if((ig >= ccd.getNColumns()) && !ccd.isEndcap()) ig-=ccd.getNColumns();
-                if((ig < 0) && !ccd.isEndcap()) ig+=ccd.getNColumns();
-                if((pxls[i][j] > 1.)&&(ig >= 0)&&(ig < ccd.getNColumns())
-                 &&(jg >= 0)&&(jg < ccd.getNRows()))
-                {   
-                  ccd.addSignal(jg,ig,pxls[i][j],parent);
-                  napad++;
-                }
-            }
-        }
-        if(napad > 0) { ccd.incPartCount(); return 1;}
-        if(napad == 0) System.out.println("Failed hit lr "+ccd.getLayerNumber()+" col. offs. "+
-                ofsx+" row offs. "+ofsy);
-        return -1;
     }
     
     
-    
     private float geteloss()
     {
         double rn = ranlan.nextLandau();

lcsim/src/org/lcsim/mc/CCDSim
RawTrackSegment.java 1.1 -> 1.2
diff -u -r1.1 -r1.2
--- RawTrackSegment.java	17 Feb 2007 02:00:18 -0000	1.1
+++ RawTrackSegment.java	3 Mar 2007 02:45:58 -0000	1.2
@@ -20,7 +20,7 @@
 
 /**
  * @author sinev U of Oregon; SLAC x2970; [log in to unmask]
- * @version $Id: RawTrackSegment.java,v 1.1 2007/02/17 02:00:18 sinev Exp $
+ * @version $Id: RawTrackSegment.java,v 1.2 2007/03/03 02:45:58 sinev Exp $
  */
 
            public class RawTrackSegment
@@ -31,10 +31,18 @@
 	    double dedx = 0.;
 	    double time = 0.;
 	    double atime = 0.;
+            double minr=999999.;
+            double maxr=0.;
+            double minz=999999.;
+            double maxz=-999999.;
+            double minphi = 999999.;
+            double maxphi = -999999.;   
             double[] cg = new double[3];
             double[] aco = new double[3];
+            double tlp = 0.;
             Hep3Vector dir = null;
             double base = 0.;
+            boolean _isblob = false;
     
 	    public RawTrackSegment(int lrn)
 	    {
@@ -57,56 +65,73 @@
 	     dedx+=hit.getdEdx();
 	     atime+=hit.getTime();
 	     nhits++;
+             double[] pnt = hit.getPoint();
+             double rht = Math.sqrt(pnt[0]*pnt[0]+pnt[1]*pnt[1]);
+             double hphi = Math.atan2(pnt[1],pnt[0]);
+             if(rht < minr) minr=rht;
+             if(rht > maxr) maxr=rht;
+             if(pnt[2]<minz) minz=pnt[2];
+             if(pnt[2]>maxz) maxz=pnt[2];
+             if(hphi < minphi) minphi = hphi;
+             if(hphi > maxphi) maxphi = hphi;  
              if(nhits > 1)
              {
               SimTrackerHit frsthit = hits.get(0);
               double[] fhc = frsthit.getPoint();
               double rfh = Math.sqrt(fhc[0]*fhc[0]+fhc[1]*fhc[1]);
-              double[] pnt = hit.getPoint();
-              double rht = Math.sqrt(pnt[0]*pnt[0]+pnt[1]*pnt[1]);
               double dx = 0.;
               double dy = 0.;
               double dz = 0.;
               double bdz = Math.abs(pnt[2]-fhc[2]);
-              if(bdz < Math.abs(rht-rfh))
-              {             
-               if(rht > rfh)
-               {
-                dx = pnt[0]-fhc[0];
-                dy = pnt[1]-fhc[1];
-                dz = pnt[2]-fhc[2];
-               }
-               else
-               {
-                dx = fhc[0]-pnt[0];
-                dy = fhc[1]-pnt[1];
-                dz = fhc[2]-pnt[2];
-               }
+              double tlc = bdz/(rht-rfh);
+              if((nhits > 2) && (Math.abs(tlc-tlp)/Math.abs(tlp) > 0.3))
+              { 
+                _isblob = true;
+//               System.out.println("For hit: "+nhits+" tlc= "+tlc+" while tlp= "+tlp+" Hit makes segment blob"); 
               }
-              else
+              if(nhits == 2) tlp=tlc;
+              if(!_isblob)
               {
-               if(Math.abs(pnt[2]) > Math.abs(fhc[2]))
-               {
+               if(bdz < Math.abs(rht-rfh))
+               {             
+                if(rht > rfh)
+                {
                 dx = pnt[0]-fhc[0];
                 dy = pnt[1]-fhc[1];
-                dz = pnt[2]-fhc[2];                 
+                dz = pnt[2]-fhc[2];
+                }
+                else
+                {
+                 dx = fhc[0]-pnt[0];
+                 dy = fhc[1]-pnt[1];
+                 dz = fhc[2]-pnt[2];
+                }
                }
                else
                {
-                dx = fhc[0]-pnt[0];
-                dy = fhc[1]-pnt[1];
-                dz = fhc[2]-pnt[2];
-               } 
-              } 
-              double nbase = Math.sqrt(dx*dx+dy*dy+dz*dz);
-              if(nbase > base)
-              {  
-               Hep3Vector hvec = new BasicHep3Vector(dx,dy,dz);
-               Hep3Vector ndir = VecOp.unit(hvec);
-               dir = new BasicHep3Vector(ndir.x(),ndir.y(),ndir.z());
-               base = nbase;
-              }  
-             }  
+                if(Math.abs(pnt[2]) > Math.abs(fhc[2]))
+                {
+                 dx = pnt[0]-fhc[0];
+                 dy = pnt[1]-fhc[1];
+                 dz = pnt[2]-fhc[2];                 
+                }
+                else
+                {
+                 dx = fhc[0]-pnt[0];
+                 dy = fhc[1]-pnt[1];
+                 dz = fhc[2]-pnt[2];
+                }
+               }  
+               double nbase = Math.sqrt(dx*dx+dy*dy+dz*dz);
+               if(nbase > base)
+               {  
+                Hep3Vector hvec = new BasicHep3Vector(dx,dy,dz);
+                Hep3Vector ndir = VecOp.unit(hvec);
+                dir = new BasicHep3Vector(ndir.x(),ndir.y(),ndir.z());
+                base = nbase;
+               }  
+              }
+             }   
              for(int j=0; j<3; j++)
              {
 	      aco[j]+=hit.getPoint()[j];
@@ -114,10 +139,17 @@
              }
 	     time = atime/nhits;
 	    }
-
+            
+            public boolean isBlob() { return _isblob; }  
 	    public double getdEdx() {return dedx;}
 	    public double getTime() {return time;}
 	    public double[] getCenter() { return cg;}
+            public double getMinR() { return minr; }
+            public double getMaxR() { return maxr; }
+            public double getMinZ() { return minz; }
+            public double getMaxZ() { return maxz; }
+            public double getMinPhi() { return minphi; }
+            public double getMaxPhi() { if((maxphi-minphi) > Math.PI) maxphi-=2.*Math.PI; return maxphi; } 
 	    public List<SimTrackerHit> getHits() { return hits; }
 	    public int getNHits() { return nhits; }
             public int getLayer() { return layer; }
CVSspam 0.2.8