Moller skimmer for 2019

recon/src/main/java/org/hps/recon/skims/MollerSkimmer.java

@@ -1,40 +1,164 @@
 package org.hps.recon.skims;
+import static java.lang.Math.abs;
 
+import java.util.List;
 import java.util.Set; 
+import java.io.BufferedReader;
+import java.io.IOException;
+import java.io.InputStream;
+import java.io.InputStreamReader;
 
-import org.lcsim.event.EventHeader;
+import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.Hep3Vector;
 
+import org.lcsim.event.EventHeader;
+import org.hps.recon.ecal.cluster.ClusterUtilities;
+import org.hps.recon.particle.ReconParticleDriver;
+import org.hps.recon.tracking.TrackData;
+import org.hps.recon.tracking.TrackType;
+import org.hps.record.epics.EpicsData;
+import org.hps.record.scalers.ScalerData;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.ReconstructedParticle;
+import org.lcsim.event.Vertex;
 
 public class MollerSkimmer extends Skimmer {
-    private String _MollerCandidateCollectionName = "UnconstrainedMollerCandidates";
-    private double _clusterTimingCut = 20.0; // only used if _tight is true
-    private double _v0Chi2Cut = 100.0;
-    private double _trackChi2Cut = 80.0;
-    private double _trackDtCut = 20.0;
-    private double _trackPMax = 0.9;
-    private double _v0PMax = 1.4;
-    private int    _nHitsMin=10;
-
+    private String _MollerCandidateCollectionName = "UnconstrainedMollerCandidates_KF";
+    private String _MollerVertexCollectionName = "UnconstrainedMollerVertices_KF";
+    private double _vtxChi2Cut = 100.0;
+    private double _trackChi2Cut = 30.0;
+    private double _trackDtCut = 20.0; 
+    private int    _nHitsMin=9;
+    private boolean _debug=false;    
+    private int totalMollers=0; 
+    private int totalMollersPassing=0; 
+
     @Override
     public boolean passSelection(EventHeader event){
-	System.out.println(this.getClass().getName()+":: in pass selection"); 
-	boolean pass=true; 
+
+	if(_debug)
+	    System.out.println(this.getClass().getName()+":: in pass selection"); 
+	incrementEventProcessed();
 
+	if (!event.hasCollection(ReconstructedParticle.class, _MollerCandidateCollectionName)) {
+	     return false; 
+        }	
+	if (!event.hasCollection(Vertex.class, _MollerVertexCollectionName)) {
+	     return false; 
+        }
+
+        List<ReconstructedParticle> V0Candidates = event.get(ReconstructedParticle.class, _MollerCandidateCollectionName);
+	List<Vertex> V0Vertexes= event.get(Vertex.class, _MollerVertexCollectionName);
+
+	if(V0Candidates.size() != V0Vertexes.size())
+	    System.out.println(this.getClass().getName()+":: Number of Vertexes = "+V0Vertexes.size()+
+			       "; number of candidates = "+V0Candidates.size());
+
+        int nMollers = 0; 
+	totalMollers += V0Candidates.size();
+        for (ReconstructedParticle v0 : V0Candidates) {
+
+            ReconstructedParticle eleTop = v0.getParticles().get(ReconParticleDriver.MOLLER_TOP);
+            ReconstructedParticle eleBot = v0.getParticles().get(ReconParticleDriver.MOLLER_BOT);
+
+            if (v0.getStartVertex().getChi2() > _vtxChi2Cut) {
+		if(_debug)System.out.println(this.getClass().getName()+"::  failed vertex chi2");
+                continue;
+            }
+	    if(eleTop.getTracks().get(0).getTrackerHits().size()<_nHitsMin
+               || eleBot.getTracks().get(0).getTrackerHits().size()<_nHitsMin){
+		if(_debug)System.out.println(this.getClass().getName()+"::  failed nHitsMin "+eleTop.getTracks().get(0).getTrackerHits().size()+"  "+eleBot.getTracks().get(0).getTrackerHits().size()+" nHitsMin = "+_nHitsMin);
+                continue;
+            }
+            if ((eleTop.getTracks().get(0).getChi2()/eleTop.getTracks().get(0).getNDF()) > _trackChi2Cut
+                    || (eleBot.getTracks().get(0).getChi2()/eleBot.getTracks().get(0).getNDF()) > _trackChi2Cut) {
+		if(_debug)System.out.println(this.getClass().getName()+"::  failed track chi2");
+                continue;
+            }
+
+            double eleTime = TrackData.getTrackTime(TrackData.getTrackData(event, eleTop.getTracks().get(0)));
+            double posTime = TrackData.getTrackTime(TrackData.getTrackData(event, eleBot.getTracks().get(0)));
+            if (Math.abs(eleTime - posTime) > _trackDtCut) {
+		if(_debug)System.out.println(this.getClass().getName()+"::  failed track dt");
+                continue;
+            }
+            nMollers++;
+	    totalMollersPassing++; 
+        }
+
+	if (nMollers>0){
+	    incrementEventPassed();
+	    return true;
+	} else
+	    return false;
 
-	return pass; 
+    }    
+
+    public MollerSkimmer(String file) {
+	super(file, null); 
+    }
+    public MollerSkimmer(String file, Set<String> ignore) {
+	super(file, ignore); 
     }
 
-
     @Override
     public void setParameters(String parsFileName){
+	String infilePreResDir = "/org/hps/recon/skims/"; 
+	String infile=infilePreResDir+parsFileName; 
+        InputStream inParamStream = this.getClass().getResourceAsStream(infile);
+        System.out.println(this.getClass().getName()+"::  reading in Moller skimming cuts from "+infile);
+        BufferedReader reader = new BufferedReader(new InputStreamReader(inParamStream));
+        String line;
+        String delims = "[ ]+";// this will split strings between one or more spaces
+	try {
+            while ((line = reader.readLine()) != null) {
+                String[] tokens = line.split(delims);
+		String parName=tokens[0].replaceAll("\\s+","");
+                System.out.println(this.getClass().getName()+"::  parameter name = " + parName + "; value = " + tokens[1]);
+		putParam(parName,tokens[1]);
+
+            }
+        } catch (IOException ex) {
+	    System.out.println(this.getClass().getName()+":: died while reading parameters");
+            return;
+        }        
 	return;
     }
 
-    public MollerSkimmer(String file) {
-	super(file, null); 
-	//        this(super.addFileExtension(file), null);
+
+    private void putParam(String parName, String var){
+	if(parName.equals("MollerCandidateCollectionName"))
+	    _MollerCandidateCollectionName=var;
+	else if(parName.equals("vtxChi2Cut"))
+	    _vtxChi2Cut=Double.parseDouble(var);
+	else if(parName.equals("trackChi2Cut"))
+	    _trackChi2Cut=Double.parseDouble(var);
+	else if(parName.equals("trackDtCut"))
+	    _trackDtCut=Double.parseDouble(var);
+	else if(parName.equals("nHitsMin"))
+	    _nHitsMin=Integer.parseInt(var);
+	else
+	    System.out.println(this.getClass().getName()+":: couldn't find "+parName+"!");  
     }
-    public MollerSkimmer(String file, Set<String> ignore) {
-	super(file, ignore); 
+
+    public int getTotalMollersPassing(){
+	return totalMollersPassing; 
+    }
+
+    public int getTotalMollers(){
+	return totalMollers; 
+    }
+
+    public void setVtxChi2Cut(double cutVal){
+	this._vtxChi2Cut=cutVal; 
+    }
+    public void setTrackChi2Cut(double cutVal){
+	this._trackChi2Cut=cutVal; 
+    }
+    public void setTrackDtCut(double cutVal){
+	this._trackDtCut=cutVal; 
+    }
+    public void setNHitsMin(int cutVal){
+	this._nHitsMin=cutVal; 
     }
 }

recon/src/main/java/org/hps/recon/skims/MultiSkimDriver.java

@@ -51,15 +51,25 @@ public class MultiSkimDriver extends Driver {
     List<Skimmer> writeSkimList=new ArrayList<Skimmer>(); 
 
     public void endOfData() {
+
         System.out.println(this.getClass().getSimpleName() + " Summary: ");
+
 	if(skimV0){
-
 	    System.out.println("V0 skim events processed = " + v0Skimmer.getNProcessed());
-	    System.out.println("events passed            = " + v0Skimmer.getNPassed());
-	    System.out.println("       pass efficiency   = " + v0Skimmer.getPassFraction());
+	    System.out.println("Events passed            = " + v0Skimmer.getNPassed());
+	    System.out.println("       Pass efficiency   = " + v0Skimmer.getPassFraction());
 	    System.out.println("Total number of V0s processed = "+((V0Skimmer)v0Skimmer).getTotalV0s()); 
 	    System.out.println("Total number of V0s passing skim = "+((V0Skimmer)v0Skimmer).getTotalV0sPassing()); 
 	}
+
+        if(skimMoller){
+            System.out.println("Moller skim events processed = " + mollerSkimmer.getNProcessed());
+            System.out.println("Events passed            = " + mollerSkimmer.getNPassed());
+            System.out.println("       Pass efficiency   = " + mollerSkimmer.getPassFraction());
+            System.out.println("Total number of Mollers processed = "+((MollerSkimmer)mollerSkimmer).getTotalMollers());
+            System.out.println("Total number of Mollers passing skim = "+((MollerSkimmer)mollerSkimmer).getTotalMollersPassing());
+        }
+
     }
 
     @Override

recon/src/main/java/org/hps/recon/skims/V0Skimmer.java

@@ -149,7 +149,7 @@ public void setParameters(String parsFileName){
 	String infilePreResDir = "/org/hps/recon/skims/"; 
 	String infile=infilePreResDir+parsFileName; 
         InputStream inParamStream = this.getClass().getResourceAsStream(infile);
-        System.out.println(this.getClass().getName()+"::  reading in per-sensor per-phase calibs from "+infile);
+        System.out.println(this.getClass().getName()+"::  reading in V0 skimming cuts from "+infile);
         BufferedReader reader = new BufferedReader(new InputStreamReader(inParamStream));
         String line;
         String delims = "[ ]+";// this will split strings between one or more spaces

recon/src/main/resources/org/hps/recon/skims/mollerskim_parameters_2019_ver0.txt

@@ -0,0 +1,5 @@
+MollerCandidateCollectionName	  UnconstrainedMollerCandidates_KF
+trackChi2Cut			  30.0
+nHitsMin			  6
+vtxChi2Cut			  30.0
+trackDtCut			  20.0

steering-files/src/main/resources/org/hps/steering/recon/PhysicsRun2019_pass1_recon_skimmed.lcsim

@@ -43,7 +43,7 @@
         <!-- Event filtering --> 
         <driver name="StripEvent"/>
 
-        <!--<driver name="LCIOWriter"/>-->
+        <driver name="LCIOWriter"/>
 	<driver name="AidaSaveDriver"/>
         <driver name="CleanupDriver"/>
 
@@ -183,11 +183,13 @@
         </driver>  
         <driver name="StripEvent" type="org.hps.recon.skims.MultiSkimDriver">
 	    <skimV0>true</skimV0>
-	    <skimThreeBody>false</skimThreeBody>
-	    <skimFEE>false</skimFEE>
-	    <skimMoller>false</skimMoller>
 	    <v0ParamFile>v0skim_parameters_2019_ver0.txt</v0ParamFile>  
 	    <v0OutputFile>${outputFile}_v0skim.slcio</v0OutputFile>
+            <skimMoller>true</skimMoller>
+            <mollerParamFile>mollerskim_parameters_2019_ver0.txt</mollerParamFile>
+            <mollerOutputFile>${outputFile}_mollerskim.slcio</mollerOutputFile>
+            <skimThreeBody>false</skimThreeBody>
+            <skimFEE>false</skimFEE>
         </driver>