7 modified files
lcsim/src/org/lcsim/contrib/Cassell/recon/analysis
diff -u -r1.1 -r1.2
--- CalClusterComposition.java 13 Sep 2007 15:10:11 -0000 1.1
+++ CalClusterComposition.java 16 Oct 2007 18:20:13 -0000 1.2
@@ -41,8 +41,13 @@
/** Creates a new instance of CalClusterComposition */
public CalClusterComposition(List<MCParticle> fsp)
{
+ this(fsp,null);
+ }
+ public CalClusterComposition(List<MCParticle> fsp,List<MCParticle> fst)
+ {
tr = new TraceOrigin(fsp);
pt = new ParticleType();
+ if(fst != null)pt.setTracks(fst);
ntypes = pt.getTypes().length;
pfrac = new double[ntypes];
np = new int[ntypes];
lcsim/src/org/lcsim/contrib/Cassell/recon/analysis
diff -u -r1.1 -r1.2
--- CalorimeterHitCollectionEnergies.java 13 Sep 2007 15:10:11 -0000 1.1
+++ CalorimeterHitCollectionEnergies.java 16 Oct 2007 18:20:13 -0000 1.2
@@ -27,7 +27,10 @@
String[] coll;
String pre;
String FSname;
+ String FSTrackedname;
String[] type;
+ String[] rtype;
+ String[] bltype;
int nmax = 10000000;
TraceOrigin to;
double[] totE;
@@ -36,23 +39,36 @@
ParticleType pt;
/** Creates a new instance of CalorimeterHitCollectionEnergies */
- public CalorimeterHitCollectionEnergies(String Prefix,String fsname,String[] collnames,boolean plot)
+ public CalorimeterHitCollectionEnergies(String Prefix,String fsname,String fstrackedname, String[] collnames,boolean plot)
{
coll = collnames;
pre = Prefix+"/";
FSname = fsname;
+ FSTrackedname = fstrackedname;
makeHists = plot;
pt = new ParticleType();
type = pt.getTypes();
+ rtype = pt.getRTypes();
+ bltype = pt.getBLTypes();
runE = new double[type.length];
}
public CalorimeterHitCollectionEnergies(String Prefix,String fsname,String[] collnames)
{
- this(Prefix,fsname,collnames,false);
+ this(Prefix,fsname,"All",collnames,false);
+ }
+ public CalorimeterHitCollectionEnergies(String Prefix,String fsname,String[] collnames,boolean plot)
+ {
+ this(Prefix,fsname,"All",collnames,plot);
}
protected void process(EventHeader event)
{
List<MCParticle> fslist = event.get(MCParticle.class,FSname);
+ List<MCParticle> fstracked = null;
+ if(FSTrackedname.compareTo("All") != 0)
+ {
+ fstracked = event.get(MCParticle.class,FSTrackedname);
+ pt.setTracks(fstracked);
+ }
to = new TraceOrigin(fslist);
List<List<CalorimeterHit>> all = event.get(CalorimeterHit.class);
totE = new double[type.length];
@@ -110,9 +126,21 @@
}
if(makeHists)
{
+ double[] rE = new double[rtype.length];
+ double[] blE = new double[bltype.length];
for(int i=0;i<type.length;i++)
{
aida.cloud1D(pre+type[i]+" cal energy per event",nmax).fill(totE[i]);
+ rE[pt.mapTypeToRType(i)] += totE[i];
+ blE[pt.mapRTypeToBLType(pt.mapTypeToRType(i))] += totE[i];
+ }
+ for(int i=0;i<rtype.length;i++)
+ {
+ aida.cloud1D(pre+rtype[i]+" cal energy per event",nmax).fill(rE[i]);
+ }
+ for(int i=0;i<bltype.length;i++)
+ {
+ aida.cloud1D(pre+bltype[i]+" cal energy per event",nmax).fill(blE[i]);
}
}
}
lcsim/src/org/lcsim/contrib/Cassell/recon/analysis
diff -u -r1.1 -r1.2
--- ClusterCollectionEnergies.java 13 Sep 2007 15:10:11 -0000 1.1
+++ ClusterCollectionEnergies.java 16 Oct 2007 18:20:13 -0000 1.2
@@ -28,7 +28,10 @@
String coll;
String pre;
String FSname;
+ String FSTrackedname;
String[] type;
+ String[] rtype;
+ String[] bltype;
int nmax = 10000000;
double[] totE;
double[] runE;
@@ -37,24 +40,37 @@
boolean makeHists;
/** Creates a new instance of ClusterCollectionEnergies */
- public ClusterCollectionEnergies(String fsname,String collname,boolean plot)
+ public ClusterCollectionEnergies(String fsname,String fstrackedname,String collname,boolean plot)
{
pt = new ParticleType();
type = pt.getTypes();
+ rtype = pt.getRTypes();
+ bltype = pt.getBLTypes();
coll = collname;
FSname = fsname;
+ FSTrackedname = fstrackedname;
makeHists = plot;
runE = new double[type.length];
pre = coll+"/";
}
public ClusterCollectionEnergies(String fsname,String collname)
{
- this(fsname,collname,false);
+ this(fsname,"All",collname,false);
+ }
+ public ClusterCollectionEnergies(String fsname,String collname,boolean plot)
+ {
+ this(fsname,"All",collname,plot);
}
protected void process(EventHeader event)
{
List<MCParticle> fslist = event.get(MCParticle.class,FSname);
- ccc = new CalClusterComposition(fslist);
+ List<MCParticle> fstracked = null;
+ if(FSTrackedname.compareTo("All") != 0)
+ {
+ fstracked = event.get(MCParticle.class,FSTrackedname);
+ ccc = new CalClusterComposition(fslist,fstracked);
+ }
+ else ccc = new CalClusterComposition(fslist);
totE = new double[type.length];
List<Cluster> cl = event.get(Cluster.class,coll);
for(Cluster c:cl)
@@ -86,12 +102,24 @@
if(makeHists)
{
double evte = 0.;
+ double[] rE = new double[rtype.length];
+ double[] blE = new double[bltype.length];
for(int i=0;i<type.length;i++)
{
evte += totE[i];
aida.cloud1D(pre+"Energy per event from "+type[i], nmax).fill(totE[i]);
+ rE[pt.mapTypeToRType(i)] += totE[i];
+ blE[pt.mapRTypeToBLType(pt.mapTypeToRType(i))] += totE[i];
}
aida.cloud1D(pre+"Total Energy per event", nmax).fill(evte);
+ for(int i=0;i<rtype.length;i++)
+ {
+ aida.cloud1D(pre+"Energy per event from "+rtype[i], nmax).fill(rE[i]);
+ }
+ for(int i=0;i<bltype.length;i++)
+ {
+ aida.cloud1D(pre+"Energy per event from "+bltype[i], nmax).fill(blE[i]);
+ }
}
}
public double[] getEnergies(){return totE;}
lcsim/src/org/lcsim/contrib/Cassell/recon/analysis
diff -u -r1.1 -r1.2
--- ExampleAnalysisDriver.java 13 Sep 2007 15:10:11 -0000 1.1
+++ ExampleAnalysisDriver.java 16 Oct 2007 18:20:13 -0000 1.2
@@ -29,6 +29,7 @@
String GenFSname = "GenFinalStateParticles";
String CheatReconRname = "ReconPerfectReconParticles";
String CheatReconFSname = "ReconFSParticles";
+ String CheatReconFSTrackedname = "ReconFSParticlesTracked";
String PPRPflowRname = "PPRReconParticles";
String[] DigisimCollectionNames =
{"EcalBarrDigiHits", "EcalEndcapDigiHits","HcalBarrDigiHits", "HcalEndcapDigiHits"};
@@ -45,7 +46,7 @@
ClusterCollectionEnergies[] PPRClusterEnergies;
ReconstructedParticleCollectionEnergies PerfectReconEnergies;
ReconstructedParticleCollectionEnergies PPRReconEnergies;
-
+
String[] type;
ParticleType pt;
@@ -62,7 +63,7 @@
crd.setCheatFSParticleOutputName(CheatReconFSname);
add(crd);
//
-// Make the perfect pattern recognition pflow reconstructed particles from the
+// Make the perfect pattern recognition pflow reconstructed particles from the
// Cheat Recon particles
//
add(new PPRDriver(CheatReconRname,PPRPflowRname));
@@ -76,16 +77,16 @@
//
GenFSEnergies = new FSParticleCollectionEnergies(GenFSname,true);
add(GenFSEnergies);
- ReconFSEnergies = new FSParticleCollectionEnergies(CheatReconFSname,true);
+ ReconFSEnergies = new FSParticleCollectionEnergies(CheatReconFSname,CheatReconFSTrackedname,true);
add(ReconFSEnergies);
AllCalHitEnergies = new CalorimeterHitCollectionEnergies(
- "AllCalHits",CheatReconFSname,DigisimCollectionNames,true);
+ "AllCalHits",CheatReconFSname,CheatReconFSTrackedname,DigisimCollectionNames,true);
add(AllCalHitEnergies);
PerfectReconEnergies = new ReconstructedParticleCollectionEnergies(
- CheatReconFSname,CheatReconRname, true);
+ CheatReconFSname,CheatReconFSTrackedname,CheatReconRname, true);
add(PerfectReconEnergies);
PPRReconEnergies = new ReconstructedParticleCollectionEnergies(
- CheatReconFSname,PPRPflowRname, true);
+ CheatReconFSname,CheatReconFSTrackedname,PPRPflowRname, true);
add(PPRReconEnergies);
//
@@ -101,10 +102,10 @@
ClusterFromCheatReconRname[i] = CheatReconRname+type[i]+"Clusters";
ClusterFromPPRPflowRname[i] = PPRPflowRname+type[i]+"Clusters";
PerfectClusterEnergies[i] = new ClusterCollectionEnergies(
- CheatReconFSname,ClusterFromCheatReconRname[i],true);
+ CheatReconFSname,CheatReconFSTrackedname,ClusterFromCheatReconRname[i],true);
add(PerfectClusterEnergies[i]);
PPRClusterEnergies[i] = new ClusterCollectionEnergies(
- CheatReconFSname,ClusterFromPPRPflowRname[i],true);
+ CheatReconFSname,CheatReconFSTrackedname,ClusterFromPPRPflowRname[i],true);
add(PPRClusterEnergies[i]);
}
ievt = 0;
@@ -119,26 +120,26 @@
double Zmass = 0.;
for(MCParticle p:mcl)
{
- int id = Math.abs(p.getPDGID());
- if( (id == 1 )||(id == 2 )||(id == 3 ) )
- {
- if(p.getParents().get(0).getPDGID() == 23 )
- {
- Zmass = p.getParents().get(0).getMass();
- Hep3Vector pp = p.getMomentum();
- double cost = Math.abs(pp.z()/pp.magnitude());
- if(cost > .8)
- {
- keepit = false;
- continue;
- }
- }
- }
+ int id = Math.abs(p.getPDGID());
+ if( (id == 1 )||(id == 2 )||(id == 3 ) )
+ {
+ if(p.getParents().get(0).getPDGID() == 23 )
+ {
+ Zmass = p.getParents().get(0).getMass();
+ Hep3Vector pp = p.getMomentum();
+ double cost = Math.abs(pp.z()/pp.magnitude());
+ if(cost > .8)
+ {
+ keepit = false;
+ continue;
+ }
+ }
+ }
}
if(!keepit)
{
- ievt++;
- return;
+ ievt++;
+ return;
}
System.out.println("Processing event "+ievt);
super.process(event);
@@ -149,8 +150,12 @@
String pre = PPRPflowRname+" to "+CheatReconRname+"/";
double[] Echeat = PerfectReconEnergies.getEnergies();
double[] Eppr = PPRReconEnergies.getEnergies();
+ double EcheatTot = 0.;
+ double EpprTot = 0.;
for(int i=0;i<type.length;i++)
{
+ EcheatTot += Echeat[i];
+ EpprTot += Eppr[i];
if(Echeat[i] > 1.)
{
double dE = Eppr[i] - Echeat[i];
@@ -159,6 +164,10 @@
aida.cloud1D(pre+"delta E over rootE "+type[i], nmax).fill(dE/Math.sqrt(Echeat[i]));
}
}
+ double dE = EpprTot - EcheatTot;
+ aida.cloud1D(pre+"delta E ", nmax).fill(dE);
+ aida.cloud1D(pre+"delta E over E ", nmax).fill(dE/EcheatTot);
+ aida.cloud1D(pre+"delta E over rootE ", nmax).fill(dE/Math.sqrt(EcheatTot));
ievt++;
}
lcsim/src/org/lcsim/contrib/Cassell/recon/analysis
diff -u -r1.1 -r1.2
--- FSParticleCollectionEnergies.java 13 Sep 2007 15:10:11 -0000 1.1
+++ FSParticleCollectionEnergies.java 16 Oct 2007 18:20:13 -0000 1.2
@@ -23,8 +23,11 @@
{
private AIDA aida = AIDA.defaultInstance();
String FSname;
+ String FSTrackedname;
String pre;
String[] type;
+ String[] rtype;
+ String[] bltype;
double[] totE;
double[] runE;
int nmax = 10000000;
@@ -32,24 +35,38 @@
boolean makeHists;
double mass;
double nonmass;
+ double nonEtot;
/** Creates a new instance of FSParticleCollectionEnergies */
- public FSParticleCollectionEnergies(String fsname, boolean hist)
+ public FSParticleCollectionEnergies(String fsname,String fstrackedname, boolean hist)
{
FSname = fsname;
+ FSTrackedname = fstrackedname;
pt = new ParticleType();
type = pt.getTypes();
+ rtype = pt.getRTypes();
+ bltype = pt.getRTypes();
pre = FSname + "/";
makeHists = hist;
runE = new double[type.length];
}
public FSParticleCollectionEnergies(String fsname)
{
- this(fsname, false);
+ this(fsname,"All", false);
+ }
+ public FSParticleCollectionEnergies(String fsname, boolean hist)
+ {
+ this(fsname,"All", hist);
}
protected void process(EventHeader event)
{
List<MCParticle> fslist = event.get(MCParticle.class,FSname);
+ List<MCParticle> fstracked = null;
+ if(FSTrackedname.compareTo("All") != 0)
+ {
+ fstracked = event.get(MCParticle.class,FSTrackedname);
+ pt.setTracks(fstracked);
+ }
totE = new double[type.length];
double px = 0.;
double py = 0.;
@@ -75,12 +92,30 @@
if(makeHists)aida.histogram1D(pre+type[itype]+" Energy distibution of particles",200,0.,100.).fill(p.getEnergy());
}
double Etot = 0.;
- double nonEtot = 0.;
+ nonEtot = 0.;
+ double[] rE = new double[rtype.length];
+ double[] blE = new double[bltype.length];
for(int i=0;i<type.length;i++)
{
if(makeHists)aida.cloud1D(pre+type[i]+" Energy distibution per event",nmax).fill(totE[i]);
Etot += totE[i];
- if(type[i].compareTo("neutrino") != 0)nonEtot += totE[i];
+ if(type[i].compareTo("neutrino") != 0)
+ {
+ nonEtot += totE[i];
+ rE[pt.mapTypeToRType(i)] += totE[i];
+ blE[pt.mapRTypeToBLType(pt.mapTypeToRType(i))] += totE[i];
+ }
+ }
+ if(makeHists)
+ {
+ for(int i=0;i<rtype.length;i++)
+ {
+ aida.cloud1D(pre+rtype[i]+" Energy distibution per event",nmax).fill(rE[i]);
+ }
+ for(int i=0;i<bltype.length;i++)
+ {
+ aida.cloud1D(pre+bltype[i]+" Energy distibution per event",nmax).fill(blE[i]);
+ }
}
mass = Math.sqrt(Etot*Etot - px*px - py*py - pz*pz);
nonmass = Math.sqrt(nonEtot*nonEtot - nonpx*nonpx - nonpy*nonpy - nonpz*nonpz);
@@ -98,6 +133,8 @@
{return mass;}
public double getNoNeutrinoMass()
{return nonmass;}
+ public double getNoNeutrinoEnergy()
+ {return nonEtot;}
public double[] getRunEnergies()
{return runE;}
lcsim/src/org/lcsim/contrib/Cassell/recon/analysis
diff -u -r1.1 -r1.2
--- ParticleType.java 13 Sep 2007 15:10:11 -0000 1.1
+++ ParticleType.java 16 Oct 2007 18:20:13 -0000 1.2
@@ -9,6 +9,8 @@
package org.lcsim.contrib.Cassell.recon.analysis;
import org.lcsim.event.MCParticle;
import org.lcsim.event.ReconstructedParticle;
+import org.lcsim.event.base.MCReconstructedParticle;
+import java.util.List;
/**
*
@@ -17,14 +19,25 @@
public class ParticleType
{
- String[] type = {"photon","ch hadron","n hadron","electron","neutrino"};
- String[] rtype = {"photon","charged","n hadron"};
- int[] ttrtype = {0,1,2,1,2};
+ String[] type = {"photon","ch hadron","n hadron","electron","neutrino","untracked electron","untracked ch hadron"};
+ String[] rtype = {"recon photon","recon track","recon n hadron"};
+ String[] bltype = {"neutral","charged"};
+ int[] ttrtype = {0,1,2,1,2,0,2};
+ int[] rtbltype = {0,1,0};
+ List<MCParticle> tracks = null;
/** Creates a new instance of ParticleType */
public ParticleType()
{
}
+ public ParticleType(List<MCParticle> tr)
+ {
+ tracks = tr;
+ }
+ public void setTracks(List<MCParticle> tr)
+ {
+ tracks = tr;
+ }
//
// Return photon=0, charged hadron=1, neutral hadron=2, electron=3, neutrino=4
//
@@ -44,13 +57,54 @@
{
itype = 1;
if(Math.abs(pdg) == 11)itype = 3;
+ if(tracks != null)
+ {
+ boolean tracked = false;
+ for(MCParticle tp:tracks)
+ {
+ if(p == tp)
+ {
+ tracked = true;
+ break;
+ }
+ }
+ if(!tracked)
+ {
+ if(itype == 1)itype = 6;
+ else itype = 5;
+ }
+ }
}
}
return itype;
}
public int getRType(ReconstructedParticle p)
{
- if(p.getCharge() != 0.)return 1;
+ if(p.getCharge() != 0.)
+ {
+ if(tracks == null)return 1;
+ if(p instanceof MCReconstructedParticle)
+ {
+ MCReconstructedParticle rp = (MCReconstructedParticle) p;
+ MCParticle pp = rp.getMCParticle();
+ boolean tracked = false;
+ for(MCParticle tp:tracks)
+ {
+ if(pp == tp)
+ {
+ tracked = true;
+ break;
+ }
+ }
+ if(!tracked)
+ {
+ if(Math.abs(pp.getPDGID()) == 11)return 0;
+ return 2;
+ }
+ else return 1;
+ }
+ else return 1;
+ }
if(p.getMass() == 0)return 0;
return 2;
}
@@ -62,9 +116,18 @@
{
return rtype;
}
+ public String[] getBLTypes()
+ {
+ return bltype;
+ }
public int mapTypeToRType(int itype)
{
if( (itype >= 0)&&(itype < type.length) )return ttrtype[itype];
return -1;
}
+ public int mapRTypeToBLType(int irtype)
+ {
+ if( (irtype >= 0)&&(irtype < rtype.length) )return rtbltype[irtype];
+ return -1;
+ }
}
lcsim/src/org/lcsim/contrib/Cassell/recon/analysis
diff -u -r1.2 -r1.3
--- ReconstructedParticleCollectionEnergies.java 14 Sep 2007 13:58:23 -0000 1.2
+++ ReconstructedParticleCollectionEnergies.java 16 Oct 2007 18:20:13 -0000 1.3
@@ -32,8 +32,10 @@
String coll;
String pre;
String FSname;
+ String FSTrackedname;
String[] type;
String[] rtype;
+ String[] bltype;
int nmax = 10000000;
TraceOrigin to;
double[] totE;
@@ -46,25 +48,37 @@
ParticleType pt;
/** Creates a new instance of ReconstructedParticleCollectionEnergies */
- public ReconstructedParticleCollectionEnergies(String fsname,String collname,boolean plot)
+ public ReconstructedParticleCollectionEnergies(String fsname,String fstrackedname, String collname,boolean plot)
{
coll = collname;
pre = coll+"/";
FSname = fsname;
+ FSTrackedname = fstrackedname;
pt = new ParticleType();
type = pt.getTypes();
rtype = pt.getRTypes();
+ bltype = pt.getBLTypes();
runE = new double[rtype.length];
- runCalE = new double[type.length][rtype.length];
+ runCalE = new double[rtype.length][rtype.length];
makeHists = plot;
}
public ReconstructedParticleCollectionEnergies(String fsname,String collname)
{
- this(fsname,collname,false);
+ this(fsname,"All",collname,false);
+ }
+ public ReconstructedParticleCollectionEnergies(String fsname,String collname,boolean plot)
+ {
+ this(fsname,"All",collname,plot);
}
protected void process(EventHeader event)
{
List<MCParticle> fslist = event.get(MCParticle.class,FSname);
+ List<MCParticle> fstracked = null;
+ if(FSTrackedname.compareTo("All") != 0)
+ {
+ fstracked = event.get(MCParticle.class,FSTrackedname);
+ pt.setTracks(fstracked);
+ }
to = new TraceOrigin(fslist);
List<ReconstructedParticle> rpl = event.get(ReconstructedParticle.class,coll);
totE = new double[rtype.length];
@@ -162,8 +176,11 @@
mass = Math.sqrt(Etot*Etot - pxtot*pxtot - pytot*pytot - pztot*pztot);
if(makeHists)
{
+ double[][] blCalE = new double[bltype.length][bltype.length];
+ double[] blE = new double[bltype.length];
for(int i=0;i<rtype.length;i++)
{
+ blE[pt.mapRTypeToBLType(i)] += totE[i];
aida.cloud1D(pre+rtype[i]+" Energy distibution per event",nmax).fill(totE[i]);
double calE = 0.;;
for(int j=0;j<rtype.length;j++)
@@ -174,12 +191,28 @@
{
for(int j=0;j<rtype.length;j++)
{
+ blCalE[pt.mapRTypeToBLType(i)][pt.mapRTypeToBLType(j)] += totCalE[i][j];
aida.histogram1D(
pre+rtype[i]+":wted fraction of calE from"+rtype[j],
51, 0., 1.02).fill(totCalE[i][j]/calE, calE);
}
}
}
+ for(int i=0;i<bltype.length;i++)
+ {
+ aida.cloud1D(pre+bltype[i]+" Energy distibution per event",nmax).fill(blE[i]);
+ double calE = 0.;;
+ for(int j=0;j<bltype.length;j++)
+ {
+ calE += blCalE[i][j];
+ }
+ for(int j=0;j<bltype.length;j++)
+ {
+ aida.histogram1D(
+ pre+bltype[i]+":wted fraction of calE from"+bltype[j],
+ 51, 0., 1.02).fill(blCalE[i][j]/calE, calE);
+ }
+ }
aida.cloud1D(pre+"Event Mass",nmax).fill(mass);
aida.cloud1D(pre+"Event Energy",nmax).fill(Etot);
}
@@ -198,11 +231,15 @@
{return Etot;}
protected void endOfData()
{
+ System.out.println();
System.out.println(coll+" run summary");
System.out.println();
System.out.println(" Energy sums");
+ double[] blE = new double[bltype.length];
+ double[][] blCalE = new double[bltype.length][bltype.length];
for(int i=0;i<rtype.length;i++)
{
+ blE[pt.mapRTypeToBLType(i)] += runE[i];
System.out.println(rtype[i]+" = "+(int)runE[i]);
}
System.out.println();
@@ -212,11 +249,29 @@
String str = rtype[i];
for(int j=0;j<rtype.length;j++)
{
+ blCalE[pt.mapRTypeToBLType(i)][pt.mapRTypeToBLType(j)] += runCalE[i][j];
str += " from "+rtype[j]+" = "+(int)runCalE[i][j]+";";
}
System.out.println(str);
}
System.out.println();
+ System.out.println(" Bottom line Energy sums");
+ for(int i=0;i<bltype.length;i++)
+ {
+ System.out.println(bltype[i]+" = "+(int)blE[i]);
+ }
+ System.out.println();
+ System.out.println(" Bottom line Cal Energy sums");
+ for(int i=0;i<bltype.length;i++)
+ {
+ String str = bltype[i];
+ for(int j=0;j<bltype.length;j++)
+ {
+ str += " from "+bltype[j]+" = "+(int)blCalE[i][j]+";";
+ }
+ System.out.println(str);
+ }
+ System.out.println();
}
}
CVSspam 0.2.8