lcsim/src/org/lcsim/contrib/RobKutschke/TRFSelfTest/fit
diff -N FitTest_sid00Driver.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ FitTest_sid00Driver.java 8 Jun 2009 06:03:17 -0000 1.1
@@ -0,0 +1,455 @@
+package org.lcsim.contrib.RobKutschke.TRFSelfTest.fit;
+
+import java.util.ArrayList;
+import java.util.Iterator;
+import java.util.List;
+import java.util.Set;
+import Jama.Matrix;
+
+import org.lcsim.event.EventHeader;
+
+import org.lcsim.geometry.Detector;
+
+import org.lcsim.util.Driver;
+import org.lcsim.util.aida.AIDA;
+import org.lcsim.constants.Constants;
+
+import org.lcsim.recon.tracking.trfbase.TrackVector;
+import org.lcsim.recon.tracking.trfbase.TrackError;
+import org.lcsim.recon.tracking.trfbase.ETrack;
+import org.lcsim.recon.tracking.trfbase.Hit;
+import org.lcsim.recon.tracking.trfbase.VTrack;
+import org.lcsim.recon.tracking.trfbase.PropStat;
+import org.lcsim.recon.tracking.trfbase.PropDir;
+import org.lcsim.recon.tracking.trfbase.Propagator;
+
+import org.lcsim.recon.tracking.trfcyl.SurfCylinder;
+import org.lcsim.recon.tracking.trfzp.SurfZPlane;
+
+import org.lcsim.recon.tracking.trffit.HTrack;
+
+import org.lcsim.event.EventHeader.LCMetaData;
+import org.lcsim.event.SimTrackerHit;
+
+import org.lcsim.math.chisq.ChisqProb;
+
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.generator.RKTrack;
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.generator.RKTrackGen;
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.generator.RKTrackGenParams;
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.generator.toTRF;
+
+import org.lcsim.contrib.RobKutschke.TRF.trffit.FullFitKalman;
+
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.util.RKDebug;
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.geom.RKGeom;
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.geom.RKSurf;
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.hits.RKHit;
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.hits.RKHitList;
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.hits.RKMakeHits;
+
+import org.lcsim.contrib.RobKutschke.ToyConfig.*;
+
+/**
+ * Driver for exercises learning how to use the TRF code.
+ *
+ *@author $Author: kutschke $
+ *@version $Id: FitTest_sid00Driver.java,v 1.1 2009/06/08 06:03:17 kutschke Exp $
+ *
+ * Date $Date: 2009/06/08 06:03:17 $
+ *
+ */
+
+
+public class FitTest_sid00Driver extends Driver
+{
+
+ private AIDA aida = AIDA.defaultInstance();
+
+ // The run time configuration system.
+ ToyConfig config;
+
+ // Geometry summary extracted from GeomConverter information.
+ RKGeom rkgeom = null;
+
+ // Internal track generate for self test and its parameter set.
+ RKTrackGen gen = null;
+ RKTrackGenParams param = null;
+
+ // Internal self test code to generate hits on tracks.
+ // Includes energy loss and multiple scattering.
+ private RKMakeHits HitMaker = null;
+
+ // Propagator and fitter.
+ Propagator prop = null;
+ FullFitKalman fitk = null;
+
+ // Starting values of covariance matrices.
+ private TrackError vstartc_in = new TrackError();
+ private TrackError vstartz_in = new TrackError();
+ private TrackError vstart_out = new TrackError();
+
+ // Instances of the diagnostic package, one for each direction of fit.
+ private RKTrackFitDiag DiagForward = null;
+ private RKTrackFitDiag DiagBackward = null;
+
+ // Information copied from the runtime configuration system.
+ private int nSelfTest = 1;
+ private boolean doms = false;
+ private boolean doeloss = false;
+ private double msfac = 1.0;
+ //private double dedxscale = 7.5;
+ private double dedxsigma = 0.00;
+
+ private int FitTestDebugLevel = 0;
+
+ // The only constructor.
+ public FitTest_sid00Driver(){
+
+ // Extract information from the run time configuration system.
+ try{
+ ToyConfig config = ToyConfig.getInstance();
+
+ // Assorted control parameters.
+ nSelfTest = config.getInt("nSelfTest");
+ doms = config.getBoolean("DoMS");
+ doeloss = config.getBoolean("DoELoss");
+ msfac = config.getDouble("MsFac");
+ //dedxscale = config.getDouble("dEdXScale");
+ dedxsigma = config.getDouble("dEdXSigma");
+
+ // Limits for generation of track parameters.
+ param = new RKTrackGenParams( config );
+ FitTestDebugLevel = config.getInt("FitTestDebugLevel", FitTestDebugLevel );
+
+ // Starting values of the covariance matrices.
+ double[] v1 = config.getArrayDouble("Vstartc_in");
+ double[] v2 = config.getArrayDouble("Vstartz_in");
+ double[] v3 = config.getArrayDouble("Vstart_out");
+ for ( int i=0; i<5; ++i ){
+ vstartc_in.set(i,i,v1[i]);
+ vstartz_in.set(i,i,v2[i]);
+ vstart_out.set(i,i,v3[i]);
+ }
+
+
+ } catch (ToyConfigException e){
+ System.out.println (e.getMessage() );
+ System.out.println ("Stopping now." );
+ System.exit(-1);
+ }
+
+ // Set flags in the debug system.
+ RKDebug.Instance().setDoMs(doms);
+ RKDebug.Instance().setDoEloss(doeloss);
+ RKDebug.Instance().setMsFac(msfac);
+ RKDebug.Instance().setPrintOutliers(false);
+
+ }
+
+ public void detectorChanged( Detector detector) {
+
+ // Build summary of the geometry information.
+ rkgeom = new RKGeom(detector );
+
+ // Update B field inside the track generation parameters.
+ param.setbz(rkgeom.getBz());
+ System.out.println (param);
+
+ // Seed for track random number generator.
+ //long seed = -987835122;
+ long seed = -398783512;
+
+ // An increment used to generate independent streams of random
+ // numbers - I have no proof that this really works.
+ long seedIncrement = 876633217;
+
+ // Instantiate a new track generator for this parameter set and geometry.
+ gen = new RKTrackGen(param,seed);
+
+ // Seed for hit generator.
+ seed += seedIncrement;
+ RKHit.setSeed(seed);
+
+ // Instantiate new making code.
+ HitMaker = new RKMakeHits( rkgeom, getRandom() );
+
+ // Instantiate and configure the kalman fitter.
+ prop = rkgeom.newPropagator();
+ fitk = new FullFitKalman(prop);
+ fitk.setDoMs(doms);
+ fitk.setDoEloss(doeloss);
+
+ }
+
+ // Do the test only once.
+ private boolean dotest = true;
+
+ public void process(EventHeader event){
+
+ if ( dotest ){
+ dotest = false;
+
+ // Execute the self test.
+ for ( int i=0; i<nSelfTest; ++i ){
+ RKDebug.Instance().setTrack(i);
+ SelfTestOneTrack(i);
+ }
+ RKDebug.Instance().setTrack(-99);
+ }
+
+ }
+
+ /**
+ * Perform the self test on one track.
+ *
+ * The argument is a track number used for debug printout.
+ *
+ */
+ private void SelfTestOneTrack( int itrk){
+
+ // Generate a track and plot monitoring info.
+ RKTrack rktrk = gen.newRKTrack();
+
+ // Heartbeat.
+ if ( (itrk < 10) || ( itrk%100 == 0) ){
+ System.out.printf ("New track heartbeat: %6d %12.6f %12.6f\n", itrk, rktrk.pt(), rktrk.cz() );
+ }
+
+ // Load track info into the debug system for use elsewhere.
+ RKDebug.Instance().setRKTrack(rktrk);
+
+ // Generate hits on this track in my own format.
+ // Save lots of generator info with each hit.
+ RKHitList hlist = HitMaker.GenerateMixedOneArcHits( rktrk );
+
+ // Check for sufficient hits.
+ if ( ! hlist.fitable() ){
+ double acz = Math.abs(rktrk.cz());
+ aida.cloud1D("Summary/abs(Cz) for tracks with too few hits").fill(acz);
+ aida.cloud1D("Summary/Pt for tracks with too few hits").fill(rktrk.pt());
+ aida.cloud2D("Summary/Pt vs cz for tracks with too few hits").fill(acz,rktrk.pt());
+ return;
+ }
+
+ // Get starting track parameters and covariance matrix for outward fit ( at PCA to z axis ).
+ toTRF trftrk = new toTRF(rktrk);
+ VTrack vtdca = trftrk.atDcaZaxis();
+ TrackError vstart = vstart_out;
+ RKDebug.Instance().setStartType(0);
+
+ // Prep for the outward fit.
+ ETrack et = new ETrack( vtdca.surface().newPureSurface(), vtdca.vector(), vstart );
+ HTrack ht = new HTrack( et );
+ for ( RKHit h: hlist.getForward() ){
+ Hit hh = h.MakeHit();
+ if ( hh != null ) ht.addHit(hh);
+ }
+
+ int nleads = hlist.nLeadingStrips();
+ int firstz = hlist.firstZType();
+
+ if ( nleads > 0 & firstz == 2){
+ aida.cloud1D("Test/cz for Type 1 bug(fixed)").fill(rktrk.cz());
+ }
+
+ // Do the outward 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();
+
+ // Diagnostics for forward fit.
+ if ( DiagForward == null ) {
+ int nbins = 50;
+ 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, rktrk);
+
+ // Prep for the inward fit.
+ if ( vtend.surface().pureType().equals(SurfCylinder.staticType()) ){
+ vstart = vstartc_in;
+ RKDebug.Instance().setStartType(1);
+ } else{
+ vstart = vstartz_in;
+ RKDebug.Instance().setStartType(2);
+ }
+
+ // Seed inward fit with the output of the outward fit or
+ // with the generated values.
+ ETrack etin = new ETrack( vtend.surface().newPureSurface(),
+ vtend.vector(),
+ //ht.newTrack().vector(),
+ vstart );
+
+ // When the starting surface is a zplane, we need to tell the ETrack
+ // which direction it is forward.
+ if ( vtend.surface().type() == SurfZPlane.staticType() ){
+ if ( rktrk.cz() > 0 ){
+ etin.setForward();
+ } else{
+ etin.setBackward();
+ }
+ }
+
+ double dz = hlist.deltaZ();
+ int nskip = 0;
+ int nskip0 = 0;
+ double askip = nskip0;
+
+ /*
+ if ( dz > 1000. ) {
+ RKSurf s0 = hlist.getBackward().get(0).getSurface();
+ if ( s0.getType() == RKSurf.type_zdisc ){
+ ++nskip;
+ }
+ RKSurf s1 = hlist.getBackward().get(0).getSurface();
+ if ( s1.getType() == RKSurf.type_zdisc ){
+ ++nskip;
+ }
+ }
+ nskip0 = nskip;
+ askip = nskip0;
+ aida.histogram1D("Test/nskip",10,0.,10.).fill(askip);
+ */
+
+ // Add hits to the track.
+ HTrack htin = new HTrack( etin );
+ for ( RKHit h: hlist.getBackward() ){
+ if ( nskip > 0 ){
+ --nskip;
+ continue;
+ }
+ Hit hh = h.MakeHit();
+ if ( hh != null ) htin.addHit(hh);
+ }
+
+ // Do inward fit.
+ int fstatb = fitk.fitBackward(htin);
+ if ( fstatb != 0 ){
+ System.out.println("Backwards status: " +
+ fstatb + " " +
+ hlist.nHits() + " " +
+ hlist.nCyl() + " " +
+ hlist.nZp() + " " +
+ rktrk.cz() );
+ double acz = Math.abs(rktrk.cz());
+ if ( fstatb < 100 ){
+ aida.cloud1D("/Bugs/FailedFit/Backward: abs(cz) for failed propagation").fill(acz);
+ } else{
+ aida.cloud1D("/Bugs/FailedFit/Backward: abs(cz) for failed addhit").fill(acz);
+ }
+ }
+ if ( fstatb != 0 ){
+ if ( FitTestDebugLevel > 0 ){
+ System.out.println ("Bad fit: "
+ + itrk + " "
+ + fstatb + " "
+ + rktrk.cz() + " "
+ + hlist.nHits() + " "
+ + hlist.nZ() + " "
+ + hlist.nDof()
+ );
+ for ( RKHit h: hlist.getBackward() ){
+ System.out.printf ( "%3d %10.4f %10.4f %14.8f\n",
+ h.getSurface().getType(),
+ h.getSurface().radius,
+ h.getSurface().zc,
+ h.getPsi()
+ );
+ }
+ }
+ aida.cloud2D("/Bugs/FailedFit/Bad fit nZ vs nCyl").fill( hlist.nCyl(), hlist.nZp() );
+ aida.histogram1D( "/Bugs/FailedFit/Bad fit: Start type",5,0.,5.).fill(RKDebug.Instance().getStartType());
+ }
+
+
+
+ // Diagnostics for inward fit.
+ if ( DiagBackward == null ) {
+ int nbins = 50;
+ double[] limits = { 0.002, 0.002, 0.0001, 0.001, 0.0005};
+ DiagBackward = new RKTrackFitDiag ( "Backward", "DCA", nbins, limits );
+ }
+ DiagBackward.fill( fstatb, htin, vtdca, hlist.nDof()-nskip0, vstart, rktrk);
+
+ // Alternate: comapare at inner most hit.
+ //VTrack vtin = hlist.innerMostHit().getVTrack();
+ //DiagBackward.fill( fstatb, htin, vtin, hlist.nDof()-nskip0, vstart, rktrk);
+
+ ChisqProb prob = new ChisqProb();
+ double cl = prob.gammq( hlist.nDof(), htin.chisquared() );
+
+ int pat = hlist.Pattern();
+ double cz = rktrk.cz();
+ double acz = Math.abs(cz);
+ if ( cl < 0.001 | dz > 1000. ) {
+ System.out.printf ( "Deltaz: %15.6f %15.6f %15.6f %4d\n",
+ dz,
+ cl,
+ rktrk.cz(),
+ pat
+ );
+ //System.out.printf("Cz: %15.6f %15.6f\n",
+ // rktrk.cz(),
+ // hlist.deltaZ()
+ // );
+ hlist.PrintBackward();
+ }
+
+ if ( dz > 1000. ) aida.histogram1D( "Test/cl for dz gt 1000.", 200, 0., 1. ).fill(cl);
+ if ( dz < 1000. & dz > 0. ) {
+ aida.histogram1D( "Test/cl for other dz gt 0.", 200, 0., 1. ).fill(cl);
+ aida.histogram1D( "Test/dz for other dz gt 0.", 200, 0., 1000. ).fill(dz);
+ }
+
+
+ if ( pat == 3 ){
+ aida.histogram1D("Test/pattern 3", 200, 0.78, 0.81).fill(acz);
+ aida.histogram1D("Test/pattern 3 or 4", 200, 0.78, 0.81).fill(acz);
+ aida.cloud2D("/Test/pattern 3: cl vs cz").fill( acz, cl);
+ if ( dz > 1000. ){
+ aida.cloud2D("/Test/pattern 3: cl vs cz dz>1000.").fill( acz, cl);
+ }
+ } else if ( pat == 4 ){
+ aida.histogram1D("Test/pattern 4", 200, 0.78, 0.81).fill(acz);
+ aida.histogram1D("Test/pattern 3 or 4", 200, 0.78, 0.81).fill(acz);
+ aida.histogram1D("Test/Dz for pattern 4", 100, 0., 2000.).fill(dz);
+ aida.cloud2D("/Test/pattern 4: cl vs cz").fill( acz, cl);
+ if ( dz > 1000. ){
+ aida.cloud2D("/Test/pattern 4: cl vs cz dz>1000.").fill( acz, cl);
+ } else if ( dz > 0 && dz < 1000 ){
+ aida.cloud2D("/Test/pattern 4: cl vs cz 0<=dz<1000.").fill( acz, cl);
+ } else{
+ aida.cloud2D("/Test/pattern 4: cl vs cz z==0.").fill( acz, cl);
+ }
+ } else if ( pat == 5 ){
+ aida.histogram1D("Test/Dz for pattern 5", 100, 0., 2000.).fill(dz);
+ aida.histogram1D("Test/pattern 5", 200, 0.78, 0.81).fill(acz);
+ }
+
+ // Check that diag values of starting covariance matrices are large enough.
+ // Compare final cov of outward fit with starting cov of inward fit, and vice versa.
+ DiagBackward.checkStartCov( ht, vstart);
+ DiagForward.checkStartCov ( htin, vstart_out);
+
+ }
+
+ protected void endOfData(){
+
+ // Check that the track parameter pull distributions are unit gaussians.
+ DiagForward.PullSummary();
+ DiagBackward.PullSummary();
+ }
+
+
+}
+
lcsim/src/org/lcsim/contrib/RobKutschke/TRFSelfTest/fit
diff -N RKTrackFitDiag.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ RKTrackFitDiag.java 8 Jun 2009 06:03:17 -0000 1.1
@@ -0,0 +1,287 @@
+package org.lcsim.contrib.RobKutschke.TRFSelfTest.fit;
+
+import org.lcsim.util.aida.AIDA;
+import hep.aida.IHistogram1D;
+import hep.aida.ICloud1D;
+import hep.aida.ICloud2D;
+
+import org.lcsim.recon.tracking.trfbase.VTrack;
+import org.lcsim.recon.tracking.trfbase.ETrack;
+import org.lcsim.recon.tracking.trfbase.TrackError;
+
+import org.lcsim.recon.tracking.trfutil.TRFMath;
+
+import org.lcsim.recon.tracking.trfcyl.SurfCylinder;
+import org.lcsim.recon.tracking.trfzp.SurfZPlane;
+import org.lcsim.recon.tracking.trfdca.SurfDCA;
+
+import org.lcsim.recon.tracking.trffit.HTrack;
+
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.generator.RKTrack;
+import org.lcsim.contrib.RobKutschke.TRFSelfTest.util.RKDebug;
+
+import org.lcsim.math.chisq.ChisqProb;
+
+
+/**
+ *
+ * Make some diagnostic histograms a fitted track.
+ *
+ *@author $Author: kutschke $
+ *@version $Id: RKTrackFitDiag.java,v 1.1 2009/06/08 06:03:17 kutschke Exp $
+ *
+ * Date $Date: 2009/06/08 06:03:17 $
+ *
+ */
+
+public class RKTrackFitDiag{
+
+ public RKTrackFitDiag( String tag, String type, int nbins, double[] limits ){
+ this.tag = tag;
+
+ // Directory into which histograms and clouds will go.
+ tagdir = "/" + tag;
+
+ String pwd = aida.tree().pwd();
+ aida.tree().mkdir(tagdir);
+ aida.tree().cd(tagdir);
+
+ // Book histograms.
+ for ( int i=0; i<5; ++i ){
+ resid[i] = aida.histogram1D(tag + " Residual " +i, nbins, -limits[i], limits[i] );
+ //sigma[i] = aida.histogram1D(tag + " Sigma " +i, nbins, 0., limits[i] );
+ sigma[i] = aida.cloud1D(tag + " Sigma " +i, -1);
+ pull[i] = aida.histogram1D(tag + " Pull " +i, nbins, -5., 5.);
+ ratio[i] = aida.histogram1D(tag + " E Ratio " +i, nbins, -10., 10.);
+ start[i] = aida.histogram1D(tag + " log_10(V Start(" +i +"))", 60, -5., 5.);
+ isAzimuth[i] = false;
+
+ // Force an entry at 0 in the sigma plots to help make a sensible scale.
+ sigma[i].fill(0.);
+ }
+ cl = aida.histogram1D( tag + " Confidence Level", nbins, 0, 1.);
+ czbadcl = aida.histogram1D( tag + " Cz for bad CL", 1000, -1., 1.);
+ nDof = aida.histogram1D( tag + " nDof",100,0.,50.);
+ errorCode = aida.histogram1D( tag + " Error codes", 20, 0., 20.);
+ sType = aida.histogram1D( tag + " Comparison Surface Type", 10, 0., 10.);
+ startType = aida.histogram1D( tag + " Start Type", 5, 0., 5.);
+
+ // Mark which parameters are azimuths.
+ if ( type.compareToIgnoreCase("Cyl") == 0 ){
+ isAzimuth[0] = true;
+ isAzimuth[2] = true;
+ } else if ( type.compareToIgnoreCase("DCA") == 0 ){
+ isAzimuth[2] = true;
+ }
+
+ aida.tree().cd(pwd);
+ }
+
+ public void fill( int fstat, HTrack ht, VTrack vt, int ndof, TrackError vstart, RKTrack rkt){
+
+ String pwd = aida.tree().pwd();
+ aida.tree().cd(tagdir);
+
+ if ( fstat != 0 ){
+ errorCode.fill(5.5);
+ aida.tree().cd(pwd);
+ return;
+ }
+
+ double cz = rkt.cz();
+
+ // Compute confidence level.
+ ChisqProb prob = new ChisqProb();
+ double chisq = ht.chisquared();
+
+ // Why does this happen?
+ if ( chisq < 0. ){
+ System.out.println( "Error negative chisquared: "
+ + RKDebug.Instance().getTrack() + " "
+ + ndof + " " + ht.chisquared() );
+
+ errorCode.fill(15.5);
+ if ( BadChisq == null ){
+ BadChisq = aida.cloud1D( "/Bugs/" + tag + ": Negative Chisquared");
+ }
+ BadChisq.fill(chisq);
+ chisq = 0.;
+
+ }
+
+ // Confidence level of the fit.
+ double c = prob.gammq( ndof, chisq );
+
+ // Tag the troublesome region in which we have bad inwards fits.
+ boolean czcut = ( Math.abs(cz)>0.65 && Math.abs(cz)<0.85);
+
+ double[] pulls = { 0., 0., 0., 0., 0.};
+
+ // Fill per parameter histograms.
+ ETrack et = ht.newTrack();
+ for ( int j=0; j<5; ++j ){
+
+ double res = et.vector(j) - vt.vector(j);
+ double sigsq = et.error(j,j);
+ if ( sigsq < 0. ){
+ errorCode.fill( j+0.5);
+ continue;
+ }
+ double sig = Math.sqrt(sigsq);
+
+ double sigsqstart = vstart.get(j,j);
+
+ // For angles, deal with the seam.
+ if ( isAzimuth[j] ){
+ res = TRFMath.fmod2( res, TRFMath.TWOPI);
+ }
+
+ double pul = res/sig;
+
+ resid[j].fill(res);
+ sigma[j].fill(sig);
+ pull [j].fill(pul);
+ start[j].fill(Math.log10(sigsqstart));
+
+ pulls[j] = pul;
+
+
+ if ( Math.abs(pul) > 5. & RKDebug.Instance().getPrintOutliers() ){
+ System.out.printf ("Outlier: track, in/out, parameter, residual, pull, cz: %3d %-9s %3d %12.6f %10.2f %10.4f\n",
+ RKDebug.Instance().getTrack(),
+ tag,
+ j,
+ res,
+ pul,
+ rkt.cz()
+ );
+ }
+
+ }
+
+
+
+
+ // Fill remaining histograms.
+ cl.fill(c);
+ nDof.fill(ndof);
+
+ if ( c < 0.005 ) {
+ czbadcl.fill( cz );
+ if ( czcut ){
+ startType.fill( RKDebug.Instance().getStartType() );
+ }
+ }
+
+ if ( c < 0.001 ) {
+ if ( tag.compareTo("Backward") == 0 ){
+ System.out.printf ("Pulls: %3d | %9.6f | %8.3f %8.3f %8.3f %8.3f %8.3f\n",
+ RKDebug.Instance().getTrack(),
+ c,
+ pulls[0],
+ pulls[1],
+ pulls[2],
+ pulls[3],
+ pulls[4]
+ );
+ }
+ }
+
+
+ // Type of the surface at which the comparison is made.
+ int stype = 0;
+ if ( vt.surface().type().compareTo(SurfCylinder.staticType()) == 0 ){
+ stype = 1;
+ } else if ( vt.surface().type().compareTo(SurfZPlane.staticType()) == 0 ){
+ stype = 2;
+ } else if ( vt.surface().type().compareTo(SurfDCA.staticType()) ==0 ){
+ stype = 3;
+ }
+ sType.fill(stype+0.5);
+
+ aida.tree().cd(pwd);
+ }
+
+ // Check that starting values of the covariance matrices are large enough
+ // compared to the final value of the other way fit.
+ public void checkStartCov( HTrack ht, TrackError vstart ){
+ ETrack et = ht.newTrack();
+ for ( int j=0; j<5; ++j ){
+
+ double sigsq = et.error(j,j);
+ if ( sigsq < 0. ) continue;
+ double sig = Math.sqrt(sigsq);
+
+ double sigsq0 = vstart.get(j,j);
+ if ( sigsq0 < 0. ) continue;
+ double sig0 = Math.sqrt(sigsq0);
+
+ double r= -9.9;
+ if ( sig0 > 0. ){
+ r = sig/sig0;
+ }
+ ratio[j].fill(Math.log10(r));
+ }
+
+ }
+
+ public void PullSummary(){
+ String pwd = aida.tree().pwd();
+ aida.tree().cd(tagdir);
+
+ IHistogram1D pulSum = aida.histogram1D(tag + ": Pull Summary",20, -5., 5.);
+ ICloud1D outSum = aida.cloud1D(tag + ": Out of Bounds Entries");
+ for ( int i=0; i<pull.length; ++i ){
+ double mean = pull[i].mean();
+ double rms = pull[i].rms();
+ int all = pull[i].allEntries();
+ int out = pull[i].extraEntries();
+ int n = all-out;
+
+ double mpull = mean*Math.sqrt(n);
+ double rpull = (rms-1.0)*Math.sqrt(2*n);
+
+ pulSum.fill(mpull);
+ pulSum.fill(rpull);
+ pulSumAll.fill(mpull);
+ pulSumAll.fill(rpull);
+ if ( out > 0 ) {
+ outSum.fill(out);
+ outSumAll.fill(out);
+ }
+
+ }
+ aida.tree().cd(pwd);
+ }
+
+ private AIDA aida = AIDA.defaultInstance();
+
+ // Copy of the idendifier; used for printout.
+ private String tag = null;
+ private String tagdir = null;
+
+ // Histograms for residuals, pulls, sigams, ndof, confidence level and error codes.
+ private IHistogram1D[] resid = new IHistogram1D[5];
+ private IHistogram1D[] pull = new IHistogram1D[5];
+ //private IHistogram1D[] sigma = new IHistogram1D[5];
+ private ICloud1D[] sigma = new ICloud1D[5];
+ private IHistogram1D[] start = new IHistogram1D[5];
+ private IHistogram1D[] ratio = new IHistogram1D[5];
+ private IHistogram1D cl = null;
+ private IHistogram1D czbadcl = null;
+ private IHistogram1D nDof = null;
+ private IHistogram1D sType = null;
+ private IHistogram1D startType = null;
+ private ICloud1D BadChisq = null;
+
+ private IHistogram1D errorCode = null;
+
+ private IHistogram1D pulSumAll = aida.histogram1D("Summary/All Fits Pull Summary",20,-5.,5.);
+ private ICloud1D outSumAll = aida.cloud1D("Summary/All Fits Out of Bounds Pulls");
+
+
+ // Identify which track parameters are azimuthal angles.
+ // For these parameters we need to deal with the seam.
+ private boolean[] isAzimuth = new boolean[5];
+
+}
