Author: [log in to unmask]
Date: Thu Apr 16 10:48:02 2015
New Revision: 2719
Log:
simple ana driver
Modified:
java/trunk/users/src/main/java/org/hps/users/phansson/TruthMomentumResolutionDriver.java
Modified: java/trunk/users/src/main/java/org/hps/users/phansson/TruthMomentumResolutionDriver.java
=============================================================================
--- java/trunk/users/src/main/java/org/hps/users/phansson/TruthMomentumResolutionDriver.java (original)
+++ java/trunk/users/src/main/java/org/hps/users/phansson/TruthMomentumResolutionDriver.java Thu Apr 16 10:48:02 2015
@@ -1,12 +1,14 @@
package org.hps.users.phansson;
import hep.aida.IAnalysisFactory;
-import hep.aida.ICloud1D;
import hep.aida.IDataPointSet;
import hep.aida.IDataPointSetFactory;
import hep.aida.IHistogram1D;
import hep.aida.IPlotter;
import hep.aida.IPlotterStyle;
+import hep.aida.ref.histogram.Histogram1D;
+import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.Hep3Vector;
import java.io.IOException;
import java.util.HashMap;
@@ -16,11 +18,13 @@
import java.util.logging.Logger;
import org.hps.analysis.ecal.HPSMCParticlePlotsDriver;
+import org.lcsim.constants.Constants;
import org.lcsim.event.EventHeader;
import org.lcsim.event.MCParticle;
import org.lcsim.event.Track;
import org.lcsim.event.base.ParticleTypeClassifier;
import org.lcsim.geometry.Detector;
+import org.lcsim.recon.tracking.seedtracker.SeedTrack;
import org.lcsim.util.Driver;
import org.lcsim.util.aida.AIDA;
@@ -36,9 +40,12 @@
int totalTracks = 0;
private boolean _debug = false;
private AIDA aida = AIDA.defaultInstance();
+ double bfield;
+ double bfac;
//private AIDAFrame pFrame;
IAnalysisFactory af = aida.analysisFactory();
IPlotter pPlotter;
+ IPlotter pPlotterProp;
IPlotter pPlotter2;
IPlotter pPlotter22;
IPlotter pPlotter3;
@@ -55,17 +62,19 @@
IHistogram1D hTruthMatchedNegTrackPdiff;
IHistogram1D[] hTruthMatchedNegTrackPdiffvsP = new IHistogram1D[5];
IDataPointSet hTruthMatchedNegTrackPres;
- ICloud1D hTruthMatchedPosTrackPdiffPrev[] = new ICloud1D[5];
- ICloud1D hTruthMatchedNegTrackPdiffPrev[] = new ICloud1D[5];
- ICloud1D hNTracks;
+ IHistogram1D hTruthMatchedPosTrackPdiffPrev[] = new Histogram1D[5];
+ IHistogram1D hTruthMatchedNegTrackPdiffPrev[] = new Histogram1D[5];
+ IHistogram1D hNTracks;
IHistogram1D trkCountVsEventPlot;
- ICloud1D hNPosTracks;
- ICloud1D hNNegTracks;
- ICloud1D hNPositronsForTrack;
- ICloud1D hNElectronsForTrack;
- ICloud1D hNPositronsForTrackInv;
- ICloud1D hNElectronsForTrackInv;
+ IHistogram1D hNPosTracks;
+ IHistogram1D hNNegTracks;
+ IHistogram1D hNPositronsForTrack;
+ IHistogram1D hNElectronsForTrack;
+ IHistogram1D hNPositronsForTrackInv;
+ IHistogram1D hNElectronsForTrackInv;
+
+ IHistogram1D hTrackChi2;
HashMap<Integer,MCParticle> mc_ele_prev = new HashMap<Integer, MCParticle>();
HashMap<Integer,MCParticle> mc_pos_prev = new HashMap<Integer, MCParticle>();
@@ -93,8 +102,12 @@
public void detectorChanged(Detector detector) {
-
-
+ bfield = detector.getFieldMap().getField(new BasicHep3Vector(0., 0., 1.)).y();
+ System.out.println("b-field = " + bfield);
+
+
+ bfac = Constants.fieldConversion;
+
//pFrame = new AIDAFrame();
//pFrame.setTitle("Truth p Plots");
makePlots();
@@ -107,9 +120,15 @@
}
- public double getMomentum(Track track) {
- double[] p_vec = track.getTrackStates().get(0).getMomentum();
- return Math.sqrt(p_vec[0]*p_vec[0] + p_vec[1]*p_vec[1] +p_vec[2]*p_vec[2]);
+ public double getPt(Track trk, double b) {
+ double R = 1.0/trk.getTrackStates().get(0).getOmega();
+ System.out.println("Omega " + trk.getTrackStates().get(0).getOmega() + " R " + R + " b " + b + " bfac " + Constants.fieldConversion);
+ return Constants.fieldConversion * b * Math.abs(R);
+ }
+
+ public double getP(Track trk, double b) {
+ double theta = Math.atan(1.0/trk.getTrackStates().get(0).getTanLambda());
+ return getPt(trk,b) / Math.sin(theta);
}
@@ -156,8 +175,11 @@
int[] ntrks = {0,0,0};
for (Track trk : tracklist) {
- double p = this.getMomentum(trk);
+ //Hep3Vector p_vec = new BasicHep3Vector(trk.getTrackStates().get(0).getMomentum());
+ double p = getP(trk,bfield);
+ System.out.printf("p %f pT %f\n", p, getPt(trk,bfield));
this.hTrackP.fill(p);
+ this.hTrackChi2.fill(trk.getChi2());
if(this.isElectronTrack(trk)) {
this.hNegTrackP.fill(p);
if(electron!=null) {
@@ -283,7 +305,15 @@
private void makePlots() {
+ pPlotterProp = af.createPlotterFactory().create("Track recon plots");
+ pPlotterProp.show();
+ pPlotterProp.setTitle("Track recon plots");
+ pPlotterProp.createRegions(1, 2);
+ hTrackChi2 = aida.histogram1D("Track chi2", 30, 0,10);
+ pPlotterProp.region(0).plot(hTrackChi2);
+
pPlotter = af.createPlotterFactory().create("Truth p Plots");
+ pPlotter.show();
pPlotter.setTitle("Truth p Plots");
//pFrame.addPlotter(pPlotter);
IPlotterStyle style0 = pPlotter.style();
@@ -307,17 +337,19 @@
IDataPointSetFactory dpsf = aida.analysisFactory().createDataPointSetFactory(aida.tree());
hTruthMatchedPosTrackPres = dpsf.create("hTruthMatchedPosTrackPres", "RMS(Track p - p(e+)) q>0 vs P",2);
hTruthMatchedNegTrackPres = dpsf.create("hTruthMatchedNegTrackPres", "RMS(Track p - p(e-)) q<0 vs P",2);
- hNTracks = aida.cloud1D("Ntrks");
+ hNTracks = aida.histogram1D("Ntrks",6, -0.5, 5.5);
+
+
trkCountVsEventPlot = aida.histogram1D("Number of Tracks vs Event Nr", 501, -0.5, 500.5);
- trkCountVsEventPlot.annotation().addItem("xAxisLabel", "Event Number");
- hNPosTracks = aida.cloud1D("Ntrks q>0");
- hNNegTracks = aida.cloud1D("Ntrks q<0");
+ trkCountVsEventPlot.annotation().addItem("xAxisLabel", "Event Number");
+ hNPosTracks = aida.histogram1D("Ntrks q>0",6, -0.5, 5.5);
+ hNNegTracks = aida.histogram1D("Ntrks q<0",6, -0.5, 5.5);
hElectronP = aida.histogram1D("Electron Momentum", 50, 0,4);
hPositronP = aida.histogram1D("Positron Momentum", 50, 0,4);
- hNPositronsForTrack = aida.cloud1D("N positrons given track with q>0");
- hNElectronsForTrack = aida.cloud1D("N electrons given track with q<0");
- hNPositronsForTrackInv = aida.cloud1D("N positrons given track with q<0");
- hNElectronsForTrackInv = aida.cloud1D("N electrons given track with q>0");
+ hNPositronsForTrack = aida.histogram1D("N positrons given track with q>0",6, -0.5, 5.5);
+ hNElectronsForTrack = aida.histogram1D("N electrons given track with q<0",6, -0.5, 5.5);
+ hNPositronsForTrackInv = aida.histogram1D("N positrons given track with q<0",6, -0.5, 5.5);
+ hNElectronsForTrackInv = aida.histogram1D("N electrons given track with q>0",6, -0.5, 5.5);
pPlotter.region(0).plot(hTrackP);
@@ -371,8 +403,8 @@
pPlotter3.createRegions(2, 5);
for(int i=0;i<5;++i) {
- hTruthMatchedPosTrackPdiffPrev[i] = aida.cloud1D("Track p - p(e+) q>0 BC=-" + i);
- hTruthMatchedNegTrackPdiffPrev[i] = aida.cloud1D("Track p - p(e-) q<0 BC=-" + i);
+ hTruthMatchedPosTrackPdiffPrev[i] = aida.histogram1D("Track p - p(e+) q>0 BC=-" + i, 100, -0.2,0.2);
+ hTruthMatchedNegTrackPdiffPrev[i] = aida.histogram1D("Track p - p(e-) q<0 BC=-" + i, 100, -0.2,0.2);
pPlotter3.region(i).plot(hTruthMatchedPosTrackPdiffPrev[i]);
pPlotter3.region(5+i).plot(hTruthMatchedNegTrackPdiffPrev[i]);
}
