hps-java/src/main/java/org/lcsim/hps/util
diff -u -r1.1 -r1.2
--- CalculateAcceptanceFromMadGraph.java 5 Sep 2012 00:17:54 -0000 1.1
+++ CalculateAcceptanceFromMadGraph.java 5 Sep 2012 19:41:04 -0000 1.2
@@ -50,7 +50,7 @@
static double sigz = 0.02; // Luminous region size in z (mm)
static double maxLen = 200; // maximum decay length (mm)
static double _ecm = 2.2; //GeV
- static int nInFiles = 1; //number of input files
+ static int nInFiles = 1000; //number of input files
static int nmax = 500000; //maximum number of events to write to 1 stdhep file (new one opens if n>nmax)
static String fileType = "lhe";
static int nread = -1; // a running total of number of events read/written to stdhep files
@@ -78,7 +78,7 @@
static int[] recoMassLyr1 = new int[nbins];
static int[] recoMassLyr50 = new int[nbins];
- public static void doAccounting(List<Boolean> passEle, List<Boolean> passPos, double invMass) {
+ public static void doAccounting(List<Boolean> passEle, List<Boolean> passPos, boolean passEvt, double invMass) {
boolean passLyr700 = false;
boolean passLyr100 = false;
boolean passLyr500 = false;
@@ -86,7 +86,9 @@
//find the bin...
int bin = (int) (invMass / maxMass * nbins);
// System.out.println("invMass = "+invMass + "...goes in bin #"+bin);
+ if(bin<nbins){
genMass[bin]++;
+
//see if pos+ele is in layer 70cm layer (required form both test & full)
if (passEle.get(5) && passPos.get(5))
passLyr700 = true;
@@ -105,7 +107,10 @@
if (passLyr700 && passLyr50) //require they go through 6 layers...
recoMassLyr50[bin]++;
-
+ }else{
+ System.out.println("Mass out of range! "+invMass);
+ }
+
}
private static Options createCommandLineOptions() {
@@ -113,6 +118,7 @@
options.addOption(new Option("e", true, "Beam Energy (GeV)"));
options.addOption(new Option("n", true, "Number of files to run."));
+ options.addOption(new Option("b", true, "B-Field"));
return options;
}
@@ -152,7 +158,7 @@
if (cl.hasOption("e")) {
ecmString = cl.getOptionValue("e");
- _ecm = Double.valueOf(ninString);
+ _ecm = Double.valueOf(ecmString);
System.out.println(ecmString);
eptString = convertDecimal(ecmString);
}
@@ -170,8 +176,7 @@
// String outDir = "/nfs/slac/g/hps/mgraham/DarkPhoton/SignalEvents/";
String outDir = "./";
-// for (int i = 0; i < nInFiles; i++) {
- for (int i = 933; i < 933 + nInFiles; i++) {
+ for (int i = 0; i < nInFiles; i++) {
int fnum = i + 1;
String snum = "_" + fnum;
if (fnum < 10000) {
@@ -517,9 +522,9 @@
double xc = getxc(R, d0, phi0);
double yc = getyc(R, d0, phi0);;
-// System.out.println(p.toString());
-// System.out.println("d0 = "+d0+"; phi0 = "+phi0+"; z0 = "+z0+"; slope = "+slope+"; R = "+R);
-// System.out.println("x0 = "+x0+"; y0 = "+y0+"; xc = "+xc+"; yc = "+yc);
+ // System.out.println(p.toString());
+ // System.out.println("d0 = "+d0+"; phi0 = "+phi0+"; z0 = "+z0+"; slope = "+slope+"; R = "+R);
+ // System.out.println("x0 = "+x0+"; y0 = "+y0+"; xc = "+xc+"; yc = "+yc);
List<Boolean> passLayer = new ArrayList<Boolean>();
for (int ii = 0; ii < nLayers; ii++) {
double pathL = PathToXPlane(x0, y0, xc, yc, R, x[ii]);
@@ -541,8 +546,8 @@
i++;
}
}
- if (eventPass(pEle, pPos))
- doAccounting(passEle, passPos, 1000.0 * getInvMass(pEle, pPos));
+
+ doAccounting(passEle, passPos, eventPass(pEle,pPos),1000.0 * getInvMass(pEle, pPos));
}
@@ -571,13 +576,18 @@
}
public static boolean eventPass(Hep3Vector p1, Hep3Vector p2) {
-
+ // System.out.println("p1.magnitude = "+p1.magnitude()+"; p2.magnitude = "+p2.magnitude()+"; 0.8*_ecm = "+0.8*_ecm);
+
if (p1.magnitude() + p2.magnitude() < 0.8 * _ecm)//trigger requires 80% of beam energy
return false;
+// System.out.println("Passed totenergy");
if (p1.magnitude() < _pCut)
return false;
if (p2.magnitude() < _pCut)
return false;
+ if(p2.z()*p1.z()>0) // this is basically the opposite quadrant cut in the trigger (B-field makes them opposite in y)
+ return false;
+ // System.out.println("Event is good!!!!!");
return true;
}
hps-java/src/main/java/org/lcsim/hps/util
diff -u -r1.5 -r1.6
--- ConvertToStdhep.java 29 Aug 2012 15:40:46 -0000 1.5
+++ ConvertToStdhep.java 5 Sep 2012 19:41:04 -0000 1.6
@@ -44,7 +44,7 @@
static double _yoff = 0.0;
static double _zoff = 0.03;
static double aMass = 99; //Aprime mass (MeV)
- static double sigx = 0.02; // Luminous region size in x (mm)
+ static double sigx = 0.2; // Luminous region size in x (mm)
static double sigy = 0.02; // Luminous region size in y (mm)
static double sigz = 0.0; // Luminous region size in z (mm)
//beam is positioned so that at first beam direction is in z, then rotated to correct orientation
@@ -60,16 +60,18 @@
static IRotation3D rot = new RotationGeant(rotx, roty, rotz);
// static String fileType="dat";
static int nread = -1; // a running total of number of events read/written to stdhep files
- static boolean _eventFilter=false;
- static int _nEleRequired=2;
-
+ static boolean _eventFilter = false;
+ static int _nEleRequired = 2;
+
private static Options createCommandLineOptions() {
Options options = new Options();
options.addOption(new Option("m", true, "A' Mass (MeV)"));
options.addOption(new Option("e", true, "Beam Energy (GeV)"));
options.addOption(new Option("n", true, "Number of files to run."));
- options.addOption(new Option("s", false, "Filter Events"));
+ options.addOption(new Option("x", true, "Beam sigma in x"));
+ options.addOption(new Option("y", true, "Beam sigma in y"));
+ options.addOption(new Option("s", false, "Filter Events"));
return options;
}
@@ -98,6 +100,8 @@
String ninString = String.valueOf(nInFiles);
String ecmString = String.valueOf(_ecm);
String eptString = String.valueOf(_ecm);
+ String sigxString = String.valueOf(sigx);
+ String sigyString = String.valueOf(sigy);
eptString = convertDecimal(eptString);
// LCIO output file.
if (cl.hasOption("m")) {
@@ -111,6 +115,17 @@
nInFiles = Integer.valueOf(ninString);
System.out.println(ninString);
}
+
+ if (cl.hasOption("x")) {
+ sigxString = cl.getOptionValue("x");
+ sigx = Double.valueOf(sigxString);
+ System.out.println(sigxString);
+ }
+ if (cl.hasOption("y")) {
+ sigyString = cl.getOptionValue("y");
+ sigy = Double.valueOf(sigyString);
+ System.out.println(sigyString);
+ }
if (cl.hasOption("e")) {
ecmString = cl.getOptionValue("e");
@@ -118,14 +133,17 @@
System.out.println(ecmString);
eptString = convertDecimal(ecmString);
}
-
- String filter="all";
+
+ String filter = "all";
if (cl.hasOption("s")) {
- _eventFilter=true;
- filter="selected";
+ _eventFilter = true;
+ filter = "selected";
}
+ sigxString=convertMicron(sigx);
+ sigyString=convertMicron(sigy);
- String postfix = "_20ux200u_beamspot_gammactau_0cm.stdhep";
+// String postfix = "_20ux200u_beamspot_gammactau_0cm.stdhep";
+ String postfix = "_"+sigxString+"x"+sigyString+"_beamspot_gammactau_0cm.stdhep";
// String fDir="/nfs/slac/g/hps/mgraham/DarkPhoton/tvm/testrun/";
// String fileLabel = "ap2.2gev40mevsel";
@@ -137,7 +155,7 @@
//String inLabel = "W2pt2GeV_Ap100MeV_";
String fDir = "/nfs/slac/g/hps/mgraham/DarkPhoton/MadGraph/Events" + eptString + "Ap" + massString + "MeV/";
- String fileLabel = "ap" + ecmString + "gev" + massString + "mev"+filter;
+ String fileLabel = "ap" + ecmString + "gev" + massString + "mev" + filter;
String inLabel = "W" + eptString + "GeV_Ap" + massString + "MeV_";
String inPost = "_unweighted_events.lhe";
@@ -584,7 +602,7 @@
double ApMass = 0;
double ApEnergy = 0;
boolean acceptEvent = false;
- int nElePass = 0;
+ int nElePass = 0;
for (int npart = 0; npart < nhep; npart++) {
List<Double> vals = getNumbersInLine(tok);
if (vals.size() != 13) {
@@ -683,11 +701,12 @@
}
}
StdhepEvent ev = new StdhepEvent(nevhep, nhep, isthep, idhep, jmohep, jdahep, phep, vhep);
- if(nElePass>=_nEleRequired)
- acceptEvent=true;
+ if (nElePass >= _nEleRequired) {
+ acceptEvent = true;
+ }
if (acceptEvent) {
sw.writeRecord(ev);
- } else{
+ } else {
// System.out.println("LHE Event failed acceptance cuts");
}
}
@@ -704,19 +723,25 @@
}
return num;
}
+
+ public static String convertMicron(double num) {
+ double mic=num*1000.0;
+ String out=Integer.toString((int)mic);
+ return out+"u";
+ }
public static boolean inAcceptance(Hep3Vector ph) {
boolean ok = false;
double[] p = {ph.x(), ph.y(), ph.z()};
double ptot = Math.sqrt(p[0] * p[0] + p[1] * p[1] + p[2] * p[2]);
-
+
double sinThx = p[0] / Math.sqrt(p[0] * p[0] + p[2] * p[2]);
double sinThy = p[1] / Math.sqrt(p[1] * p[1] + p[2] * p[2]);
// for now, just use thetay (non-bend direction)
// System.out.println("px = "+p[0]+"; py = "+p[1]+"; pz = "+p[2]);
// System.out.println(sinThy+" "+sinThx+" "+ptot);
-
- if (Math.abs(sinThy) > 0.012&&ptot>0.1) {
+
+ if (Math.abs(sinThy) > 0.012 && ptot > 0.1) {
ok = true;
}
return ok;