Print

Print


Commit in lcsim/src/org/lcsim/contrib/SteveMagill on MAIN
ClusterToReconstructedParticleDriver.java+85-431.2 -> 1.3


lcsim/src/org/lcsim/contrib/SteveMagill
ClusterToReconstructedParticleDriver.java 1.2 -> 1.3
diff -u -r1.2 -r1.3
--- ClusterToReconstructedParticleDriver.java	18 May 2007 19:58:10 -0000	1.2
+++ ClusterToReconstructedParticleDriver.java	29 May 2008 21:31:15 -0000	1.3
@@ -1,4 +1,4 @@
-package org.lcsim.contrib.SteveMagill;
+package org.lcsim.contrib.compile.SteveMagill;
 
 import hep.aida.ITree;
 import hep.physics.vec.BasicHep3Vector;
@@ -6,7 +6,9 @@
 import hep.physics.particle.properties.ParticlePropertyManager;
 import hep.physics.particle.properties.ParticleType;
 import java.util.ArrayList;
-import java.util.List;
+// import java.util.List;
+import java.util.*;
+import org.lcsim.event.Track;
 import org.lcsim.event.CalorimeterHit;
 import org.lcsim.event.Cluster;
 import org.lcsim.event.EventHeader;
@@ -26,7 +28,7 @@
  *
  * Created on March 29, 2007, 1:52 PM
  *
- * $Id: ClusterToReconstructedParticleDriver.java,v 1.2 2007/05/18 19:58:10 magill Exp $
+ * $Id: ClusterToReconstructedParticleDriver.java,v 1.3 2008/05/29 21:31:15 magill Exp $
  */
 
 /**
@@ -41,7 +43,7 @@
 //    private FixedConeClusterer _fcc;
 //    private boolean _initialized;
     
-    private String[] _photonnames;
+    private String _photonnames;
     private String _tracknames;
     private String _neutralhadnames;
     private String _allrecoparticles;
@@ -49,6 +51,9 @@
     private double _neutralMass;
     private ParticleType _neutralParticleType;
     
+    ClusterEnergyCalculator phcalc;
+    ClusterEnergyCalculator nhcalc;
+    
     /** Creates a new instance of SimpleClusterToReconstructedParticleExampleDriver */
     public ClusterToReconstructedParticleDriver()
     {
@@ -75,60 +80,93 @@
 //            _initialized=true;
 //        }
         
+        // This rp list is for cheat clusters only (perfect all cal result)
+       // List<ReconstructedParticle> chrpList = new ArrayList<ReconstructedParticle>();
+       // List<BasicCluster> chclusters = event.get(BasicCluster.class,"PerfectCheatClusters");
+       // for (BasicCluster chclus : chclusters)
+       // {
+       //     double ClE = chclus.getEnergy();
+       //     double[] pos = chclus.getPosition();
+       //     double norm = sqrt(pos[0]*pos[0]+pos[1]*pos[1]+pos[2]*pos[2]);
+       //     double p = ClE/norm;
+       //     Hep3Vector mom = new BasicHep3Vector(p*pos[0], p*pos[1], p*pos[2]);
+       //     BaseReconstructedParticle rp = new BaseReconstructedParticle(ClE, mom);
+       //     rp.setMass(0); // all clusters have zero mass
+       //     ParticleType photype = ParticlePropertyManager.getParticlePropertyProvider().get(22);
+       //     ParticleID id = new BaseParticleID(photype);
+       //     rp.addParticleID(id);
+       //     rp.setParticleIdUsed(id);
+       //     rp.addCluster(chclus);
+       //     chrpList.add(rp);
+       // }
+        
         // array of ReconstructedParticles to add to the event...
         List<ReconstructedParticle> rpList = new ArrayList<ReconstructedParticle>();
         
-        // first, get tracks and add them to rpList
-        List<PerfectTrack> evtracks = event.get(PerfectTrack.class,_tracknames);
-        for (PerfectTrack itrack : evtracks)
+        //  define cluster energy calculators
+        phcalc = new QPhotonClusterEnergyCalculator();
+        nhcalc = new QNeutralHadronClusterEnergyCalculator();
+        
+        // first, get tracks and add them to rpList, use cluster E for energy, track P for momentum
+        List<Track> evtracks = event.get(Track.class,_tracknames);
+        Map<Track, BasicCluster> trclmap = (Map<Track, BasicCluster>) event.get("TrackClusMap");
+        for (Track itrack : evtracks)
         {
-            double TrE = itrack.getMCParticle().getEnergy();
+            BasicCluster trmclus = trclmap.get(itrack);
+            double TrclE = trmclus.getEnergy();
             double[] TrP = itrack.getMomentum();
+            double Trpmag = Math.sqrt(TrP[0]*TrP[0]+TrP[1]*TrP[1]+TrP[2]*TrP[2]);
+            double TrE = Math.sqrt(Trpmag*Trpmag+0.01948);
+            if (TrclE/Trpmag>2.0)
+            {
+                TrP[0] = TrclE*TrP[0]/Trpmag;
+                TrP[1] = TrclE*TrP[1]/Trpmag;
+                TrP[2] = TrclE*TrP[2]/Trpmag;
+                TrE = Math.sqrt(TrclE*TrclE+0.01948);
+            }
             Hep3Vector mom = new BasicHep3Vector(TrP[0],TrP[1],TrP[2]);
             BaseReconstructedParticle rp = new BaseReconstructedParticle(TrE, mom);
             rp.setMass(0.13957);
-            // rp.setMass(itrack.getMCParticle().getMass());
+//            rp.setMass(itrack.getMCParticle().getMass());
             double trq = (double) itrack.getCharge();
             rp.setCharge(trq);
-            ParticleID id = new BaseParticleID(itrack.getMCParticle().getType());
-            rp.addParticleID(id);
+            if (trq>0)
+            {
+                ParticleType pitype = ParticlePropertyManager.getParticlePropertyProvider().get(211);
+                ParticleID id = new BaseParticleID(pitype);
+                rp.addParticleID(id);
+                rp.setParticleIdUsed(id);
+            } 
+            if (trq<0)
+            {
+                ParticleType pitype = ParticlePropertyManager.getParticlePropertyProvider().get(-211);
+                ParticleID id = new BaseParticleID(pitype);
+                rp.addParticleID(id);
+                rp.setParticleIdUsed(id);
+            }
             rp.addTrack(itrack);
+            rp.addCluster(trmclus); 
             rpList.add(rp);
         }
         
         // now get photon clusters
-        List<BasicCluster> phobclusters = event.get(BasicCluster.class,"PhotonBClusters");
-        List<BasicCluster> phoecclusters = event.get(BasicCluster.class,"PhotonECClusters");
-        for (BasicCluster phobclus : phobclusters)
+        List<BasicCluster> phoclusters = event.get(BasicCluster.class,_photonnames);
+        for (BasicCluster phoclus : phoclusters)
         {
-            double ClE = phobclus.getEnergy()*1.015*1.013;
-            double[] pos = phobclus.getPosition();
+//            double ClE = phoclus.getEnergy();
+            double ClE = phcalc.getEnergy(phoclus);
+            double[] pos = phoclus.getPosition();
             double norm = sqrt(pos[0]*pos[0]+pos[1]*pos[1]+pos[2]*pos[2]);
             double p = ClE/norm;
             Hep3Vector mom = new BasicHep3Vector(p*pos[0], p*pos[1], p*pos[2]);
             BaseReconstructedParticle rp = new BaseReconstructedParticle(ClE, mom);
-            rp.setMass(0);
+            rp.setMass(0.);
+            rp.setCharge(0.);
             ParticleType photype = ParticlePropertyManager.getParticlePropertyProvider().get(22);
             ParticleID id = new BaseParticleID(photype);
             rp.addParticleID(id);
             rp.setParticleIdUsed(id);
-            rp.addCluster(phobclus);
-            rpList.add(rp);
-        }
-        for (BasicCluster phoecclus : phoecclusters)
-        {
-            double ClE = phoecclus.getEnergy()*1.015*1.013;
-            double[] pos = phoecclus.getPosition();
-            double norm = sqrt(pos[0]*pos[0]+pos[1]*pos[1]+pos[2]*pos[2]);
-            double p = ClE/norm;
-            Hep3Vector mom = new BasicHep3Vector(p*pos[0], p*pos[1], p*pos[2]);
-            BaseReconstructedParticle rp = new BaseReconstructedParticle(ClE, mom);
-            rp.setMass(0);
-            ParticleType photype = ParticlePropertyManager.getParticlePropertyProvider().get(22);
-            ParticleID id = new BaseParticleID(photype);
-            rp.addParticleID(id);
-            rp.setParticleIdUsed(id);
-            rp.addCluster(phoecclus);
+            rp.addCluster(phoclus);
             rpList.add(rp);
         }
         
@@ -138,13 +176,21 @@
         List<BasicCluster> neuhclusters = event.get(BasicCluster.class,_neutralhadnames);
         for (BasicCluster neuhclus : neuhclusters)
         {
-            double ClE = neuhclus.getEnergy()*1.3;
+            double ClEa = neuhclus.getEnergy();
+            double ClE = nhcalc.getEnergy(neuhclus);
+//            if (ClE<0.5) ClE = 0.5;
+            double p2mag = ClE*ClE-0.4977*0.4977;
+            if (p2mag<0) p2mag = 0;
+            double pmag = Math.sqrt(p2mag);
             double[] pos = neuhclus.getPosition();
             double norm = sqrt(pos[0]*pos[0]+pos[1]*pos[1]+pos[2]*pos[2]);
-            double p = ClE/norm;
+//            double p = ClE/norm;
+            double p = pmag/norm;
             Hep3Vector mom = new BasicHep3Vector(p*pos[0], p*pos[1], p*pos[2]);
             BaseReconstructedParticle rp = new BaseReconstructedParticle(ClE, mom);
-            rp.setMass(718.5);
+//            rp.setMass(0.);
+            rp.setMass(0.4977);  // use klong value
+            rp.setCharge(0.);
             ParticleType neutype = ParticlePropertyManager.getParticlePropertyProvider().get(130);
             ParticleID id = new BaseParticleID(neutype);
             rp.addParticleID(id);
@@ -160,11 +206,7 @@
         // add the list of rps to the event...
         event.put(_allrecoparticles,rpList, ReconstructedParticle.class,0);
         //some simple histograms...
-        aida.cloud1D("Number of ReconstructedParticles found").fill(rpList.size());
-        for(ReconstructedParticle rp : rpList)
-        {
-            aida.cloud1D("rp energy").fill(rp.getEnergy());
-        }
+        aida.cloud1D("Number of AllRecoParticles found").fill(rpList.size());
         
     }
     
@@ -172,7 +214,7 @@
   {
       _tracknames = tnames;
   }
-  public void setPhotonNames(String[] pnames)
+  public void setPhotonNames(String pnames)
   {
       _photonnames = pnames;
   }
CVSspam 0.2.8