macros/leptoquarks_analysis/MachineLearning_CSV.C

0001 #include "TROOT.h"
0002 #include "TClass.h"
0003 #include "TGraph.h"
0004 #include "TF1.h"
0005 #include "TH1F.h"
0006 #include "TH2F.h"
0007 #include "TH3F.h"
0008 #include "TH3D.h"
0009 #include "TTree.h"
0010 #include "TCanvas.h"
0011 #include "TBranch.h"
0012 #include "Riostream.h"
0013 #include "TStyle.h"
0014 #include "TFile.h"
0015 #include "TString.h"
0016 #include "TLegend.h"
0017 #include "TRandom3.h"
0018 #include "TMath.h"
0019 #include "math.h"
0020 #include "TColor.h"
0021 #include <vector>
0022 #include <sstream>
0023 #include <algorithm>
0024
0025 #include <cstdlib>
0026 #include "TMath.h"
0027 #include <iostream>
0028 #include <fstream>
0029 #include <string>
0030 #include <math.h>
0031 #include <cmath>
0032 #include "TGraph.h"
0033 #include "TGraph2D.h"
0034 #include <algorithm>
0035
0036
0037 /*
0038    Write out jet variables to csv file
0039    to be used for machine learning
0040
0041    written by sean.jeffas@stonybrook.edu
0042 */
0043
0044 int MachineLearning_CSV()
0045 {
0046
0047
0048   ofstream myfile;
0049
0050   string filename;
0051
0052   // Names to loop over different types of files
0053   string seed[10] = {"1","2","3","4","5","6","7","8","9","10"};
0054   string type[2] = {"3pion","SM"};
0055
0056   filename = "./data/JetSummary_p250_e20_1000events_r05.csv";
0057   myfile.open(filename.c_str());
0058
0059   // Loop over all LQ, NC, and CC geant files
0060   for(int a = 0; a<10; a++){
0061     for(int b=0; b<2; b++){
0062
0063       const std::string inFile = "LeptoAna_p250_e20_1000events_"+seed[a]+"seed_"+type[b]+"_r05.root";
0064       const std::string inDirectory = "/gpfs/mnt/gpfs02/phenix/scratch/spjeffas/data/";
0065       std::string inputFile = inDirectory+inFile;
0066
0067       TFile *f = TFile::Open(inputFile.c_str());
0068       TTree *t = (TTree*)f->Get("event");
0069
0070       const int Nevent = t->GetEntries();
0071
0072       // Variables for different jet characteristics
0073       vector<float> * tracks_rmax;
0074       vector<float> * tracks_count;
0075       vector<float> * tracks_chargesum;
0076       vector<float> * tracks_vertex;
0077       vector<float> * jetshape_radius;
0078       vector<float> * jetshape_econe_1;
0079       vector<float> * jetshape_econe_2;
0080       vector<float> * jetshape_econe_5;
0081       vector<float> * jet_eta;
0082       vector<float> * jet_minv;
0083       vector<float> * jet_etotal;
0084       vector<float> * jet_ptrans;
0085       vector<int> * evtgen_pid;
0086
0087       //point to variables in tree
0088       t->SetBranchAddress("tracks_rmax_R",&tracks_rmax);
0089       t->SetBranchAddress("tracks_count_R",&tracks_count);
0090       t->SetBranchAddress("tracks_chargesum_R",&tracks_chargesum);
0091       t->SetBranchAddress("tracks_vertex",&tracks_vertex);
0092       t->SetBranchAddress("jetshape_radius",&jetshape_radius);
0093       t->SetBranchAddress("jetshape_econe_r01",&jetshape_econe_1);
0094       t->SetBranchAddress("jetshape_econe_r02",&jetshape_econe_2);
0095       t->SetBranchAddress("jetshape_econe_r05",&jetshape_econe_5);
0096       t->SetBranchAddress("jet_eta",&jet_eta);
0097       t->SetBranchAddress("jet_minv",&jet_minv);
0098       t->SetBranchAddress("jet_etotal",&jet_etotal);
0099       t->SetBranchAddress("evtgen_pid",&evtgen_pid);
0100       t->SetBranchAddress("jet_ptrans",&jet_ptrans);
0101
0102
0103       //loop over all events
0104       for(int i = 0; i < Nevent; i++)
0105     {
0106       //Get entry for each event
0107       t->GetEntry(i);
0108
0109       for(int l=0; l < tracks_rmax->size(); l++){
0110
0111         // Get variables for each jet in event
0112         double rmax = (*tracks_rmax)[l];
0113         int count = (*tracks_count)[l];
0114         int chargesum = (*tracks_chargesum)[l];
0115         double vertex = (*tracks_vertex)[l];
0116         double radius = (*jetshape_radius)[l];
0117         double econe_1 = (*jetshape_econe_1)[l];
0118         double econe_2 = (*jetshape_econe_2)[l];
0119         double econe_5 = (*jetshape_econe_5)[l];
0120         double eta = (*jet_eta)[l];
0121         double minv = (*jet_minv)[l];
0122         double etotal = (*jet_etotal)[l];
0123         double ptrans = (*jet_ptrans)[l];
0124         int pid = (*evtgen_pid)[l];
0125
0126         if(ptrans < 5) continue;
0127
0128         // If LQ then tag as LQ in csv file
0129         if(b == 0 && pid == 15 && vertex == vertex) myfile << count << "," << chargesum << "," << eta << "," << vertex << "," << "tau" << endl;
0130
0131         //If SM then tag as SM in csv file
0132         if(b != 0 && pid != 15 && pid != 11 && vertex == vertex) myfile << count << "," << chargesum << "," << eta << "," << vertex << "," << "DIS" << endl;
0133       }
0134     }
0135     }
0136   }
0137
0138   myfile.close();
0139
0140
0141
0142
0143   return 0;
0144 }