hps-java/src/main/java/org/lcsim/hps/users/phansson
diff -u -r1.1 -r1.2
--- FastTrackResidualDriver.java 4 May 2012 15:45:44 -0000 1.1
+++ FastTrackResidualDriver.java 4 May 2012 18:23:35 -0000 1.2
@@ -5,6 +5,7 @@
package org.lcsim.hps.users.phansson;
import hep.aida.*;
+import hep.aida.ref.plotter.PlotterRegion;
import hep.physics.vec.BasicHep3Vector;
import hep.physics.vec.Hep3Vector;
import java.io.IOException;
@@ -29,20 +30,21 @@
public class FastTrackResidualDriver extends Driver {
int nevents = 0;
- private boolean debug = false;
- private boolean doFastTracking = false;
- private boolean doDefaultResiduals = false;
+ boolean debug = false;
+ boolean saveFile = false;
protected IDDecoder dec = null;
protected Subdetector ecal;
private String ecalName = "Ecal";
+ String sides[] = {"up","down"};
double crystalX;
double crystalY;
double beamGap;
- private double EcalZPosition;
+ double EcalZPosition;
+ double conversionZ;
private static int crystalCols;
private static int crystalRows;
- private String outputPlotFileName;
+ private String outputPlotFileName = "";
private AIDA aida = AIDA.defaultInstance();
private IAnalysisFactory af = aida.analysisFactory();
@@ -52,6 +54,7 @@
private List< IHistogram1D > nhits_tracker = new ArrayList<IHistogram1D>();
private List< IHistogram1D > ncl_ecal = new ArrayList<IHistogram1D>();
private List< IHistogram2D > ncl_ecal_map = new ArrayList<IHistogram2D>();
+ private List< IHistogram2D > nselcl_ecal_map = new ArrayList<IHistogram2D>();
public void startOfData() {
@@ -72,14 +75,9 @@
crystalCols = 46;
crystalRows = 5;
- displayFastTrackingPlots();
+ // Position of the conversion
+ conversionZ = -1000.0;
-
-
- }
-
- public FastTrackResidualDriver() {
-
IHistogramFactory hf = aida.histogramFactory();
String side;
for (int iSide=0;iSide<2;++iSide) {
@@ -92,17 +90,19 @@
List<IHistogram1D> listLay1Hit = new ArrayList<IHistogram1D>();
resy_org_lay1hit.add(listLay1Hit);
- nhits_tracker.add(hf.createHistogram1D("FT_nhits_tracker_" + side , 20, 0, 20));
- ncl_ecal.add(hf.createHistogram1D("FT_ncl_ecal_" + side , 20, 0, 20));
+ nhits_tracker.add(hf.createHistogram1D("FT_nhits_tracker_" + side , 40, 0, 40));
+ ncl_ecal.add(hf.createHistogram1D("FT_ncl_ecal_" + side , 60, 0, 60));
//Setup the ecal 2D plot
- //IHistogram2D tmp = hf.createHistogram2D("htmp", 10, 0.0, 10., 10, 0.0, 10.0);
- if (side == "up") ncl_ecal_map.add(hf.createHistogram2D("FT_ecal_crystalclusterhitmap_" + side, 92, -46.0, 46.0, 5, 0.0, 5.0));
- else ncl_ecal_map.add(hf.createHistogram2D("FT_ecal_crystalclusterhitmap_" + side, 92, -46.0, 46.0, 5, -5.0, 0.0));
- //hf.createHistogr
- //hf.createHistogram2D(side, iSide, beamGap, beamGap, iSide, beamGap, beamGap)
- //hf.createHistogram2D
+ if (side == "up") {
+ ncl_ecal_map.add(hf.createHistogram2D("FT_ecal_hitmap_" + side, 92, -46.0, 46.0, 5, 0.0, 5.0));
+ nselcl_ecal_map.add(hf.createHistogram2D("FT_sel_ecal_hitmap_" + side, 92, -46.0, 46.0, 5, 0.0, 5.0));
+ } else {
+ ncl_ecal_map.add(hf.createHistogram2D("FT_ecal_hitmap_" + side, 92, -46.0, 46.0, 5, -5.0, 0.0));
+ nselcl_ecal_map.add(hf.createHistogram2D("FT_sel_ecal_hitmap_" + side, 92, -46.0, 46.0, 5, -5.0, 0.0));
+ }
+
for (int iLayer=1;iLayer<6;++iLayer) {
IHistogram1D h = hf.createHistogram1D("FT_res_" + side + "_l"+iLayer, 50, -25, 25);
@@ -117,6 +117,30 @@
}
}
+
+ displayFastTrackingPlots();
+
+
+
+ }
+
+ public FastTrackResidualDriver() {}
+
+
+ public void setConversionZ(double z) {
+ this.conversionZ = z;
+ }
+
+ public void setDebug(boolean flag) {
+ this.debug = flag;
+ }
+
+ public void setOutputPlotFileName( String name ) {
+ this.outputPlotFileName = name;
+ }
+
+ public void setEcalZPosition(double val) {
+ this.EcalZPosition = val;
}
public void process(EventHeader event) {
@@ -149,7 +173,7 @@
//if ( 1==1 ) return;
- String sides[] = {"up","down"};
+
int nhits;
int nhitsInTracker;
int nhitsInLayer1;
@@ -164,25 +188,25 @@
int layer;
String si_side;
double res;
+ int layerIndex = -1;
for (int iSide=0;iSide<2;++iSide) {
- if(debug) System.out.println("Side: " + sides[iSide]);
-
+
//if ( 1==1 ) return;
ecal_cls = getEcalClustersForFastTracking(ecal_all_clusters, sides[iSide]);
- if (debug) System.out.println("This side has " + ecal_cls.size() + " clusters");
+
+ if (debug) System.out.println("Found " + ecal_cls.size() +" Ecal clusters on the " + sides[iSide]);
//if ( 1==1 ) return;
ncl_ecal.get(iSide).fill(ecal_cls.size());
if( ecal_cls.size() ==0 ) {
- System.out.println("No clusters...");
+ System.out.println("No clusters on this side...");
continue;
}
- //if ( 1==1 ) return;
//Fill map of Ecal hits
for( Hep3Vector cl : ecal_cls) {
@@ -191,25 +215,38 @@
}
ecal_cl = selectCluster(ecal_cls);
-
+
+ if (debug) System.out.println("Selected " + ecal_cl.toString() +" cluster to use as pointer of fast track ");
+
+
+ int selclpos[] = getCrystalPair(ecal_cl);
+ nselcl_ecal_map.get(iSide).fill(selclpos[0], selclpos[1]);
+
//Get "target" position i.e. the origin of the radiation
origin_pos = getFastTrackOrigin(event);
+ if (debug) System.out.println("Conversion started at " + origin_pos.toString());
+
//Create the fast track
fastTrack = new FastTrack(origin_pos,ecal_cl,debug);
if ( debug ) System.out.println(fastTrack.toString());
- if ( 1==1 ) return;
+
nhits = 0;
nhitsInTracker = getNLayersWithAxialHit(trackerHits,sides[iSide]);
+
+ if( debug ) System.out.println("There are " + nhitsInTracker + " hits on this side of the tracker");
+
//most upstream layer nr is different for top and bottom
int firstLayer = 1;
if (sides[iSide] == "down") firstLayer=2;
nhitsInLayer1 = getNAxialHitsInLayers(trackerHits,sides[iSide],firstLayer);
-
+
+ if( debug ) System.out.println("There are " + nhitsInLayer1 + " hits on the first layer for this side (layer=" + firstLayer +")");
+
nhits_tracker.get(iSide).fill(nhitsInTracker);
@@ -223,6 +260,7 @@
siSensor = stripCluster.getSensor();
name = siSensor.getName();
+
if ( name.length() < 14) {
System.err.println("This name is too short!!");
throw new RuntimeException("SiSensor name " + name + " is invalid?");
@@ -234,18 +272,17 @@
si_side = getSideFromSiCluster(stripCluster);
if ( debug ) {
- System.out.println("hit " + nhits + " on " + si_side + " side of detector");
- System.out.println("SiSensor layer " + layer + "(" + name + ")");
- System.out.println("isAxial? " + isaxial);
+ System.out.println("hit " + nhits + " in layer " + layer + " on " + si_side + " side of detector");
}
- if( sides[iSide] == si_side) {
- System.out.println("This hit is same side as Ecal cluster side");
- } else {
- System.out.println("This hit is opposite side side as Ecal cluster side -> skip!");
+ if( sides[iSide] != si_side) {
+ if (debug) System.out.println("This hit is opposite side side as Ecal cluster side -> skip!");
continue;
}
+
+ if (debug) System.out.println("Hit and cluster both on " + sides[iSide] + " side -> use it!");
+
res = fastTrack.getFastTrackResidual(stripCluster);
@@ -253,9 +290,9 @@
// 0 1 2 3 4
// 0 1 2 3 4
- // Fix the layer for the list index
- int layerIndex = layer;
- if ( layer % 2 ==0 ) layerIndex=layer/2-1;
+ // Fix the layer for the list index -> ugly
+ layerIndex = layer;
+ if ( layer % 2 == 0 ) layerIndex=layer/2-1;
else layerIndex = (layer-1)/2;
resy_org.get(iSide).get(layerIndex).fill(res);
@@ -288,17 +325,13 @@
//First a workaround - use hits from track..
//List<TrackerHit> hits = getAllHitsFromTracks(event);
List<SiTrackerHitStrip1D> hits = getAllHits(event);
- if( debug ) {
- System.err.println("Found " + hits.size() + " hits");
- }
return hits;
}
private Hep3Vector getFastTrackOrigin(EventHeader event) {
double tx = 0.0;
double ty = 0.0;
- double tz = 0.0;
- Hep3Vector pos = new BasicHep3Vector(tx,ty,tz);
+ Hep3Vector pos = new BasicHep3Vector(tx,ty,conversionZ);
return pos;
}
@@ -314,18 +347,7 @@
// return clusters;
// }
- private double[] getClusterPosition(HPSEcalCluster cluster) {
- CalorimeterHit hit = cluster.getSeedHit();
-
- //IDDecoder dec = dec.getSubdetector("Ecal").getIDDecoder();
- dec.setID(hit.getCellID());
- int ix = dec.getValue("ix");
- int iy = dec.getValue("iy");
- double position[] = new double[2];
- position[0] = crystalX * ix;
- position[1] = crystalY * iy + beamGap * Math.signum(iy);
- return position;
- }
+
private int[] getCrystalPair(Hep3Vector cluster) {
double x = cluster.x();
@@ -335,7 +357,7 @@
//position[0] = crystalX * ix;
position[1] = (int) Math.floor( (y - beamGap * Math.signum(y) ) / crystalY);
//position[1] = crystalY * iy + beamGap * Math.signum(iy);
- System.out.println("x " + x + " y " + y + " --> ix " + position[0] + " iy " + position[1]);
+ //System.out.println("x " + x + " y " + y + " --> ix " + position[0] + " iy " + position[1]);
return position;
//IDDecoder dec = dec.getSubdetector("Ecal").getIDDecoder();
@@ -345,7 +367,19 @@
}
-
+ private double[] getClusterPosition(HPSEcalCluster cluster) {
+ CalorimeterHit hit = cluster.getSeedHit();
+ //return cluster.getPosition();
+
+ //IDDecoder dec = dec.getSubdetector("Ecal").getIDDecoder();
+ dec.setID(hit.getCellID());
+ int ix = dec.getValue("ix");
+ int iy = dec.getValue("iy");
+ double position[] = new double[2];
+ position[0] = crystalX * ix;
+ position[1] = crystalY * iy + beamGap * Math.signum(iy);
+ return position;
+ }
private List<Hep3Vector> getEcalClustersForFastTracking(List<HPSEcalCluster> clusters, String side) {
if(!side.equalsIgnoreCase("up") && !side.equalsIgnoreCase("down")) {
@@ -353,25 +387,21 @@
}
List<Hep3Vector> cls = new ArrayList<Hep3Vector>();
boolean save;
- double [] pos_xy;
- //if ( 1==1 ) return cls;
+ double [] pos;
for ( HPSEcalCluster cl : clusters) {
- //for ( : clusters) {
-
- //if (1==1) return cls;
save=false;
- if (1==1) return cls;
- pos_xy = getClusterPosition(cl);
- if (1==1) return cls;
+ pos = getClusterPosition(cl);
- if(pos_xy[1]>=0 && side.equalsIgnoreCase("up")) {
+ if(pos[1]>=0 && side.equalsIgnoreCase("up")) {
save=true;
- } else if(pos_xy[1]<0 && side.equalsIgnoreCase("down")) {
+ } else if(pos[1]<0 && side.equalsIgnoreCase("down")) {
save=true;
}
- if(save) cls.add(new BasicHep3Vector(pos_xy[0],pos_xy[1],EcalZPosition));
- //if (1==1) return cls;
+ if(save) {
+ cls.add(new BasicHep3Vector(pos[0],pos[1],EcalZPosition));
+ //cls.add(new BasicHep3Vector(pos[0],pos[1],pos[2]));
+ }
}
return cls;
@@ -391,10 +421,11 @@
private List<HPSEcalCluster> getAllEcalClustersForFastTracking(EventHeader event) {
//List<HPSEcalCluster> clusters = new ArrayList<HPSEcalCluster>();
- //List<HPSEcalCluster> clusters = event.get(HPSEcalCluster.class, "EcalClusters");
- List<HPSEcalCluster> clusters = event.get(HPSEcalCluster.class, "EcalReadoutHits");
+ List<HPSEcalCluster> clusters = event.get(HPSEcalCluster.class, "EcalClusters");
+ //List<HPSEcalCluster> clusters = event.get(HPSEcalCluster.class, "EcalReadoutHits");
+ //List<HPSEcalCluster> clusters = event.get(HPSEcalCluster.class, "EcalCalHits");
if ( debug) {
- System.out.println("Found " + clusters.size() + " clusters");
+ System.out.println("Found " + clusters.size() + " EcalClusters");
}
@@ -406,10 +437,10 @@
private boolean isAxialHit(SiTrackerHitStrip1D strip) {
boolean axial = false;
Hep3Vector m = strip.getMeasuredCoordinate();
- System.out.println("y " + m.y() + " -> abs(y) " + Math.abs(m.y()) );
+ //System.out.println("y " + m.y() + " -> abs(y) " + Math.abs(m.y()) );
if ( Math.abs((Math.abs(m.y())-1.0))<0.0000000001 ) {
axial = true;
- System.out.println( " ===> " + axial );
+ //System.out.println( " ===> " + axial );
}
return axial;
}
@@ -476,13 +507,10 @@
String l = name.substring(13);
l = l.substring(0, l.indexOf("_"));
- System.out.println("layer: " + l + " (length: " + l.length() + ")");
- //String ll = l.subs
- //System.out.println("ll: " + ll);
ilayer = Integer.parseInt(l);
- System.out.println("il: " + ilayer);
+
if ( ilayer < 1 || ilayer > 10 ) {
- System.err.println("This layer doesn't exist?");
+ System.out.println("This layer doesn't exist?");
throw new RuntimeException("SiSensor name " + name + " is invalid?");
}
@@ -490,23 +518,21 @@
throw new RuntimeException("This sensor name do not have a layer!");
}
-
return ilayer;
}
-
-
-
+
+
public void displayFastTrackingPlots() {
//IPlotter plotter = af.createPlotterFactory(af.createTreeFactory().create());
- IPlotter plotter_org = af.createPlotterFactory().create();
+ IPlotter plotter_org = af.createPlotterFactory().create("HPS SVT Fast Track Residuals (All tracks)");
plotter_org.createRegions(5,2,0);
- IPlotter plotter_org_LayAllHit = af.createPlotterFactory().create();
+ IPlotter plotter_org_LayAllHit = af.createPlotterFactory().create("HPS SVT Fast Track Residuals (All layers has hits)");
plotter_org_LayAllHit.createRegions(5,2,0);
- IPlotter plotter_org_Lay1Hit = af.createPlotterFactory().create();
+ IPlotter plotter_org_Lay1Hit = af.createPlotterFactory().create("HPS SVT Fast Track Residuals (1st layer has hits)");
plotter_org_Lay1Hit.createRegions(5,2,0);
for ( int iSide=0;iSide<2;++iSide) {
@@ -520,29 +546,42 @@
plotter_org.region(idx).plot(resy_org.get(iSide).get(iLayer));
plotter_org_LayAllHit.region(idx).plot(resy_org_layallhit.get(iSide).get(iLayer));
plotter_org_Lay1Hit.region(idx).plot(resy_org_lay1hit.get(iSide).get(iLayer));
+
}
}
plotter_org.show();
plotter_org_LayAllHit.show();
plotter_org_Lay1Hit.show();
+ //Hit mupltiplicity
IPlotter plotter_hitmult = af.createPlotterFactory().create();
plotter_hitmult.createRegions(1,2,0);
plotter_hitmult.region(0).plot(nhits_tracker.get(0));
plotter_hitmult.region(1).plot(nhits_tracker.get(1));
-
plotter_hitmult.show();
+ //ECal plots
IPlotter plotter_ecalhitmult = af.createPlotterFactory().create();
- plotter_ecalhitmult.createRegions(2,2,0);
- plotter_ecalhitmult.region(0).plot(ncl_ecal_map.get(0));
- plotter_ecalhitmult.region(1).plot(ncl_ecal_map.get(1));
- plotter_ecalhitmult.region(2).plot(ncl_ecal.get(0));
- plotter_ecalhitmult.region(3).plot(ncl_ecal.get(1));
+ plotter_ecalhitmult.createRegions(3,2,0);
+ plotter_ecalhitmult.region(0).plot(ncl_ecal.get(0));
+ plotter_ecalhitmult.region(1).plot(ncl_ecal.get(1));
+ plotter_ecalhitmult.region(2).plot(ncl_ecal_map.get(0));
+ plotter_ecalhitmult.region(3).plot(ncl_ecal_map.get(1));
+ plotter_ecalhitmult.region(4).plot(nselcl_ecal_map.get(0));
+ plotter_ecalhitmult.region(5).plot(nselcl_ecal_map.get(1));
plotter_ecalhitmult.show();
+ for (int idx=2;idx<6;++idx) {
+ IPlotterStyle style = plotter_ecalhitmult.region(idx).style();
+ style.setParameter("hist2DStyle", "colorMap");
+ style.dataStyle().fillStyle().setParameter("colorMapScheme", "rainbow");
+ ((PlotterRegion) plotter_ecalhitmult.region(idx)).getPlot().setAllowUserInteraction(false);
+ ((PlotterRegion) plotter_ecalhitmult.region(idx)).getPlot().setAllowPopupMenus(false);
+ }
+
+
}
@@ -550,25 +589,16 @@
public void endOfData() {
- //try {
- // aida.saveAs(outputPlotFileName);
- //} catch (IOException ex) {
- // Logger.getLogger(FastTrackResidualDriver.class.getName()).log(Level.SEVERE, "Couldn't save aida plots to file " + outputPlotFileName, ex);
- //}
- //displayFastTrackingPlots();
+
+ if (outputPlotFileName != "")
+ try {
+ aida.saveAs(outputPlotFileName);
+ } catch (IOException ex) {
+ Logger.getLogger(FastTrackResidualDriver.class.getName()).log(Level.SEVERE, "Couldn't save aida plots to file " + outputPlotFileName, ex);
+ }
+ displayFastTrackingPlots();
}
- public void setDebug(boolean debug) {
- this.debug = true;
- }
-
- public void setOutputPlotFileName( String name ) {
- this.outputPlotFileName = name;
- }
-
- public void setEcalZPosition(double val) {
- this.EcalZPosition = val;
- }
}