lcsim/src/org/lcsim/contrib/uiowa
diff -u -r1.2 -r1.3
--- ShowerPointFinder.java 10 Jul 2008 22:18:43 -0000 1.2
+++ ShowerPointFinder.java 12 Jul 2008 01:20:34 -0000 1.3
@@ -27,32 +27,23 @@
import org.lcsim.geometry.*;
import org.lcsim.util.decision.*;
- public class ShowerPointFinder{
+ public class ShowerPointFinder{
protected LocalHelixExtrapolator m_findCluster;
protected Map<Track, Set<Cluster>> m_newMapTrackToShowerComponents;
protected List<Cluster> m_mips;
- protected Set<CalorimeterHit> m_allhits;
- protected int exam = 0;
- protected boolean nomip = false;
- protected boolean debug = true;
- protected double m_distance = 0;
- protected double m_angle = 0;
- protected double m_degree = 0;
- protected double m_costheta = 0;
- protected double m_trans = 0;
+ protected Set<CalorimeterHit> m_allhits;
+ protected int exam = 0; //to test neighbour
+ protected boolean debug = false;
- //sid01
+ //sid01
protected double ecalRs=1270; //ECAL inner R
- protected double ecalRe=1410; //ECAL outer 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
-
- protected EventHeader m_event;
+ 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;
@@ -62,190 +53,189 @@
}
- public Map<Track, BasicCluster> findMips(List<Track> tracksSortedByMomentum) {
- Map<Track, BasicCluster> MapTrkToMIP = new HashMap<Track, BasicCluster>();
- for(Track tr : tracksSortedByMomentum){
- Set<Cluster> clusoftrk = m_newMapTrackToShowerComponents.get(tr);
- Set<Cluster> mipsoftrk = new HashSet<Cluster>();
- for(Cluster cl : clusoftrk)
- if(m_mips.contains(cl))
- mipsoftrk.add(cl);
- if(mipsoftrk.size() == 0) continue; //no any mip found
-
- if(debug) System.out.println("Given track");
-
- //new cluster for finding shower point
- BasicCluster mipclus = new BasicCluster();
-
- Hep3Vector interceptPoint = m_findCluster.performExtrapolation(tr);
- Hep3Vector tangent = m_findCluster.getTangent();
- Hep3Vector tangentUnit = VecOp.unit(tangent);
+ 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>();
+ for(Cluster cl : clusoftrk){
+ if(m_mips.contains(cl)) {
+ mipsoftrk.add(cl);
+ }
+ }
+
+ if(mipsoftrk.size() == 0) {//No any mip found
+ MapTrkToMIP.put(tr,mipclus);
+ continue; //go to next track
+ }
+
+ //new cluster for finding shower point
+ Hep3Vector interceptPoint = m_findCluster.performExtrapolation(tr); //not need except debuging
+ //Hep3Vector tangent = m_findCluster.getTangent();
+ //Hep3Vector tangentUnit = VecOp.unit(tangent);
- Hep3Vector trackThreeMomentum = new BasicHep3Vector(tr.getMomentum());
- double trackMomentum = trackThreeMomentum.magnitude();
+ if(debug){
+ Hep3Vector v = new BasicHep3Vector(tr.getMomentum());
+ System.out.println("Given track");
+ System.out.println("Debug: Size of Cluster of track= " + clusoftrk.size());
+ System.out.println("Debug: Size of Cluster of Mip of track= " + mipsoftrk.size());
+ System.out.println("Debug: trk momentum= " + v.magnitude());
+ if(interceptPoint != null) System.out.println("Debug: extra trk: pos= " + interceptPoint.magnitude() + " ("+ interceptPoint.x() + " " + interceptPoint.y() + " " + interceptPoint.z() + ")" );
+ }
+
+ //find closest mip to track
+ SortedMap<Double, Cluster> mipsSortedByDistance = new TreeMap();
+ for(Cluster cl : mipsoftrk){
+ CalorimeterHit firstHit = cl.getCalorimeterHits().get(0);
+ Hep3Vector first = hitPosition(firstHit);
+ double distance = first.magnitude();
+ mipsSortedByDistance.put(distance, cl);
+ }
- if(debug){
- System.out.println("Debug: Size of Cluster of track= " + clusoftrk.size());
- System.out.println("Debug: Size of Cluster of Mip of track= " + mipsoftrk.size());
- System.out.println("Debug: trk momentum= " + trackMomentum);
- if(interceptPoint != null) System.out.println("Debug: extra trk: pos= " + interceptPoint.magnitude() + " ("+ interceptPoint.x() + " " + interceptPoint.y() + " " + interceptPoint.z() + ")" );
- }
-
- //find closest mip to track
- SortedMap<Double, Cluster> mipsSortedByDistance = new TreeMap();
- for(Cluster cl : mipsoftrk){
- CalorimeterHit firstHit = cl.getCalorimeterHits().get(0);
- Hep3Vector first = hitPosition(firstHit);
- if (interceptPoint != null) {
- Hep3Vector displacement = VecOp.sub(first, interceptPoint);
- double distance = displacement.magnitude();
- mipsSortedByDistance.put(distance, cl);
- }
- }
-
- //arrange by the length of the hit position
- SortedMap<Double, CalorimeterHit> sortedMIPbyPos = new TreeMap();
- Cluster cl = mipsSortedByDistance.get(mipsSortedByDistance.firstKey());
+ //arrange by the length of the hit position
+ SortedMap<Double, CalorimeterHit> sortedMIPbyPos = new TreeMap();
+ Cluster cl = mipsSortedByDistance.get(mipsSortedByDistance.firstKey());
for( int i = 0 ; i < cl.getCalorimeterHits().size() ; i++){
- CalorimeterHit hit=cl.getCalorimeterHits().get(i);
- Hep3Vector v = hitPosition(hit);
- double pmag = v.magnitude();
+ CalorimeterHit hit=cl.getCalorimeterHits().get(i);
+ Hep3Vector v = hitPosition(hit);
+ double pmag = v.magnitude();
sortedMIPbyPos.put(pmag, hit);
- }
+ }
-
- Set<Long> allHitIDs = new HashSet<Long>();
- for (CalorimeterHit hit : m_allhits) allHitIDs.add(hit.getCellID());
-
- //String lastDet = new String();
-
- if(debug) System.out.println("Debug: Starting mip");
- Iterator itermip = sortedMIPbyPos.entrySet().iterator();
- for(int i=0 ; i < sortedMIPbyPos.size() ;i++){
- Map.Entry entry = (Map.Entry) itermip.next();
- Double key = (Double)entry.getKey();
- CalorimeterHit hit = (CalorimeterHit)entry.getValue();
- examPosition(hit,allHitIDs);
- mipclus.addHit(hit);
- }
- if(debug) System.out.println("Debug: The end of the mip" );
-
- for( int i = 0 ; i < 100 ; i++){ // 100? 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();
- id.setID(hit.getCellID());
- double x = id.getPosition()[0];
- double y = id.getPosition()[1];
- double z = id.getPosition()[2];
- double p = Math.sqrt(x*x+y*y+z*z);
- double r = Math.sqrt(x*x+y*y);
- double lastlayer = id.getLayer();
- //lastDet = hit.getSubdetector().getName();
- long[] lastlayer5x5x2 = id.getNeighbourIDs(1, 2, 2);
- 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
- 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));
- boolean EndToBarrelHcal = ((r < ecalRs && r > ecalRs-4) && hitmag > p && ((hcalZs < z) || (z < hcalZe)) && (idall.getLayer()==0));
+ //test if the mip is enough closet to track otherwise take it as no mip.
+ CalorimeterHit firstmip = sortedMIPbyPos.get(sortedMIPbyPos.firstKey());
+ IDDecoder firstmipid = firstmip.getIDDecoder();
+ firstmipid.setID(firstmip.getCellID());
+ if(firstmipid.getLayer() > 3) { // if the mip is not in first 3 layers, put 0 in mipclus
+ MapTrkToMIP.put(tr,mipclus);
+ continue; //go to next track
+ }
- if(debug){
- 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(debug) System.out.println("Debug: Starting mip");
+ Iterator itermip = sortedMIPbyPos.entrySet().iterator();
+ for(int i=0 ; i < sortedMIPbyPos.size() ;i++){
+ Map.Entry entry = (Map.Entry) itermip.next();
+ Double key = (Double)entry.getKey();
+ CalorimeterHit hit = (CalorimeterHit)entry.getValue();
+ examPosition(hit);
+ mipclus.addHit(hit);
+ }
+ if(debug) System.out.println("Debug: The end of the mip" );
- }
-
- if(sortedHitbyPos.size() == 0) continue; // add next hit if eligible hit in next layer.
- CalorimeterHit hitAdd = sortedHitbyPos.get(sortedHitbyPos.firstKey()); //add the hit with a smallest angle
- examPosition(hitAdd,allHitIDs);
- mipclus.addHit(hitAdd);
-
- if (exam > 2) {
- if(debug) System.out.println("Debug: showering is already started " + exam + " layers ago" );
- if(mipclus.getCalorimeterHits().size() < 5) break; // if hit size is fewer than 5, don't remove
- else for(int k=0; k < 3 ; k++) {
- CalorimeterHit toRemove = mipclus.getCalorimeterHits().get(mipclus.getCalorimeterHits().size()-1);
- mipclus.removeHit(toRemove);
- }
- break;
- }
- }
-
- if(debug){
- int k=0;
- if(mipclus.getCalorimeterHits().size() < 4) k = mipclus.getCalorimeterHits().size();
- else k = 4;
- for(int i =0; i < k ;i++){
- CalorimeterHit hit = mipclus.getCalorimeterHits().get( mipclus.getCalorimeterHits().size()-1-i);
- Hep3Vector v = hitPosition(hit);
- System.out.println("Debug: last "+i+" pos= "+v.magnitude()+" (" +v.x() + " " + v.y() + " " + v.z()+ ")" );
- }
- }
-
- if(debug) System.out.println("MIP hit size= " + mipclus.getCalorimeterHits().size());
-
- CalorimeterHit hit = mipclus.getCalorimeterHits().get(0);//first fit
- IDDecoder id = hit.getIDDecoder();
- id.setID(hit.getCellID());
- if(id.getLayer() <= 3) MapTrkToMIP.put(tr,mipclus); //remove possible secondary mip.
- }
- return MapTrkToMIP;
- }
-
- //give the hit position vector and examine neighbour hits
- protected void examPosition(CalorimeterHit hit, Set<Long> allHitIDs){
+ for( int i = 0 ; i < 70 ; i++){ // 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();
id.setID(hit.getCellID());
- int lay = id.getLayer();
- Hep3Vector v = hitPosition(hit);
- int p = (int) (v.magnitude());
- int count[] = new int[4];
- for (int i =0 ; i < 4 ; i++){
- long[] neighbours = id.getNeighbourIDs(0,i+1,i+1);
- for(long neighbourID : neighbours)
- if (allHitIDs.contains(neighbourID)) count[i]++;
+ double x = id.getPosition()[0];
+ double y = id.getPosition()[1];
+ double z = id.getPosition()[2];
+ 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);
+ 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
+ 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));
+
+ if(debug){
+ 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(sortedHitbyPos.size() == 0) continue; // go to next layer
+ CalorimeterHit hitAdd = sortedHitbyPos.get(sortedHitbyPos.firstKey()); //add the hit with a smallest angle
+ examPosition(hitAdd);
+ mipclus.addHit(hitAdd);
+
+ if (exam > 2) {
+ if(debug) System.out.println("Debug: showering is already started " + exam + " layers ago" );
+ if(mipclus.getCalorimeterHits().size() < 5) break; // if hit size is fewer than 5, don't remove
+ else for(int k=0; k < 3 ; k++) {
+ CalorimeterHit toRemove = mipclus.getCalorimeterHits().get(mipclus.getCalorimeterHits().size()-1);
+ mipclus.removeHit(toRemove);
+ }
+ break;
}
- if ( count[3] >=2 ) exam++ ; // if number of neighbour hits in 9x9 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){
- System.out.println(" hit: pos=" + p + " , " + count[0] + " hits in 3x3, " + count[1] + " hits in 5x5 " + count[2] + " hits in 7x7 " + count[3] + " hits in 9x9 " + lay + " layer of " + hit.getSubdetector().getName());
+ } //layer loop after mip
+
+ if(debug){
+ int k=0;
+ if(mipclus.getCalorimeterHits().size() < 4) k = mipclus.getCalorimeterHits().size();
+ else k = 4;
+ for(int i =0; i < k ;i++){
+ CalorimeterHit hit = mipclus.getCalorimeterHits().get( mipclus.getCalorimeterHits().size()-1-i);
+ Hep3Vector v = hitPosition(hit);
+ System.out.println("Debug: last "+i+" pos= "+v.magnitude()+" (" +v.x() + " " + v.y() + " " + v.z()+ ")" );
}
+ System.out.println("MIP hit size= " + mipclus.getCalorimeterHits().size());
+ }
+
+ MapTrkToMIP.put(tr,mipclus);
+ } //track loop
+ return MapTrkToMIP;
+ }
+
+ //examine neighbour hits
+ protected void examPosition(CalorimeterHit hit){
+ Set<Long> allHitIDs = new HashSet<Long>();
+ for (CalorimeterHit allhit : m_allhits) allHitIDs.add(allhit.getCellID());
+
+ IDDecoder id = hit.getIDDecoder();
+ id.setID(hit.getCellID());
+ int count[] = new int[4];
+ for (int i =0 ; i < 4 ; i++){
+ long[] neighbours = id.getNeighbourIDs(0,i+1,i+1);
+ 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
+ else exam = 0; // if number of neighbour is none, set to zero and start again from 0
- protected Hep3Vector hitPosition(CalorimeterHit hit){
- IDDecoder id = hit.getIDDecoder();
- id.setID(hit.getCellID());
- double x = id.getPosition()[0];
- double y = id.getPosition()[1];
- double z = id.getPosition()[2];
- Hep3Vector v = new BasicHep3Vector(x,y,z);
- return v;
+ if(debug){
+ Hep3Vector v = hitPosition(hit);
+ int p = (int) (v.magnitude());
+ System.out.println(" hit: pos=" + p + " , " + count[0] + " hits in 3x3, " + count[1] + " hits in 5x5 " + count[2] + " hits in 7x7 " + count[3] + " hits in 9x9 " + id.getLayer() + " layer of " + hit.getSubdetector().getName());
}
+ }
+
+ //give the vector of hit
+ protected Hep3Vector hitPosition(CalorimeterHit hit){
+ IDDecoder id = hit.getIDDecoder();
+ id.setID(hit.getCellID());
+ double x = id.getPosition()[0];
+ double y = id.getPosition()[1];
+ double z = id.getPosition()[2];
+ Hep3Vector v = new BasicHep3Vector(x,y,z);
+ return v;
+ }
}