vtxreco.cxx

Go to the documentation of this file.

#include <algorithm>
#include <cfloat>
#include <iostream>
#include <string>

#include "TChain.h"
#include "TClonesArray.h"
#include "TFile.h"
#include "TString.h"

#include "St_db_Maker/St_db_Maker.h"
#include "StEmcRawMaker/StBemcTables.h"
#include "StMuDSTMaker/COMMON/StMuDst.h"
#include "StMuDSTMaker/COMMON/StMuDstMaker.h"
#include "StMuDSTMaker/COMMON/StMuEvent.h"

#include "StGenericVertexMaker/StGenericVertexFinder.h"
#include "StGenericVertexMaker/StiPPVertex/StPPVertexFinder.h"

#include "travex/ProgramOptions.h"


bool switched_file(const StMuDstMaker& maker);
void process_muDst(tvx::ProgramOptions& prgOpts);



int main(int argc, char **argv)
{
   tvx::ProgramOptions prgOpts(argc, argv);

   prgOpts.ProcessOptions();
   prgOpts.Print();

   process_muDst(prgOpts);

   return EXIT_SUCCESS;
}



void process_muDst(tvx::ProgramOptions& prgOpts)
{
   StMuDstMaker *muDstMaker = new StMuDstMaker(0, 0, "", prgOpts.PathToInputFile().c_str(), "st:MuDst.root", 1e9); // set up maker in read mode
   //                                     0, 0                        this mean read mode
   //                                           dir                    read all files in this directory
   //                                               file               bla.lis real all file in this list, if (file!="") dir is ignored
   //                                                    filter        apply filter to filenames, multiple filters are separated by ':'
   //                                                           10      maximum number of file to read

   // Specify (active) branches to read but first disable all branches
   muDstMaker->SetStatus("*", 0);

   std::vector<std::string> activeBranchNames = {
      "MuEvent",
      "GlobalTracks",
      "CovGlobTrack",
      "EmcTow",
      "BTofHit",
      "BTofHeader",
      "StStMuMcVertex",
      "StStMuMcTrack"
   };

   // Enable selected branches
   for (const auto& branchName : activeBranchNames)
      muDstMaker->SetStatus(branchName.c_str(), 1);

   TChain *muDstChain = muDstMaker->chain();

   unsigned int nentries = muDstChain->GetEntries();
   unsigned int nEventsUser = prgOpts.GetMaxEventsUser();
   unsigned int nEventsToRead = nEventsUser > 0 ? std::min(nEventsUser, nentries) : nentries;

   std::cout << nentries << " events in chain, " << nEventsToRead << " will be read." << std::endl;

   // Output
   // Create new branch
   TClonesArray  *verticesRefitted = new TClonesArray("StMuPrimaryVertex", 1000);

   TFile*   outFile      = nullptr;
   TTree*   muDstTreeOut = nullptr;
   TBranch* new_branch   = nullptr;

   TString new_name = TString("out_all.MuDst.root");
   outFile = new TFile(new_name, "RECREATE");
   //muDstTreeOut = muDstMaker->chain()->GetTree()->CloneTree(0);
   muDstTreeOut = muDstMaker->chain()->CloneTree(0);
   new_branch = muDstTreeOut->Branch("PrimaryVertices", &verticesRefitted, 65536, 99);

   StGenericVertexFinder* vertexFinder = new StPPVertexFinder();//VertexFit_t::Beamline1D);
   //StGenericVertexFinder* vertexFinder = new StPPVertexFinder(VertexFit_t::Beamline3D);


   St_db_Maker* st_db_maker = new St_db_Maker("db", "StarDb", "MySQL:StarDb", "$STAR/StarDb");
   st_db_maker->SetMaxEntryTime(20140222, 0);
   //st_db_maker->SetMaxEntryTime(20120102, 0);
   st_db_maker->Init();

   //StBemcTables* bemcTables = new StBemcTables();

   // Keep track of the number of events with 0 reconstructed verticies
   int nEventsNoRecoVertex = 0;
   int currRunNumber = 0;
   int file_id = 0;

   // Main loop over events
   for (unsigned int iEvent = 0; iEvent < nEventsToRead; iEvent++)
   {
      if ( muDstMaker->Make() ) break;

      // Create new output file
      if ( switched_file(*muDstMaker) )
      {
      //   if (muDstTreeOut) {
      //      muDstTreeOut->Write();
      //      outFile->Close();
      //      delete outFile;
      //      file_id++;
      //   }
      //
      //   TString new_name = TString("out_");
      //   new_name += file_id;
      //   new_name += ".MuDst.root";
      //   outFile = new TFile(new_name, "RECREATE");
      //   muDstTreeOut = muDstMaker->chain()->GetTree()->CloneTree();
      //   new_branch = muDstTreeOut->Branch("PrimaryVertices", &verticesRefitted, 65536, 99);
      }

      // Access all event data
      StMuDst *muDst = muDstMaker->muDst();

      // Access event header-wise information
      StMuEvent *muEvent = muDst->event();

      // Perform new run initialization such as database access and info retrieval
      if (currRunNumber != muEvent->runNumber())
      {
         currRunNumber = muEvent->runNumber();

         st_db_maker->InitRun(currRunNumber);
         vertexFinder->InitRun(currRunNumber);//, st_db_maker);
      }

      //vertexFinder->Fit(*muDst);
      //vertexFinder->result(*verticesRefitted); // container refilled

      muDstTreeOut->Fill();

      vertexFinder->Clear();
   }

   //if (muDstTreeOut) {
      muDstTreeOut->Write();
      outFile->Close();
      delete outFile;
   //}

   std::cout << "Number of events: " <<  nEventsToRead
             << ", with zero reconstructed verticies: " << nEventsNoRecoVertex << std::endl;
}


bool switched_file(const StMuDstMaker& maker)
{
   static std::string prev_file_name("");
   static std::string curr_file_name("");

   curr_file_name = maker.GetFileName();

   if (curr_file_name != prev_file_name) {
      std::cout << "Processing new file: " << curr_file_name << std::endl;
      prev_file_name = curr_file_name;

      return true;
   }

   return false;
}