coldQCDTracking/SingleTrackAN/SingleTrackAN.cc

0001 #include "SingleTrackAN.h"
0002
0003 #include <uspin/SpinDBContent.h>
0004 #include <uspin/SpinDBContentv1.h>
0005
0006 #include <calotrigger/TriggerRunInfo.h>
0007 #include <ffarawobjects/Gl1Packet.h>
0008
0009 #include <trackbase_historic/SvtxTrack.h>
0010 #include <trackbase_historic/SvtxTrackMap.h>
0011
0012 #include <fun4all/Fun4AllReturnCodes.h>
0013
0014 #include <phool/PHCompositeNode.h>
0015 #include <phool/getClass.h>
0016
0017 #include <ffaobjects/RunHeader.h>
0018 #include <ffaobjects/EventHeader.h>
0019 #include <ffarawobjects/Gl1Packet.h>
0020
0021 #include <phool/PHNodeIterator.h>
0022
0023 #include <bcolumicount/StreamingBcoInfo.h>
0024 #include <bcolumicount/StreamingLumiInfo.h>
0025
0026
0027 #include <TFile.h>
0028 #include <TH1F.h>
0029 #include <TH2F.h>
0030 #include <TString.h>
0031 #include <TTree.h>
0032 #include <TVector3.h>
0033
0034 #include <algorithm>
0035 #include <cassert>
0036 #include <cmath>
0037 #include <iostream>
0038 #include <limits>
0039 #include <stdexcept>
0040
0041 //____________________________________________________________________________..
0042 SingleTrackAN::SingleTrackAN(const std::string &name):
0043  SubsysReco(name)
0044 {
0045   // TODO: Add Verbosity switch
0046   std::cout << "SingleTrackAN::SingleTrackAN(const std::string &name) Calling ctor" << std::endl;
0047 }
0048
0049 //____________________________________________________________________________..
0050 SingleTrackAN::~SingleTrackAN()
0051 {
0052   // TODO: Add Verbosity switch
0053   std::cout << "SingleTrackAN::~SingleTrackAN() Calling dtor" << std::endl;
0054 }
0055
0056 //____________________________________________________________________________..
0057 int SingleTrackAN::Init([[maybe_unused]] PHCompositeNode *topNode)
0058 {
0059   // TODO: Add Verbosity switch
0060   std::cout << "SingleTrackAN::Init(PHCompositeNode *topNode) Initializing" << std::endl;
0061   return Fun4AllReturnCodes::EVENT_OK;
0062 }
0063
0064
0065
0066
0067
0068
0069
0070 //____________________________________________________________________________..
0071 int SingleTrackAN::InitRun([[maybe_unused]] PHCompositeNode *topNode)
0072 {
0073   // TODO: Add Verbosity switch
0074   //std::cout << "SingleTrackAN::InitRun(PHCompositeNode *topNode) Initializing for Run XXX" << std::endl;
0075
0076   outfile = new TFile(m_output_filename.c_str(), "RECREATE");
0077   outfile->cd();
0078
0079   SpinDBContent* spindb = findNode::getClass<SpinDBContent>(topNode, "SpinDBContent");
0080   if (!spindb)
0081   {
0082     std::cout << PHWHERE << ":: SpinDB node missing! Skipping run" << std::endl;
0083     return Fun4AllReturnCodes::ABORTRUN;
0084   }
0085
0086   StreamingLumiInfo *streaming_lumi_info = findNode::getClass<StreamingLumiInfo>(topNode, "STREAMINGLUMIINFO");
0087   if (streaming_lumi_info)
0088   {
0089     std::cout << " lumi : " << streaming_lumi_info->get_lumi() << std::endl;
0090   }
0091
0092   for (int i = 0; i < 120; i++)
0093   {
0094     spinPatternBlue[i] = -1*spindb->GetSpinPatternBlue(i); //-1 for IP12 -> IP8
0095     spinPatternYellow[i] = -1*spindb->GetSpinPatternYellow(i); //-1 for IP12 -> IP8
0096
0097     std::cout << spinPatternBlue[i] << " " << spinPatternYellow[i] << std::endl;
0098   }
0099
0100   hbunchnum = new TH1F("hbunchnum", "Bunch Number", 120, 0, 120);
0101   hsphnxbunchnum = new TH1F("hsphnxbunchnum", "Sphenix Bunch Number", 120, 0, 120);
0102   hmultiplicity = new TH1F("hmultiplicity", "Multiplicity", 100, 0, 100);
0103   hcrossing = new TH1F("hcrossing", "Crossing", 550, -150, 400);
0104   htrackBunch = new TH1F("htrackBunch", "Track Bunch", 120, 0, 120);
0105   hbspin = new TH1F("hbspin", "Blue Spin", 3, -1, 2);
0106   hyspin = new TH1F("hyspin", "Yellow Spin", 3, -1, 2);
0107
0108   singleTrackTree = new TTree("singleTrackTree", "Single Track Analysis Tree");
0109   singleTrackTree->SetDirectory(outfile);
0110   singleTrackTree->Branch("bunchnum", &bunchnum, "bunchnum/I");
0111   singleTrackTree->Branch("sphnxbunchnum", &sphnxbunchnum, "sphnxbunchnum/I");
0112   singleTrackTree->Branch("multiplicity", &multiplicity, "multiplicity/I");
0113   singleTrackTree->Branch("crossing", &crossing, "crossing/I");
0114   singleTrackTree->Branch("trackBunch", &trackBunch, "trackBunch/I");
0115   singleTrackTree->Branch("bspin", &bspin, "bspin/I");
0116   singleTrackTree->Branch("yspin", &yspin, "yspin/I");
0117   singleTrackTree->Branch("bco", &bco, "bco/l");
0118   singleTrackTree->Branch("bcowindowlow", &bcowindowlow, "bcowindowlow/l");
0119   singleTrackTree->Branch("bcowindowhigh", &bcowindowhigh, "bcowindowhigh/l");
0120   singleTrackTree->Branch("usable_bco_tag", &usable_bco_tag, "usable_bco_tag/O");
0121   singleTrackTree->Branch("doublecount", &doublecount, "doublecount/O");
0122   singleTrackTree->Branch("p", &p, "p/F");
0123   singleTrackTree->Branch("pt", &pt, "pt/F");
0124   singleTrackTree->Branch("pz", &pz, "pz/F");
0125   singleTrackTree->Branch("xf", &xf, "xf/F");
0126   singleTrackTree->Branch("eta", &eta, "eta/F");
0127   singleTrackTree->Branch("phi", &phi, "phi/F");
0128   singleTrackTree->Branch("quality", &quality, "quality/F");
0129   singleTrackTree->Branch("charge", &charge, "charge/F");
0130
0131   runTree = new TTree("runTree", "Tree for quantities needed only once per run");
0132   runTree->SetDirectory(outfile);
0133   runTree->Branch("bpol", &bpol, "bpol/F");
0134   runTree->Branch("bpolerr", &bpolerr, "bpolerr/F");
0135   runTree->Branch("ypol", &ypol, "ypol/F");
0136   runTree->Branch("ypolerr", &ypolerr, "ypolerr/F");
0137   runTree->Branch("lumi", &lumi, "lumi/D");
0138   runTree->Branch("bunchnumber_lumi", &bunchnumber_lumi);
0139
0140   // Double check that this is working as expected... It seems like the sphnxbunchnum and bunchnum distributions are identical, which makes me think xingshift is 0
0141   xingshift = spindb->GetCrossingShift();
0142
0143   spindb->GetPolarizationBlue(0,bpol,bpolerr);
0144   spindb->GetPolarizationYellow(0,ypol,ypolerr);
0145
0146   lumi = streaming_lumi_info->get_lumi();
0147   bunchnumber_lumi = streaming_lumi_info->get_bunchnumber_lumi();
0148
0149   runTree->Fill();
0150
0151   return Fun4AllReturnCodes::EVENT_OK;
0152 }
0153
0154 //____________________________________________________________________________..
0155 int SingleTrackAN::process_event([[maybe_unused]] PHCompositeNode *topNode)
0156 {
0157
0158
0159   Gl1Packet* gl1packet = findNode::getClass<Gl1Packet>(topNode, "GL1RAWHIT");
0160   bunchnum = gl1packet->getBunchNumber(); // this will get me triggered crossing bunch number
0161   sphnxbunchnum = (bunchnum + xingshift) % 120;
0162
0163   StreamingBcoInfo* bcoinfo = findNode::getClass<StreamingBcoInfo>(topNode, "STREAMINGBCOINFO");
0164   bco = bcoinfo->get_bco();
0165   usable_bco_tag = bcoinfo->get_usable_bco_tag();
0166   auto window = bcoinfo->get_bco_streaming_window();
0167   bcowindowlow = window.first;
0168   bcowindowhigh = window.second;
0169
0170
0171   // Check BCO output -- TODO: add Verbosity switch //
0172
0173   //std::cout << "BCO: " << bco << std::endl;
0174   //std::cout << "Usable BCO Tag: " << usable_bco_tag << std::endl;
0175   //std::cout << "BCO Streaming Window: " << bcowindowlow << " to " << bcowindowhigh << std::endl;
0176
0177
0178   SvtxTrackMap* trackmap = findNode::getClass<SvtxTrackMap>(topNode, "SvtxTrackMap");
0179   if (!trackmap)
0180   {
0181     trackmap = findNode::getClass<SvtxTrackMap>(topNode, "TrackMap");
0182     if (!trackmap)
0183     {
0184       std::cout
0185           << "TrackSpin::process_event - Error can not find DST trackmap node SvtxTrackMap" << std::endl;
0186       exit(-1);
0187     }
0188   }
0189   multiplicity = trackmap->size();
0190   for (SvtxTrackMap::Iter iter = trackmap->begin(); iter != trackmap->end(); ++iter)
0191   {
0192
0193     SvtxTrack* track = iter->second;
0194     crossing = track->get_crossing();
0195     if (crossing > int(bcowindowhigh - bco))
0196     {
0197       doublecount = true;
0198     }
0199     else
0200     {
0201       doublecount = false;
0202     }
0203     charge = track->get_charge();
0204     quality = track->get_quality();
0205     p = track->get_p();
0206     pt = track->get_pt();
0207     pz = track->get_pz();
0208     xf = (2*pz)/200;
0209     eta = track->get_eta();
0210     phi = track->get_phi();
0211
0212     // It seems like sphnxbunchnum and bunchnum have identical distributions. Is this expected?
0213     trackBunch = crossing + sphnxbunchnum;
0214     while (trackBunch < 0)
0215     {
0216       trackBunch += NBUNCHES;
0217     }
0218     trackBunch %= NBUNCHES;
0219
0220     bspin = spinPatternBlue[trackBunch];
0221     yspin = spinPatternYellow[trackBunch];
0222
0223
0224     singleTrackTree->Fill();
0225
0226   }
0227
0228     // TODO: Add Verbosity switch
0229     //std::cout << "SingleTrackAN::process_event(PHCompositeNode *topNode) Processing Event" << std::endl;
0230     return Fun4AllReturnCodes::EVENT_OK;
0231 }
0232
0233 //____________________________________________________________________________..
0234 int SingleTrackAN::ResetEvent([[maybe_unused]] PHCompositeNode *topNode)
0235 {
0236   // TODO: Add Verbosity switch
0237   //std::cout << "SingleTrackAN::ResetEvent(PHCompositeNode *topNode) Resetting internal structures, prepare for next event" << std::endl;
0238   return Fun4AllReturnCodes::EVENT_OK;
0239 }
0240
0241 //____________________________________________________________________________..
0242 int SingleTrackAN::EndRun(const int runnumber)
0243 {
0244   // TODO: Add Verbosity switch
0245   std::cout << "SingleTrackAN::EndRun(const int runnumber) Ending Run for Run " << runnumber << std::endl;
0246   return Fun4AllReturnCodes::EVENT_OK;
0247 }
0248
0249 //____________________________________________________________________________..
0250 int SingleTrackAN::End([[maybe_unused]] PHCompositeNode *topNode)
0251 {
0252   // TODO: Add Verbosity switch
0253   std::cout << "SingleTrackAN::End(PHCompositeNode *topNode) This is the End..." << std::endl;
0254
0255
0256   outfile->cd();
0257   singleTrackTree->Write();
0258   runTree->Write();
0259   outfile->Close();
0260
0261
0262
0263   return Fun4AllReturnCodes::EVENT_OK;
0264 }
0265
0266 //____________________________________________________________________________..
0267 int SingleTrackAN::Reset([[maybe_unused]] PHCompositeNode *topNode)
0268 {
0269  std::cout << "SingleTrackAN::Reset(PHCompositeNode *topNode) being Reset" << std::endl;
0270   return Fun4AllReturnCodes::EVENT_OK;
0271 }
0272
0273 //____________________________________________________________________________..
0274 void SingleTrackAN::Print(const std::string &what) const
0275 {
0276   std::cout << "SingleTrackAN::Print(const std::string &what) const Printing info for " << what << std::endl;
0277 }
0278
0279
0280