lcsim/sandbox/NickSinev/Examples
diff -N FastMC_tgen.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ FastMC_tgen.java 8 Jun 2010 18:01:22 -0000 1.1
@@ -0,0 +1,291 @@
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.*;
+import org.lcsim.event.util.*;
+import org.lcsim.event.base.*;
+import org.lcsim.geometry.Detector;
+import org.lcsim.contrib.NickSinev.tracking.wmfitter.*;
+import org.lcsim.contrib.NickSinev.tracking.util.*;
+import org.lcsim.util.Driver;
+import org.lcsim.contrib.NickSinev.ztracking.*;
+import org.lcsim.contrib.NickSinev.event.util.*;
+import org.lcsim.constants.Constants;
+import hep.physics.vec.Hep3Vector;
+import hep.physics.vec.VecOp;
+import hep.physics.vec.BasicHep3Vector;
+import hep.physics.matrix.*;
+import org.lcsim.fit.circle.*;
+import java.text.*;
+import java.util.*;
+import org.lcsim.util.aida.AIDA;
+import Jama.*;
+import java.io.*;
+import org.lcsim.util.cache.*;
+
+public class FastMC_tgen extends Driver
+{
+ private boolean masklayer = false;
+ private int maskedlayer = 24;
+ private String ofnam="Bruce_finebns_nomskl_nbc.txt";
+
+ int firstev = 0;
+ int lastev=1;
+ double[] curMom = new double[51];
+ double mstep1 = 1.0838;
+ double mstep2 = 1.2589254;
+ double mmin = 0.2;
+ double[] curCT = {.000, .109, .206, .292, .369, .438, .499, .553, .602, .645,
+ .684, .718, .749, .776, .800, .822, .842, .859, .874, .888,
+ .900, .911, .921, .929, .937, .944, .950, .955, .960, .965,
+ .968, .972, .975, .978, .980, .982, .984, .986, .987, .989, .990 };
+ double[][][] table;
+
+ double curlam = 0.05;
+ boolean real = true;
+ double minPt = 0.1;
+ double maxPt =500.;
+ double curTl = 0.;
+ private boolean incair = false;
+ private boolean doeloss = true;
+ private boolean inccor = true; //*****************CORRELATIONS*********
+ private boolean IPconstrained = false;
+ private boolean incBP = true;
+ private boolean planarfit = false;
+ private boolean equalwts = false;
+ private double vtxbrphires = 0.0035;
+ private double vtxbzres = 0.0035;
+ private double vtxecrpres = 0.0035;
+ private double vtxecrres = 0.0035;
+ private double trkbrpres = 0.007;
+ private double trkbzres = 100.;
+ private double trkecrpres = 0.007;
+ private double trkecrres = 0.035;
+ private double trkfrpres = 0.007;
+ private double trkfrres = 0.007;
+ private SLDFitTrack fitter = new SLDFitTrack();
+ private CircleFitter ciftr = new CircleFitter();
+ private DecimalFormat df = new DecimalFormat(" ###0.000");
+ private DecimalFormat dfd = new DecimalFormat(".0000E00;-.0000E00");
+ private SmearTrackerHits smTh = SmearTrackerHits.getInstance();
+ private TrackingCheaterSmHits cheat = new TrackingCheaterSmHits();
+ private WeightMatrix Wm = null;
+ private WMTrackPropagator trkPr= null;
+ private int evnum = 0;
+ double rsmax = 1500.;
+ double dphi = 0.;
+ int nall = 0;
+ int ttn = 0;
+ int ngdft = 0;
+ int nalspe = 0;
+ double field = 5.*Constants.fieldConversion;
+ String dnam = " ";
+ int nctbns = 0;
+ int nmombns=0;
+ double[][] conversionFromMmToCm = { { 0.01, 0.1, 1., 0.01, 0.1},
+ { 0.1, 1.00, 10., 0.1, 1. },
+ { 1., 10.0, 100., 1., 10.},
+ { 0.01, 0.1, 1., 0.01, 0.1},
+ { 0.1, 1.0, 10., 0.1, 1. } };
+ int[][] pind = { { 0, 1, 3, 6, 10},
+ { 1, 2, 4, 7, 11},
+ { 3, 4, 5, 8, 12},
+ { 6, 7, 8, 9, 13},
+ {10, 11, 12, 13, 14}};
+
+ public FastMC_tgen()
+ {
+ curMom[0]=mmin;
+ for(int i=0; i<50; i++)
+ {
+ if(i<20) curMom[i+1]=curMom[i]*mstep1;
+ else curMom[i+1]=curMom[i]*mstep2;
+ }
+ nctbns = curCT.length;
+ nmombns = curMom.length;
+ table = new double[15][nctbns][nmombns];
+ smTh.setVtxBarrRPResolution(vtxbrphires);
+ smTh.setVtxBarrZResolution(vtxbzres);
+ smTh.setVtxECRPResolution(vtxecrpres);
+ smTh.setVtxECRResolution(vtxecrres);
+ smTh.setTrkBarrRPResolution(trkbrpres);
+ smTh.setTrkBarrZResolution(trkbzres);
+ smTh.setTrkECRPResolution(trkecrpres);
+ smTh.setTrkECRResolution(trkecrres);
+ cheat.useChargedOnly = true;
+ smTh.setDoSmear(true);
+ Wm = fitter.getWeightMatrixObject();
+ add(smTh);
+ add(cheat);
+ for(int k=0; k<nctbns; k++)
+ {
+ for(int m=0; m<nmombns; m++)
+ {
+ for(int i=0; i<5; i++)
+ {
+ for(int j=0; j<=i; j++)
+ {
+ int ind = pind[i][j];
+ if(i==j) table[ind][k][m] = 1.;
+ else table[ind][k][m] = 0.;
+ }
+ }
+ }
+ }
+ }
+
+ public void process(EventHeader event)
+ { //1
+ if((evnum > firstev-1)&&(evnum<lastev))
+ { //2
+ super.process(event);
+ if(evnum == 0)
+ { //3
+ Detector det = event.getDetector();
+ dnam = det.getName();
+ trkPr = fitter.getTrackPropagator();
+ trkPr.setIncludeBeamPipe(incBP);
+ fitter.setPlanarFit(planarfit);
+ fitter.setConstrainedFit(IPconstrained);
+ fitter.setIncludeAir(incair);
+ fitter.setDoELoss(doeloss);
+ fitter.setIncludeCorrelations(inccor);
+ fitter.setDetector(det);
+ fitter.setDoHist(false);
+ Wm.setEqualWeights(equalwts);
+ Wm.setDebugLevel(0);
+ if(masklayer) Wm.maskLayer(maskedlayer);
+ Wm.setDoHist(false);
+ } //2
+ List<CheatedBaseTrack> tracks = null;
+ try
+ { //3
+ tracks = event.get(CheatedBaseTrack.class,"CheatedBaseTracks");
+ } //2
+ catch(Exception trnotf)
+ { //3
+ System.out.println("No CheatedBaseTracks in the event ! \n");
+ } //2
+ finally
+ { //3
+ } //2
+ if(tracks != null)
+ {//3
+ int nprts = tracks.size();
+ System.out.println("Tracking cheater has found "+nprts+"BaseTracks in event "+evnum);
+ int tn = 0;
+ double[] ortpar=new double[5];
+ double[] mcpar=new double[5];
+ for (CheatedBaseTrack track : tracks )
+ {//4
+ double[] opar = track.getTrackParameters();
+ for(int i=0; i<5; i++) mcpar[i]=opar[i];
+ MCParticle mcp = track.getMCParticle();
+ if(mcp != null)
+ {//5
+ Hep3Vector pm = mcp.getMomentum();
+ double mPt = Math.sqrt(pm.x()*pm.x()+pm.y()*pm.y());
+ for(int m=0; m<nmombns; m++)
+ {//6
+ System.out.println("Simulating tracks with total momentum: "+curMom[m]);
+ for(int k=0; k<nctbns; k++)
+ {//7
+ System.out.print(" .");
+ for(int i=0; i<5; i++)
+ {//8
+ ortpar[i]=mcpar[i];
+ }//7
+ double curPt = curMom[m]*Math.sqrt(1.-curCT[k]*curCT[k]);
+ curTl=curCT[k]/Math.sqrt(1.-curCT[k]*curCT[k])+0.01;
+ double omsc = mPt/curPt;
+ ortpar[2]=Math.abs(mcpar[2])*omsc;
+ ortpar[4]=curTl;
+ track.setTrackParameters(ortpar,5.);
+ fitter.setTrack(track);
+ fitter.simulateFit(mcp);
+ int ndf = track.getNDF();
+ boolean badFit = false;
+ if(track.fitSuccess())
+ {//8
+ SymmetricMatrix trcovm = track.getErrorMatrix();
+ for(int i=0; i<5; i++)
+ {//9
+ for(int j=0; j<=i; j++)
+ {//10
+ double e = trcovm.e(j,i);
+ e*=conversionFromMmToCm[j][i];
+ trcovm.setElement(j,i,e);
+ int ind = pind[i][j];
+ table[ind][k][m]=e;
+ } //9
+ } //8
+ } //7
+ } //6
+ System.out.println();
+ } //5
+ } //4
+ tn++;
+ } //3
+ } //2
+ } //1
+ evnum++;
+ if(evnum == lastev) writeTable(ofnam);
+ } //0
+
+ private void writeTable(String fnam)
+ {
+ String cacheDir = System.getProperty("user.home");
+ File cachedir = new File(cacheDir);
+ if(cachedir == null)
+ {
+ System.out.println("CovMatrix is unable to find user home directory!");
+ }
+ else
+ {
+ File home = new File(cachedir,".cache");
+ File cache = new File(home,fnam);
+ try
+ {
+ cache.createNewFile();
+ FileWriter ofw = new FileWriter(cache);
+ int ind = 0;
+ for(int i=0; i<5; i++)
+ {
+ for(int j=0; j<=i; j++)
+ {
+ ind = pind[i][j];
+ int inf = i+1;
+ int jnf = j+1;
+ ofw.write("Cov matrix entry "+inf+" , "+jnf+" \r\n \r\n");
+ ofw.write(" "+nctbns+" \r\n");
+ ofw.write(" "+nmombns+" \r\n");
+ for(int k=0; k<nctbns; k++)
+ {
+ ofw.write(" "+df.format(curCT[k])+" \r\n");
+ for(int m=0; m<nmombns; m++)
+ {
+ double val=curMom[m];
+ if(val<10.) ofw.write(" "+df.format(val)+" "+dfd.format(table[ind][k][m])+" \r\n");
+ if((val<100.) && (val>=10.)) ofw.write(" "+df.format(val)+" "+dfd.format(table[ind][k][m])+" \r\n");
+ if((val<1000.) && (val>=100.)) ofw.write(" "+df.format(val)+" "+dfd.format(table[ind][k][m])+" \r\n");
+ if(val>=1000.) ofw.write(" "+df.format(val)+" "+dfd.format(table[ind][k][m])+" \r\n");
+ }
+ if(k != (nctbns-1))ofw.write("\r\n");
+ }
+ ofw.write(" end\r\n\r\n");
+ }
+ }
+ ofw.write("\r\n");
+ ofw.flush();
+ ofw.close();
+ System.out.println("Output file has been written\n");
+ }
+ catch(IOException e)
+ {
+ System.out.println("IOException caught: "+e.getMessage());
+ }
+ finally
+ {
+ }
+ }
+ }
+
+}
\ No newline at end of file
