| Commit in lcsim/src/org/lcsim/contrib/JanStrube/standalone on MAIN | |||
| RaveProcessor.java | +117 | added 1.1 | |
| raveTest.py | +44 | added 1.1 | |
| MainLoop.py | +1 | -1 | 1.4 -> 1.5 |
| Tracks.py | +1 | -1 | 1.2 -> 1.3 |
| at/hephy/Converter/LcioToRaveObjects.java | +208 | added 1.1 | |
| +371 | -2 | ||
adding a modified version of the RaveProcessor and the converter temporarily to this repo
diff -N RaveProcessor.java --- /dev/null 1 Jan 1970 00:00:00 -0000 +++ RaveProcessor.java 29 Nov 2006 02:32:14 -0000 1.1 @@ -0,0 +1,117 @@
+import java.util.List;
+import java.util.Vector;
+import at.hephy.Rave.RaveConstantMagneticField;
+import at.hephy.Rave.RaveFactory;
+import at.hephy.Rave.RavePoint3D;
+import at.hephy.Rave.RavePropagator;
+import at.hephy.Rave.RaveVacuumPropagator;
+import at.hephy.Rave.RaveVertexContainer;
+import at.hephy.Rave.RaveVertex;
+import at.hephy.Rave.RaveTrackContainer;
+import at.hephy.Rave.RaveTrack;
+import at.hephy.Rave.RaveCovariance3D;
+import at.hephy.Converter.LcioToRaveObjects;
+import org.lcsim.util.Driver;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.EventHeader.LCMetaData;
+import org.lcsim.event.MCParticle;
+import hep.physics.particle.Particle;
+import org.lcsim.util.aida.AIDA;
+//import org.lcsim.contrib.JanStrube.tracking.Track;
+import org.lcsim.event.Track;
+import static java.lang.Math.sqrt;
+
+
+public class RaveProcessor extends Driver {
+
+ static {
+ System.loadLibrary("rave");
+ }
+
+ private LcioToRaveObjects converter_;
+ private RaveFactory factory_;
+ private AIDA aida = AIDA.defaultInstance();
+ private double _x = 0;
+ private double _y = 0;
+ private double _z = 0;
+
+ public class MyMagneticField extends RaveConstantMagneticField {
+ public RavePoint3D inTesla(RavePoint3D p) {
+ System.out.println( "[MyMagneticField] cool!!! the right ::inTesla called!!" );
+ return new RavePoint3D(0.,0.,0. );
+ }
+
+ public MyMagneticField () {
+ super(0.,0.,5.);
+ }
+ }
+
+public class MyPropagator extends RaveVacuumPropagator {
+ public MyPropagator () {
+ super ();
+ }
+
+ public RavePropagator clone() {
+ // System.out.println( "Cloning MyPropagator");
+ return new MyPropagator();
+ }
+}
+
+public RaveProcessor() {
+ converter_ = new LcioToRaveObjects();
+ factory_ = new RaveFactory( new MyMagneticField(), new MyPropagator(), "avf" );
+}
+
+
+public void processTracks(List<Track> trks) {
+ System.out.println("Processing collection.");
+ RaveTrackContainer rtc = new RaveTrackContainer();
+ for (Track trk : trks) {
+ RaveTrack t = converter_.convert(trk);
+ rtc.add(t);
+ }
+ RaveVertexContainer vtces = factory_.create(rtc);
+ long size = vtces.size();
+ System.out.println("Rave found " + size + " vertices.");
+ aida.cloud1D("number of vtx").fill(size);
+ for (int i=0; i<size; ++i) {
+ RaveVertex vtx = vtces.get(i);
+ RavePoint3D pos = vtx.position();
+ aida.cloud1D("vtx_x").fill(pos.x());
+ aida.cloud1D("vtx_y").fill(pos.y());
+ aida.cloud1D("vtx_z").fill(pos.z());
+ RaveCovariance3D error = vtx.error();
+ aida.cloud1D("vtx_x_err").fill(sqrt(error.dxx()));
+ aida.cloud1D("vtx_y_err").fill(sqrt(error.dyy()));
+ aida.cloud1D("vtx_z_err").fill(sqrt(error.dzz()));
+ aida.cloud1D("vtx_chi2").fill(vtx.chiSquared());
+ aida.cloud1D("vtx_nDOF").fill(vtx.ndf());
+ aida.cloud1D("vtx_x_pull").fill((_x - pos.x())/sqrt(error.dxx()));
+ aida.cloud1D("vtx_y_pull").fill((_y - pos.y())/sqrt(error.dyy()));
+ aida.cloud1D("vtx_z_pull").fill((_z - pos.z())/sqrt(error.dzz()));
+ }
+}
+
+public void process(EventHeader event) {
+ _x = 0;
+ _y = 0;
+ _z = 0;
+ for (MCParticle p : event.getMCParticles()) {
+ if (p.getGeneratorStatus() == Particle.INTERMEDIATE) {
+ _x = p.getEndPoint().x();
+ _y = p.getEndPoint().y();
+ _z = p.getEndPoint().z();
+ break;
+ }
+ }
+ List<List<Track>> trkcolls = event.get(Track.class);
+ for (List<Track> tracks : trkcolls) {
+ LCMetaData meta = event.getMetaData(tracks);
+ String name = meta.getName();
+ if (name=="TPCTracks" || name=="Tracks") {
+ processTracks(tracks);
+ }
+ }
+}
+}
+
diff -N raveTest.py --- /dev/null 1 Jan 1970 00:00:00 -0000 +++ raveTest.py 29 Nov 2006 02:32:14 -0000 1.1 @@ -0,0 +1,44 @@
+#! /usr/bin/env jython
+from org.lcsim.util import Driver
+from org.lcsim.util.aida import AIDA
+from org.lcsim.util.loop import LCIODriver
+from org.lcsim.util.loop import LCSimLoop
+from org.lcsim.mc.fast import MCFast
+from java.net import URL
+from org.lcsim.util.cache import FileCache
+from java.io import File
+#from ZvTubePlotter import ZvTubePlotter
+from java.lang import System
+#from FitterTestDriver import FitterTestDriver
+from org.lcsim.contrib.JanStrube.tracking import NewMCFastTrackDriver
+from java.lang import Boolean
+import RaveProcessor
+
+True = Boolean("true")
+False = Boolean("false")
+
+def main():
+ loop = LCSimLoop()
+# location = URL("ftp://ftp-lcd.slac.stanford.edu/lcd/ILC/ILC500/Zgamma/stdhep/pythia/pythiaZgamma.stdhep")
+# location = URL("ftp://ftp-lcd.slac.stanford.edu/lcd/ILC/singleParticle/stdhep/psi_mumu_Theta4-176_5-100GeV.stdhep")
+ location = URL("ftp://ftp-lcd.slac.stanford.edu/lcd/ILC/singleParticle/stdhep/K0S_pipi_Theta45-135_5-25Gev.stdhep")
+# location = URL("ftp://ftp-lcd.slac.stanford.edu/lcd/ILC/singleParticle/sidaug05/slcio/slic/tau_5pi_Theta20-90_10-200GeV_SLIC_v1r9p3_sidaug05.slcio")
+# location = URL("ftp://ftp-lcd.slac.stanford.edu/lcd/ILC/singleParticle/sidaug05/slcio/slic/tau_3pi_Theta20-90_10-200GeV_SLIC_v1r9p3_sidaug05.slcio")
+# location = URL("ftp://ftp-lcd.slac.stanford.edu/lcd/ILC/singleParticle/sidaug05/slcio/slic/psi_mumu_Theta4-176_5-100GeV_SLIC_v1r9p3_sidaug05.slcio")
+ cache = FileCache()
+ trackFile = cache.getCachedFile(location)
+ loop.setStdhepRecordSource(trackFile, "sidaug05")
+# loop.setLCIORecordSource(trackFile)
+ # no beamspotConstraint, no simple smearing
+ loop.add(MCFast(False, False))
+# loop.add(NewMCFastTrackDriver())
+ loop.add(RaveProcessor())
+# output = File(System.getProperty("user.home"),"fastmc.slcio")
+# loop.add(LCIODriver(output))
+ loop.loop(1000)
+ loop.dispose()
+
+ AIDA.defaultInstance().saveAs("fastmc.aida")
+
+if __name__ == "__main__":
+ main()
diff -u -r1.4 -r1.5 --- MainLoop.py 28 Oct 2006 01:11:33 -0000 1.4 +++ MainLoop.py 29 Nov 2006 02:32:14 -0000 1.5 @@ -43,7 +43,7 @@
# loop.add(ReconParticleTestDriver())
# output = File(System.getProperty("user.home"),"fastmc.slcio")
# loop.add(LCIODriver(output))
- loop.loop(-1)
+ loop.loop(1000)
loop.dispose()
AIDA.defaultInstance().saveAs("fastmc.aida")
diff -u -r1.2 -r1.3 --- Tracks.py 28 Oct 2006 01:11:33 -0000 1.2 +++ Tracks.py 29 Nov 2006 02:32:14 -0000 1.3 @@ -108,7 +108,7 @@
#self.aida.cloud1D('OldTrack/diffMomentum x').fill(momentum.x() - m*momentum_old.x())
#self.aida.cloud1D('OldTrack/diffMomentum y').fill(momentum.y() - m*momentum_old.y())
#self.aida.cloud1D('OldTrack/diffMomentum z').fill(momentum.z() - m*momentum_old.z())
- #self.aida.cloud1D('OldTrack distance to Origin').fill(dist_old)
+ self.aida.cloud1D('OldTrack distance to Origin').fill(dist_old)
swimmer_new.setTrack(newT)
# should be 0, because point is instantiated there
diff -N LcioToRaveObjects.java --- /dev/null 1 Jan 1970 00:00:00 -0000 +++ LcioToRaveObjects.java 29 Nov 2006 02:32:15 -0000 1.1 @@ -0,0 +1,208 @@
+package at.hephy.Converter;
+import at.hephy.Rave.RaveTrack;
+import at.hephy.Rave.RaveVector6D;
+import at.hephy.Rave.RaveCovariance6D;
+// import hep.lcio.event.Track;
+import org.lcsim.event.Track;
+import hep.physics.matrix.SymmetricMatrix;
+import org.lcsim.spacegeom.SpacePoint;
+import org.lcsim.spacegeom.SpaceVector;
+import Jama.Matrix;
+
+import static org.lcsim.contrib.JanStrube.tracking.LCIOTrackParameters.ParameterName.d0;
+import static org.lcsim.contrib.JanStrube.tracking.LCIOTrackParameters.ParameterName.phi0;
+import static org.lcsim.contrib.JanStrube.tracking.LCIOTrackParameters.ParameterName.omega;
+import static org.lcsim.contrib.JanStrube.tracking.LCIOTrackParameters.ParameterName.tanLambda;
+import static org.lcsim.contrib.JanStrube.tracking.LCIOTrackParameters.ParameterName.z0;
+
+import static org.lcsim.contrib.JanStrube.vtxFitter.Fitter.getSpatialDerivativeMatrix;
+import static org.lcsim.contrib.JanStrube.vtxFitter.Fitter.getMomentumDerivativeMatrix;
+import static org.lcsim.constants.Constants.fieldConversion;
+
+import static java.lang.Math.sin;
+import static java.lang.Math.cos;
+import static java.lang.Math.abs;
+
+public class LcioToRaveObjects
+{
+ public RaveTrack convert( org.lcsim.contrib.JanStrube.tracking.Track lct )
+ {
+ SpacePoint refPoint = lct.getReferencePoint();
+ SpacePoint r = lct.getPosition();
+ SpaceVector p = lct.getMomentum();
+ double _d0 = lct.getParameter(d0);
+ double _omega = lct.getParameter(omega);
+ double _phi0 = lct.getParameter(phi0);
+ double _tanLambda = lct.getParameter(tanLambda);
+ double _z0 = lct.getParameter(z0);
+ RaveVector6D s = new RaveVector6D (r.x(), r.y(), r.z(), p.x(), p.y(), p.z());
+ double signOmega = _omega / abs(_omega) ;
+ int charge = lct.getCharge();
+ double Bz = p.rxy() * _omega / charge / fieldConversion;
+ // assume that the Track doesn't have to be swum for now
+ double l = 0;
+ double _phi = _phi0 - _omega*l;
+ double pt = lct.getPt();
+
+ double [][] brJacob = {
+ { Math.sin(_phi0),
+ -Math.cos(_phi0), 0.0, 0.0, 0.0, 0.0 },
+ {(1/_omega-_d0)*cos(_phi0) - cos(_phi)/_omega, (1/_omega-_d0)*sin(_phi0) - sin(_phi)/_omega, 0
+ , -pt*sin(_phi), pt*cos(_phi), 0},
+ {1/_omega/_omega*(sin(_phi0)-sin(_phi)), 1/_omega/_omega*(cos(_phi0)-cos(_phi)), 0
+ , -pt/_omega*cos(_phi)+pt*l*sin(_phi), -pt/_omega*sin(_phi)-pt*l*cos(_phi), -pt/_omega*_tanLambda},
+ {0, 0, 1, 0, 0, 0},
+ {-cos(_phi)*l*_tanLambda/(1+_tanLambda*_tanLambda), -sin(_phi)*l*_tanLambda/(1+_tanLambda*_tanLambda), l
+ , -pt*sin(_phi)*_omega*l*_tanLambda/(1+_tanLambda*_tanLambda), pt*cos(_phi)*_omega*l*_tanLambda/(1+_tanLambda*_tanLambda), pt}
+ };
+
+ SymmetricMatrix bCov = lct.getErrorMatrix();
+ double[][] rCov= { {0.,0.,0.,0.,0.,0.}, {0.,0.,0.,0.,0.,0.},{0.,0.,0.,0.,0.,0.},{0.,0.,0.,0.,0.,0.},{0.,0.,0.,0.,0.,0.},{0.,0.,0.,0.,0.,0.} };
+ for ( int x=0; x< 6; x++ )
+ {
+ for ( int y=0; y<6; y++ )
+ {
+ for ( int z=0; z<5; z++ )
+ {
+ for ( int z2=0; z2<5; z2++ )
+ {
+ rCov[x][y]+=brJacob[z][x]*bCov.e(z, z2)*brJacob[z2][y];
+ };
+ };
+ };
+ };
+
+ RaveCovariance6D err = new RaveCovariance6D ( rCov[0][0], rCov[0][1], rCov[0][2],
+ rCov[1][1], rCov[1][2], rCov[2][2],
+ rCov[0][3], rCov[0][4], rCov[0][5],
+ rCov[1][3], rCov[1][4], rCov[1][5],
+ rCov[2][3], rCov[2][4], rCov[2][5],
+ rCov[3][3], rCov[3][4], rCov[3][5],
+ rCov[4][4], rCov[4][5], rCov[5][5] ) ;
+
+ RaveTrack trk = new RaveTrack ( s, err, charge, 0, "id" );
+ return trk;
+ }
+
+
+
+ public RaveTrack convert( Track lct )
+ {
+ /*
+ float bInvP = lct.getOmega() ;
+ float bTheta = lct.getTanLambda() ;
+ float bPhi = lct.getPhi() ; */
+ double bInvP = lct.getTrackParameter ( 0 );
+ double bTheta = lct.getTrackParameter ( 1 );
+ double bPhi = lct.getTrackParameter ( 2 );
+ double bRRef = lct.getReferencePoint()[0] ;
+ double bRPhiRef = lct.getReferencePoint()[1] ;
+ double bZRef = lct.getReferencePoint()[2] ;
+
+ RaveVector6D s = new RaveVector6D ( bRRef * Math.cos( bRPhiRef / bRRef ) ,
+ bRRef * Math.sin( bRPhiRef / bRRef ) ,
+ bZRef ,
+ Math.sin( bTheta ) * Math.cos( bPhi ) / Math.abs( bInvP ) ,
+ Math.sin( bTheta ) * Math.sin( bPhi ) / Math.abs( bInvP ) ,
+ Math.cos( bTheta ) / Math.abs( bInvP ) ) ;
+ double signInvP = bInvP / Math.abs( bInvP ) ;
+
+ double [][] brJacob = {
+ { (-1) * Math.sin( bRPhiRef / bRRef ),
+ Math.cos( bRPhiRef / bRRef ), 0.0, 0.0, 0.0, 0.0 },
+ { 0.0, 0.0, 1.0, 0.0, 0.0, 0.0 },
+ { 0., 0., 0., Math.cos( bTheta ) * Math.cos( bPhi ) / Math.abs( bInvP ) ,
+ Math.cos( bTheta ) * Math.cos( bPhi ) / Math.abs( bInvP ) ,
+ (-1.0) * Math.sin( bTheta ) / Math.abs( bInvP ) },
+ { 0.0, 0.0, 0.0,
+ (-1.0) * Math.sin( bTheta ) * Math.sin( bPhi ) / Math.abs( bInvP ),
+ Math.sin( bTheta ) * Math.cos( bPhi ) / Math.abs( bInvP ), 0. },
+ { 0., 0., 0.,
+ (-1.0) * signInvP * Math.sin( bTheta ) * Math.cos( bPhi ) / Math.pow( bInvP, 2 ),
+ (-1.0) * signInvP * Math.sin( bTheta ) * Math.sin( bPhi ) / Math.pow( bInvP, 2 ),
+ (-1.0) * signInvP * Math.cos( bTheta ) / Math.pow( bInvP, 2 ) }
+ };
+ /*
+ HepMatrix brJacob = new HepMatrix ( 6, 5 ) ;
+ brJacob[0][1] = 0.0 ;
+ brJacob[1][1] = 0.0 ;
+ brJacob[2][1] = 1.0 ;
+ brJacob[3][1] = 0.0 ;
+ brJacob[4][1] = 0.0 ;
+ brJacob[5][1] = 0.0 ;
+ brJacob[0][2] = 0.0 ;
+ brJacob[1][2] = 0.0 ;
+ brJacob[2][2] = 0.0 ;
+ brJacob[3][2] = Math.cos( bTheta ) * Math.cos( bPhi ) / Math.abs( bInvP ) ;
+ brJacob[4][2] = Math.cos( bTheta ) * Math.sin( bPhi ) / Math.abs( bInvP ) ;
+ brJacob[5][2] = (-1.0) * Math.sin( bTheta ) / Math.abs( bInvP ) ;
+ brJacob[0][3] = 0.0 ;
+ brJacob[1][3] = 0.0 ;
+ brJacob[2][3] = 0.0 ;
+ brJacob[3][3] = (-1.0) * Math.sin( bTheta ) * Math.sin( bPhi ) / Math.abs( bInvP ) ;
+ brJacob[4][3] = Math.sin( bTheta ) * Math.cos( bPhi ) / Math.abs( bInvP ) ;
+ brJacob[5][3] = 0.0 ;
+ brJacob[0][4] = 0.0 ;
+ brJacob[1][4] = 0.0 ;
+ brJacob[2][4] = 0.0 ;
+ brJacob[3][4] = (-1.0) * signInvP * Math.sin( bTheta ) * Math.cos( bPhi ) / Math.pow( bInvP, 2 ) ;
+ brJacob[4][4] = (-1.0) * signInvP * Math.sin( bTheta ) * Math.sin( bPhi ) / Math.pow( bInvP, 2 ) ;
+ brJacob[5][4] = (-1.0) * signInvP * Math.cos( bTheta ) / Math.pow( bInvP, 2 ) ;
+ */
+
+ SymmetricMatrix bCov = lct.getErrorMatrix();
+ /*
+ double[][] bCov={{ m[0], m[1], m[3], m[6], m[10]}, { m[1], m[2], m[4], m[7], m[11] },
+ { m[3], m[4], m[5], m[8], m[12] }, { m[6], m[7], m[8], m[9], m[13] },
+ { m[10], m[11], m[12], m[13], m[14] } };
+ HepSymMatrix bCov = new HepSymMatrix ( 5 ) ;
+ bCov[0][0] = lct.getCovMatrix()[ 0] ;
+ bCov[0][1] = lct.getCovMatrix()[ 1] ;
+ bCov[1][1] = lct.getCovMatrix()[ 2] ;
+ bCov[0][2] = lct.getCovMatrix()[ 3] ;
+ bCov[1][2] = lct.getCovMatrix()[ 4] ;
+ bCov[2][2] = lct.getCovMatrix()[ 5] ;
+ bCov[0][3] = lct.getCovMatrix()[ 6] ;
+ bCov[1][3] = lct.getCovMatrix()[ 7] ;
+ bCov[2][3] = lct.getCovMatrix()[ 8] ;
+ bCov[3][3] = lct.getCovMatrix()[ 9] ;
+ bCov[0][4] = lct.getCovMatrix()[10] ;
+ bCov[1][4] = lct.getCovMatrix()[11] ;
+ bCov[2][4] = lct.getCovMatrix()[12] ;
+ bCov[3][4] = lct.getCovMatrix()[13] ;
+ bCov[4][4] = lct.getCovMatrix()[14] ;
+ */
+
+ /*
+ HepMatrix rCov = new HepMatrix ( 6, 6 ) ;
+ rCov = brJacob * bCov * brJacob.T() ;
+ */
+ double[][] rCov= { {0.,0.,0.,0.,0.,0.}, {0.,0.,0.,0.,0.,0.},{0.,0.,0.,0.,0.,0.},{0.,0.,0.,0.,0.,0.},{0.,0.,0.,0.,0.,0.},{0.,0.,0.,0.,0.,0.} };
+ for ( int x=0; x< 6; x++ )
+ {
+ for ( int y=0; y<6; y++ )
+ {
+ for ( int z=0; z<5; z++ )
+ {
+ for ( int z2=0; z2<5; z2++ )
+ {
+ rCov[x][y]+=brJacob[z][x]*bCov.e(z, z2)*brJacob[z2][y];
+ };
+ };
+ };
+ };
+
+ RaveCovariance6D err = new RaveCovariance6D ( rCov[0][0], rCov[0][1], rCov[0][2],
+ rCov[1][1], rCov[1][2], rCov[2][2],
+ rCov[0][3], rCov[0][4], rCov[0][5],
+ rCov[1][3], rCov[1][4], rCov[1][5],
+ rCov[2][3], rCov[2][4], rCov[2][5],
+ rCov[3][3], rCov[3][4], rCov[3][5],
+ rCov[4][4], rCov[4][5], rCov[5][5] ) ;
+
+ // int charge = signed ( copysign ( 1.0, (-1.0)*signInvP ) );
+ int charge = (int) ( -1.0 * signInvP / Math.abs ( signInvP ) );
+ RaveTrack trk = new RaveTrack ( s, err, charge, 0, "id" );
+ return trk;
+ }
+}