lcsim/sandbox/RobKutschke/TRFTests/v1
diff -u -r1.1 -r1.2
--- FitTest_sid00Driver.java 30 Aug 2007 14:16:33 -0000 1.1
+++ FitTest_sid00Driver.java 7 Sep 2007 19:04:15 -0000 1.2
@@ -29,9 +29,9 @@
* no multiple scattering or energy loss.
*
*@author $Author: kutschke $
- *@version $Id: FitTest_sid00Driver.java,v 1.1 2007/08/30 14:16:33 kutschke Exp $
+ *@version $Id: FitTest_sid00Driver.java,v 1.2 2007/09/07 19:04:15 kutschke Exp $
*
- * Date $Date: 2007/08/30 14:16:33 $
+ * Date $Date: 2007/09/07 19:04:15 $
*
*/
@@ -49,12 +49,15 @@
RKTrackGenParams param = null;
RKTrackGen gen = null;
- private TrackError vstartc = new TrackError();
- private TrackError vstartz = new TrackError();
+ private TrackError vstartc_in = new TrackError();
+ private TrackError vstartz_in = new TrackError();
+ private TrackError vstart_out = new TrackError();
private RKTrackFitDiag DiagForward = null;
private RKTrackFitDiag DiagBackward = null;
+ private RKCovSmear covsmear = null;
+
public FitTest_sid00Driver(){
}
@@ -66,7 +69,7 @@
// Limits for track generation parameters.
param = new RKTrackGenParams();
param.setbz(rkgeom.getBz());
- param.setptmin( 10.);
+ param.setptmin( 10.0);
//param.setczmin( 0.93);
//param.setczmax( 0.993);
param.setczmin( -0.993);
@@ -88,19 +91,34 @@
seed += 876633217;
RKHit.setSeed(seed);
- // Intial covariance matrix.
- double sig=10.0;
- vstartc.set(0,0,sig);
- vstartc.set(1,1,500.);
- vstartc.set(2,2,sig);
- vstartc.set(3,3,sig);
- vstartc.set(4,4,sig);
-
- vstartz.set(0,0,sig);
- vstartz.set(1,1,sig);
- vstartz.set(2,2,sig);
- vstartz.set(3,3,sig);
- vstartz.set(4,4,sig);
+ // Seed for smearing covariance matrices.
+ seed += 876633217;
+ covsmear = new RKCovSmear( seed );
+
+ // Intial covariance matrices.
+ double sig=0.0625;
+ double fac1=10.;
+ double fac2=100.;
+ vstartc_in.set(0,0,sig/fac1);
+ vstartc_in.set(1,1,sig/fac1);
+ vstartc_in.set(2,2,sig/fac2);
+ vstartc_in.set(3,3,sig/fac2);
+ vstartc_in.set(4,4,sig/fac2);
+
+ vstartz_in.set(0,0,sig);
+ vstartz_in.set(1,1,sig);
+ vstartz_in.set(2,2,sig);
+ vstartz_in.set(3,3,sig);
+ vstartz_in.set(4,4,sig);
+
+ sig=1.0;
+ vstart_out.set(0,0,sig);
+ vstart_out.set(1,1,sig);
+ vstart_out.set(2,2,sig);
+ vstart_out.set(3,3,sig);
+ vstart_out.set(4,4,sig);
+
+
}
@@ -112,8 +130,9 @@
+ event.getDetector().getName() );
//int ndo = 1000;
- int ndo = 1000;
+ int ndo = 10000;
for ( int i=0; i<ndo; ++i ){
+ RKDebug.Instance().setTrack(i);
DoIt(i);
}
@@ -122,7 +141,7 @@
private void DoIt( int itrk){
- System.out.println ("New track: " + itrk );
+ // System.out.println ("New track: " + itrk );
// Generate a track.
RKTrack rktrk = gen.newRKTrack();
@@ -130,6 +149,9 @@
// Monitor generated track info.
plotGenTrk(rktrk);
+ // Load track info into the debug-info container.
+ RKDebug.Instance().setRKTrack(rktrk);
+
// Generate hits on this track in my own format.
RKHitList hlist = GenerateOneArcHits( rktrk );
@@ -147,7 +169,8 @@
toTRF trftrk = new toTRF(rktrk);
VTrack vtdca = trftrk.atDcaZaxis();
- TrackError vstart = vstartz;
+ TrackError vstart = vstart_out;
+ RKDebug.Instance().setStartType(0);
// Prep for the outward fit.
ETrack et = new ETrack( vtdca.surface().newPureSurface(), vtdca.vector(), vstart );
@@ -160,6 +183,14 @@
// Do the outward fit.
//System.out.println ("Start forward fit.");
int fstatf = fitk.fitForward(ht);
+ if ( fstatf != 0 ){
+ System.out.println("Forwards status: " +
+ fstatf + " " +
+ hlist.nHits() + " " +
+ hlist.nCyl() + " " +
+ hlist.nZp() + " " +
+ rktrk.cz() );
+ }
// Final state reference parameters for the outward fit.
VTrack vtend = hlist.outerMostHit().getVTrack();
@@ -170,14 +201,32 @@
double[] limits = { 0.00005, 0.15, 0.0003, 0.001, 0.0003};
DiagForward = new RKTrackFitDiag ( "Forward", "Cyl", nbins, limits );
}
- DiagForward.fill( fstatf, ht, vtend, hlist.nDof(), vstart);
+ DiagForward.fill( fstatf, ht, vtend, hlist.nDof(), vstart, rktrk);
+ // Check code for smearing tracks with covariance matrix.
+ TrackError te = ht.newTrack().error();
+ TrackVector vv = covsmear.Smear( ht.newTrack().error(), ht.newTrack().vector() );
+ for ( int i=0; i<5; ++i ){
+ double xx = te.get(i,i);
+ double sig = Math.sqrt(Math.abs(xx));
+ double pull = (vv.get(i)-ht.newTrack().vector(i))/sig;
+ aida.histogram1D("Smear Check",100, -5., 5.).fill(pull);
+ }
+
+
if ( vtend.surface().pureType().equals(SurfCylinder.staticType()) ){
- vstart = vstartc;
+ vstart = vstartc_in;
+ RKDebug.Instance().setStartType(1);
+ } else{
+ vstart = vstartz_in;
+ RKDebug.Instance().setStartType(2);
}
// Prep for inward fit.
- ETrack etin = new ETrack( vtend.surface().newPureSurface(), vtend.vector(), vstart );
+ ETrack etin = new ETrack( vtend.surface().newPureSurface(),
+ ht.newTrack().vector(),
+ vstart );
+ //vtend.vector(), vstart );
// Not sure why this is needed.
if ( vtend.surface().type() == SurfZPlane.staticType() ){
@@ -204,6 +253,11 @@
hlist.nCyl() + " " +
hlist.nZp() + " " +
rktrk.cz() );
+ if ( fstatb < 100 ){
+ aida.histogram1D("Backward: cz for failed propagation",50, -1., 1. ).fill(rktrk.cz());
+ } else{
+ aida.histogram1D("Backward: cz for failed addhit",50, -1., 1. ).fill(rktrk.cz());
+ }
}
if ( fstatb != 0 ){
System.out.println ("Bad fit: "
@@ -222,7 +276,11 @@
h.getPsi()
);
}
-
+ aida.cloud2D("Bad fit nZ vs nCyl").fill( hlist.nCyl(), hlist.nZp() );
+ aida.histogram1D( "Bad fit: Start type",5,0.,5.).fill(RKDebug.Instance().getStartType());
+ } else{
+ aida.cloud2D("Good fit nZ vs nCyl").fill( hlist.nCyl(), hlist.nZp() );
+ aida.histogram1D( "Good fit: Start type",5,0.,5.).fill(RKDebug.Instance().getStartType());
}
// Diagnostics for inward fit.
@@ -231,7 +289,7 @@
double[] limits = { 0.002, 0.002, 0.0001, 0.0001, 0.0005};
DiagBackward = new RKTrackFitDiag ( "Backward", "DCA", nbins, limits );
}
- DiagBackward.fill( fstatb, htin, vtdca, hlist.nDof(), vstart);
+ DiagBackward.fill( fstatb, htin, vtdca, hlist.nDof(), vstart, rktrk);
}
@@ -260,8 +318,8 @@
if ( sum > pathlimit ) break;
if ( s.inBounds(vtdca.vector(1)) ){
rkl.addHit( new RKHit( rktrk, s, vtdca, ps, sum ) );
- aida.cloud2D( "C r vs z").fill( vtdca.vector(1), s.radius*mmTocm );
- aida.cloud2D( "Both r vs z").fill( vtdca.vector(1), s.radius*mmTocm );
+ aida.cloud2D( "C r vs z",-1).fill( vtdca.vector(1), s.radius*mmTocm );
+ aida.cloud2D( "Both r vs z",-1).fill( vtdca.vector(1), s.radius*mmTocm );
/*
double r = s.radius*mmTocm;
System.out.printf ("Adding Cyl: %10.4f %10.4f %10.4f %10.4f %10.4f %10.4f\n",
@@ -288,9 +346,9 @@
} else{
System.out.println ("Can't get here ...");
}
- aida.cloud2D("Path length vs cz").fill(rktrk.cz(),xsum);
+ aida.cloud2D("Path length vs cz",-1).fill(rktrk.cz(),xsum);
aida.cloud1D("Sum").fill(xsum);
- aida.cloud2D("Hits vs cz").fill( rktrk.cz(), ncmeas );
+ aida.cloud2D("Hits vs cz",-1).fill( rktrk.cz(), ncmeas );
sum = 0;
for ( RKSurf s: rkgeom.getZ( rktrk.z0(), rktrk.cz() ) ){
@@ -301,8 +359,8 @@
double r = Math.sqrt(vtz.vector(0)*vtz.vector(0) + vtz.vector(1)*vtz.vector(1) );
if ( s.inBounds(r) ){
rkl.addHit( new RKHit( rktrk, s, vtz, ps, sum ) );
- aida.cloud2D( "Z r vs z").fill( s.zc*mmTocm, r );
- aida.cloud2D( "Both r vs z").fill( s.zc*mmTocm, r );
+ aida.cloud2D( "Z r vs z",-1).fill( s.zc*mmTocm, r );
+ aida.cloud2D( "Both r vs z",-1).fill( s.zc*mmTocm, r );
nmeas += s.rodim;
nzmeas += s.rodim;
/*
@@ -316,9 +374,9 @@
}
}
}
- aida.cloud2D("Z Hits vs cz").fill( rktrk.cz(), nzmeas );
- aida.cloud2D("Nmeas vs cz").fill( rktrk.cz(), nmeas );
- aida.cloud2D("Nz vs Nc").fill( ncmeas, nzmeas );
+ aida.cloud2D("Z Hits vs cz",-1).fill( rktrk.cz(), nzmeas );
+ aida.cloud2D("Nmeas vs cz",-1).fill( rktrk.cz(), nmeas );
+ aida.cloud2D("Nz vs Nc",-1).fill( ncmeas, nzmeas );
return rkl;
}
@@ -339,8 +397,8 @@
aida.histogram1D("d0", nbins, d0min, d0max).fill(t.d0());
aida.histogram1D("z0", nbins, z0min, z0max).fill(t.z0());
aida.histogram1D("p", nbins, 0., 50.).fill(t.p());
- aida.cloud2D("PCA0").fill(t.x0(),t.y0());
- aida.cloud2D("D0 vs Z0").fill(t.z0(),t.d0());
+ aida.cloud2D("PCA0",-1).fill(t.x0(),t.y0());
+ aida.cloud2D("D0 vs Z0",-1).fill(t.z0(),t.d0());
}
}