lcsim/src/org/lcsim/contrib/uiowa
diff -u -r1.3 -r1.4
--- ShowerPointFinder.java 12 Jul 2008 01:20:34 -0000 1.3
+++ ShowerPointFinder.java 15 Jul 2008 23:57:15 -0000 1.4
@@ -1,4 +1,4 @@
-package org.lcsim.contrib.uiowa;
+package org.lcsim.contib.uiowa;
import java.util.*;
import org.lcsim.util.*;
@@ -35,16 +35,6 @@
protected int exam = 0; //to test neighbour
protected boolean debug = false;
- //sid01
- protected double ecalRs=1270; //ECAL inner R
- protected double ecalRe=1410; //ECAL outer R
- protected double hcalRs=1410; //HCAL inner R
- protected double hcalRe=2362; //HCAL outer R
- protected double ecalZs=1680; //ECAL inner Z
- protected double ecalZe=1820; //ECAL outer Z
- protected double hcalZs=1820; //HCAL inner Z
- protected double hcalZe=2772; //HCAL inner Z
-
public ShowerPointFinder(LocalHelixExtrapolator findCluster, Map<Track, Set<Cluster>> newMapTrackToShowerComponents, List<Cluster> mips, Set<CalorimeterHit> allHits) {
m_findCluster = findCluster;
m_newMapTrackToShowerComponents = newMapTrackToShowerComponents;
@@ -55,7 +45,6 @@
public Map<Track, BasicCluster> findMips(List<Track> tracksSortedByMomentum) {
Map<Track, BasicCluster> MapTrkToMIP = new HashMap<Track, BasicCluster>();
- BasicCluster mipclus = new BasicCluster();
for(Track tr : tracksSortedByMomentum){
Set<Cluster> clusoftrk = m_newMapTrackToShowerComponents.get(tr);
Set<Cluster> mipsoftrk = new HashSet<Cluster>();
@@ -64,7 +53,8 @@
mipsoftrk.add(cl);
}
}
-
+
+ BasicCluster mipclus = new BasicCluster();
if(mipsoftrk.size() == 0) {//No any mip found
MapTrkToMIP.put(tr,mipclus);
continue; //go to next track
@@ -123,7 +113,9 @@
}
if(debug) System.out.println("Debug: The end of the mip" );
- for( int i = 0 ; i < 70 ; i++){ // 70? until finding shower starting point
+ int lastIterationWithFoundHit = -1;
+ int maxSkippedLayers = 3;
+ for( int iIteration = 0 ; iIteration < 70 ; iIteration++){ // 70? until finding shower starting point
//To get the information from last hit
CalorimeterHit hit = mipclus.getCalorimeterHits().get( mipclus.getCalorimeterHits().size()-1);
IDDecoder id = hit.getIDDecoder();
@@ -134,45 +126,69 @@
double p = Math.sqrt(x*x+y*y+z*z);
double r = Math.sqrt(x*x+y*y);
double lastlayer = id.getLayer();
- long[] lastlayer5x5x2 = id.getNeighbourIDs(1, 2, 2);
+ String lastsubdetName = hit.getSubdetector().getName();
+ // Check for neighbouring hits in a 5x5x(2n+1) grid
+ // where n is the number of iterations since we last saw a hit.
+ // For example, if we saw a hit in the previous layer/iteration
+ // it'll be a 5x5x3 grid.
+ Set<Long> nearIDs = new HashSet<Long>();
+ int countIterationsSinceLastHit = (iIteration - lastIterationWithFoundHit);
+ if (countIterationsSinceLastHit <= 0) { throw new AssertionError("Book-keeping error"); }
+ if (countIterationsSinceLastHit > maxSkippedLayers+1) { throw new AssertionError("Book-keeping error"); }
+ long[] nearbyHitArray = id.getNeighbourIDs(countIterationsSinceLastHit, 2, 2);
+ for (long lastnearID : nearbyHitArray ) { nearIDs.add(lastnearID); }
+
Hep3Vector last0 = hitPosition(hit);
CalorimeterHit hit1 = mipclus.getCalorimeterHits().get( mipclus.getCalorimeterHits().size()-2);
Hep3Vector last1 = hitPosition(hit1);
- //Searh next hit in next layer and arrange in angular order
+ //Searh cur hit in cur layer and arrange in angular order
SortedMap<Double, CalorimeterHit> sortedHitbyPos= new TreeMap();
- for (CalorimeterHit hitall : m_allhits){
- IDDecoder idall = hitall.getIDDecoder();
- idall.setID(hitall.getCellID());
- double xall = idall.getPosition()[0];
- double yall = idall.getPosition()[1];
- double zall = idall.getPosition()[2];
- double hitmag = Math.sqrt(xall*xall+yall*yall+zall*zall);
- Hep3Vector hitpos = new BasicHep3Vector(xall,yall,zall);
- Set<Long> nearIDs = new HashSet<Long>();
- for (long lastnearID : lastlayer5x5x2 ) nearIDs.add(lastnearID);
- //use geometric information and layer number to find next hit
- boolean EcalToHcal= lastlayer == 30 && hitmag > p && idall.getLayer() == 0;
- boolean InSide = hitmag > p && nearIDs.contains(hitall.getCellID());
- boolean EndToBarrelEcal = ((r < ecalRs && r > ecalRs-4) && hitmag > p && ((ecalZs < z) || (z < ecalZe)) && (idall.getLayer()==0));//try
- boolean EndToBarrelHcal = ((r < ecalRs && r > ecalRs-4) && hitmag > p && ((hcalZs < z) || (z < hcalZe)) && (idall.getLayer()==0));
-
+ for (CalorimeterHit curhit : m_allhits){
+ IDDecoder curid = curhit.getIDDecoder();
+ curid.setID(curhit.getCellID());
+ double curx = curid.getPosition()[0];
+ double cury = curid.getPosition()[1];
+ double curz = curid.getPosition()[2];
+ double curp = Math.sqrt(curx*curx+cury*cury+curz*curz);
+ double curr = Math.sqrt(curx*curx+cury*cury);
+ Hep3Vector curpos = new BasicHep3Vector(curx,cury,curz);
+ double curlayer = curid.getLayer();
+ String cursubdetName = curhit.getSubdetector().getName();
+ //use geometric information and layer number to find cur hit
+ boolean EcalToHcal= lastlayer == 30 && curp > p && curid.getLayer() == 0;
+ boolean InSide = (curp > p || curlayer > lastlayer) && nearIDs.contains(curhit.getCellID());
+ boolean EndToBarrelEcal = (curid.getLayer()==0 && cursubdetName == "EMBarrel") && lastsubdetName=="EMEndcap";
+ boolean EndToBarrelHcal = (curid.getLayer()==0 && cursubdetName == "HADBarrel") && lastsubdetName=="HADEndcap";
if(debug){
- if(EndToBarrelEcal || EndToBarrelHcal ) System.out.println("This is a boundary condition");
+ if(EndToBarrelEcal || EndToBarrelHcal) System.out.println("This is a boundary condition");
}
- if(EcalToHcal || InSide || EndToBarrelEcal || EndToBarrelHcal){
- Hep3Vector test = VecOp.sub(hitpos, last0);
- Hep3Vector test_unit = VecOp.unit(test);
- Hep3Vector lastmip = VecOp.sub(last0, last1);
- Hep3Vector lastmip_unit = VecOp.unit(lastmip);
- double a = Math.acos(VecOp.dot(lastmip_unit, test_unit));
- sortedHitbyPos.put(a, hitall);
+ if(EcalToHcal || InSide || EndToBarrelEcal || EndToBarrelHcal ){
+ Hep3Vector cur = VecOp.sub(curpos, last0);
+ Hep3Vector curUnit = VecOp.unit(cur);
+ Hep3Vector last = VecOp.sub(last0, last1);
+ Hep3Vector lastUnit = VecOp.unit(last);
+ double a = Math.acos(VecOp.dot(lastUnit, curUnit));
+ double d = VecOp.sub(cur,last).magnitude();
+ if(a < 1.57 && d < 100){ //next hit in the range of angle 90 degree and less than 100mm(?)
+ sortedHitbyPos.put(a, curhit);
+ }
}
}
- if(sortedHitbyPos.size() == 0) continue; // go to next layer
+ if(sortedHitbyPos.size() == 0) {
+ // No hit found in this layer.
+ if (countIterationsSinceLastHit <= maxSkippedLayers) {
+ continue; // go to next layer
+ } else {
+ break; // stop -- ran out of hits and haven't seen any for several iterations
+ }
+ } else {
+ // Found one or more hits
+ lastIterationWithFoundHit = iIteration;
+ }
CalorimeterHit hitAdd = sortedHitbyPos.get(sortedHitbyPos.firstKey()); //add the hit with a smallest angle
examPosition(hitAdd);
mipclus.addHit(hitAdd);
@@ -199,9 +215,9 @@
}
System.out.println("MIP hit size= " + mipclus.getCalorimeterHits().size());
}
-
- MapTrkToMIP.put(tr,mipclus);
+ MapTrkToMIP.put(tr, mipclus);
} //track loop
+
return MapTrkToMIP;
}
@@ -218,7 +234,7 @@
for(long neighbourID : neighbours)
if (allHitIDs.contains(neighbourID)) count[i]++;
}
- if ( count[3] >=2 ) exam++ ; // if number of neighbour hits in 9x9 is larger than 1, count 1
+ if ( count[1] >=2 ) exam++ ; // if number of neighbour hits in 5x5 is larger than 1, count 1
else exam = 0; // if number of neighbour is none, set to zero and start again from 0
if(debug){