StRoot: StTpcRSMaker/TpcRS.C Source File

00001 #ifndef __CINT__
00002 #include "Riostream.h"
00003 #include "TROOT.h"
00004 #include "TSystem.h"
00005 #include "TH2.h"
00006 #include "TH3.h"
00007 #include "TF1.h"
00008 #include "TProfile.h"
00009 #include "TString.h"
00010 #include "TTree.h"
00011 #include "TChain.h"
00012 #include "TCanvas.h"
00013 #include "TClassTable.h"
00014 #include "TFileSet.h"
00015 #include "TDataSetIter.h"
00016 #include "StBFChain.h"
00017 void bfc (const Int_t Last, 
00018           const Char_t *Chain,
00019           const Char_t *infile, 
00020           const Char_t *outfile, 
00021           const Char_t *TreeFile);
00022 //R__EXTERN StBFChain *chain;
00023 #else
00024 #define SETBIT(n,i)  ((n) |=  (1 << i))
00025 #define CLRBIT(n,i)  ((n) &= ~(1 << i))
00026 class StBFChain;
00027 StBFChain *chain;
00028 class St_db_Maker;
00029 St_db_Maker *dbMk = 0;
00030 #endif
00031 //________________________________________________________________________________
00032 void TpcRS(Int_t First, Int_t Last, const Char_t *Run = "y2011,TpcRS",  
00033            const Char_t *fileIn = 0, const Char_t *opt = "Bichsel", const Char_t *kuip = 0) {
00034   gROOT->LoadMacro("bfc.C"); 
00035   TString ChainOpt("");
00036   TString RootFile("");
00037   TString Opt(opt);
00038   TString RunOpt(Run);
00039   RunOpt.ToLower();
00040   //  ChainOpt = "MakeEvent,ITTF,ForceGeometry,NoSsdIt,NoSvtIt,Idst,VFMinuit,analysis,dEdxY2,";
00041   ChainOpt = "MakeEvent,ITTF,NoSsdIt,NoSvtIt,Idst,VFMinuit,analysis,dEdxY2,";
00042   ChainOpt += "Corr4";// no dynamical distortion ! ,OSpaceZ2,OGridLeak3D,"; // check that StTpcRSMaker::kDistortion bit is set
00043   //  ChainOpt += "EvOut,MuDST,MiniMcMk,McTpcAna,IdTruth,useInTracker,-hitfilt,";
00044   ChainOpt += ",CMuDst,MiniMcMk,IdTruth,useInTracker,tree,";
00045   if (TString(gSystem->Getenv("STAR_VERSION")) == ".DEV2") ChainOpt += "McTpcAna,NoHistos,NoRunco,";
00046   else                                                     ChainOpt += "tags,";
00047   // ChainOpt += "MiniMcMk,IdTruth,useInTracker,-hitfilt,CMuDst,Tree,tags,evout,";
00048   if (RunOpt.Contains("fcf",TString::kIgnoreCase)) {
00049     ChainOpt += "tpl,tpcI,";
00050     RunOpt.ReplaceAll("TpcRS,","");
00051     RunOpt.ReplaceAll("trs,","");
00052   } else {
00053     ChainOpt += "tpcDB,TpcHitMover,TpxClu,";
00054   }
00055   //  Bool_t needAlias = kFALSE;
00056   TString FileIn(fileIn);
00057   if (FileIn == "") {
00058     ChainOpt += "gstar,"; RootFile += "gstar";
00059     if (RunOpt.Contains("hadr_of",TString::kIgnoreCase) ||
00060            Opt.Contains("hadr_of",TString::kIgnoreCase)) ChainOpt += "hadr_off,";
00061     if (RunOpt.Contains("phys_of",TString::kIgnoreCase) ||
00062            Opt.Contains("phys_of",TString::kIgnoreCase)) ChainOpt += "phys_off,";
00063     if (RunOpt.Contains("PhysicsOff",TString::kIgnoreCase) ||
00064            Opt.Contains("PhysicsOff",TString::kIgnoreCase)) ChainOpt += "phys_off,";
00065     if (! RunOpt.Contains("Y200",TString::kIgnoreCase) &&
00066         ! RunOpt.Contains("dev",TString::kIgnoreCase)
00067         ) ChainOpt += "Y2011,";
00068     if      (Opt.Contains("FieldOff" ,TString::kIgnoreCase)) ChainOpt += "FieldOff,";
00069     else if (Opt.Contains("HalfField",TString::kIgnoreCase)) ChainOpt += "HalfField,";
00070     else                                                     ChainOpt += "FieldOn,";
00071   } else {
00072     RootFile = Form("%s",gSystem->BaseName(FileIn.Data())); 
00073     if (FileIn.Contains(".daq",TString::kIgnoreCase)) {
00074       ChainOpt += "in,TpxRaw,";
00075       RootFile.ReplaceAll(".daq","");
00076     } else {
00077       if (FileIn.Contains(".fz",TString::kIgnoreCase)) {
00078         ChainOpt += "fzin,";
00079         RootFile.ReplaceAll(".fzd","");
00080         RootFile.ReplaceAll(".fz","");
00081       } else {
00082         if (FileIn.Contains(".geant.root",TString::kIgnoreCase)) {
00083           ChainOpt += "in,";
00084           RootFile.ReplaceAll(".geant.root","");
00085         } else {
00086           if (FileIn.Contains(".MuDst",TString::kIgnoreCase)) {
00087             ChainOpt += "mtin,";
00088             RootFile.ReplaceAll(".MuDst.root","");
00089           }
00090         }
00091       }
00092     }
00093   }
00094   ChainOpt += RunOpt;
00095   RootFile += Form("_%s_%i_%i",Opt.Data(),First,Last);
00096   RootFile.ReplaceAll(",","_");
00097   if (RootFile.Contains(";")) {
00098     Int_t index = RootFile.Index(";");
00099     RootFile = RootFile(0,index);
00100   }
00101   RootFile += ".root";
00102   RootFile.ReplaceAll(" ","");
00103   cout << "ChainOpt : " << ChainOpt.Data() << "\tOuput file " << RootFile.Data() << endl;
00104   
00105 
00106   TString output = RootFile;
00107   output.ReplaceAll(".root","O.root");
00108   output.ReplaceAll("*","");
00109   if (Last < 0) {
00110     bfc(-1,ChainOpt.Data(),0,0,0);
00111     return;
00112   }
00113   bfc(-1,ChainOpt.Data(),fileIn,output.Data(),RootFile.Data());
00114   if (ChainOpt.Contains("TpcRS",TString::kIgnoreCase)) {
00115     StTpcRSMaker *tpcRS = (StTpcRSMaker *) chain->Maker("TpcRS");
00116     if (tpcRS) {
00117       //      if (needAlias) tpcRS->SetInput("geant","bfc/.make/inputStream/.make/inputStream_Root/.data/bfcTree/geantBranch");
00118       Int_t m_Mode = tpcRS->GetMode();
00119       if (Opt.Contains("pai",TString::kIgnoreCase))         {SETBIT(m_Mode,StTpcRSMaker::kPAI); CLRBIT(m_Mode,StTpcRSMaker::kBICHSEL);}
00120       if (Opt.Contains("bichsel",TString::kIgnoreCase))     {SETBIT(m_Mode,StTpcRSMaker::kBICHSEL); CLRBIT(m_Mode,StTpcRSMaker::kPAI);}
00121       if (! ChainOpt.Contains("Corr",TString::kIgnoreCase)) {CLRBIT(m_Mode,StTpcRSMaker::kDistortion);}  // Check that distorton are IN chain
00122       tpcRS->SetMode(m_Mode);
00123       //      tpcRS->SetDebug(13);
00124     }
00125   }
00126   else {
00127     if (ChainOpt.Contains("trs",TString::kIgnoreCase)) {
00128       StMaker *mk = chain->Maker("tpcDB"); // simulation mode
00129       if (! mk) return;
00130       mk->SetMode(1); 
00131     }
00132   }
00133   StMaker *dEdxY2 = chain->GetMaker("dEdxY2"); 
00134   if (dEdxY2) {
00135     StdEdxY2Maker *dEdx = (StdEdxY2Maker *) dEdxY2;
00136     Int_t mask = 0;
00137 #if 0
00138     //     SETBIT(mask,StTpcdEdxCorrection::ktpcPressure); 
00139     //     SETBIT(mask,StTpcdEdxCorrection::kAdcCorrection); 
00140     //     SETBIT(mask,StTpcdEdxCorrection::kTpcSecRow); 
00141     //     SETBIT(mask,StTpcdEdxCorrection::kDrift);
00142     //     SETBIT(mask,StTpcdEdxCorrection::kzCorrection);
00143     //     SETBIT(mask,StTpcdEdxCorrection::kdXCorrection);
00144     //     SETBIT(mask,StTpcdEdxCorrection::kTpcdEdxCor);
00145     //     SETBIT(mask,StTpcdEdxCorrection::kTpcLengthCorrection);
00146     // from dEdx 06/20/04
00147     SETBIT(mask,StTpcdEdxCorrection::ktpcPressure); 
00148     //  SETBIT(mask,StTpcdEdxCorrection::ktpcMethaneIn); 
00149     //  SETBIT(mask,StTpcdEdxCorrection::ktpcGasTemperature); 
00150     //  SETBIT(mask,StTpcdEdxCorrection::ktpcWaterOut); 
00151     SETBIT(mask,StTpcdEdxCorrection::kAdcCorrection); 
00152     SETBIT(mask,StTpcdEdxCorrection::kTpcSecRow); 
00153     SETBIT(mask,StTpcdEdxCorrection::kDrift);
00154     SETBIT(mask,StTpcdEdxCorrection::kzCorrection);
00155     SETBIT(mask,StTpcdEdxCorrection::kdXCorrection);
00156     //  SETBIT(mask,StTpcdEdxCorrection::kTpcdEdxCor);
00157     //  SETBIT(mask,StTpcdEdxCorrection::kTpcLengthCorrection);
00158     SETBIT(mask,StTpcdEdxCorrection::kAdcCorrection);
00159     SETBIT(mask,StTpcdEdxCorrection::kTpcLast);
00160     //    SETBIT(Mode,StdEdxY2Maker::kOldClusterFinder); 
00161     //    SETBIT(Mode,StdEdxY2Maker::kDoNotCorrectdEdx);
00162 #endif    
00163     Int_t Mode = 2;
00164     
00165     SETBIT(Mode,StdEdxY2Maker::kPadSelection); 
00166     SETBIT(Mode,StdEdxY2Maker::kCalibration);
00167     if (TString(gSystem->Getenv("STAR_VERSION")) == ".DEV2") 
00168       SETBIT(Mode,StdEdxY2Maker::kZBGX);
00169     SETBIT(Mode,StdEdxY2Maker::kGASHISTOGRAMS);
00170     if (Mode) {
00171       cout << " set dEdxY2 Mode" << Mode << " =======================================" << endl;
00172       dEdx->SetMode(Mode); 
00173     }
00174     if (mask) {
00175       cout << " set dEdxY2 mask " << mask << " =======================================" << endl;
00176       dEdx->SetMask(mask); 
00177     }
00178   }
00179   Int_t initStat = chain->Init(); // This should call the Init() method in ALL makers
00180   if (initStat) {
00181     cout << "Chain initiation has failed" << endl;
00182     chain->Fatal(initStat, "during Init()");
00183   }
00184   if (FileIn == "" && gClassTable->GetID("TGiant3") >= 0) {
00185     St_geant_Maker *geant = (St_geant_Maker *) chain->GetMakerInheritsFrom("St_geant_Maker");
00186     if (Opt.Contains("debug",TString::kIgnoreCase)) {
00187       geant->Do("debug on");
00188       geant->SetDebug(1);
00189       geant->Do("swit 1 2");
00190       geant->Do("swit 2 2");
00191     }
00192     //                   NTRACK  ID PTLOW PTHIGH YLOW YHIGH PHILOW PHIHIGH ZLOW ZHIGH
00193     //    geant->Do("gkine 100  14   0.1    10.  -1     1      0    6.28    0.    0.;");
00194     cout << Opt << endl;
00195     if (Opt.Contains(":")) {
00196       Int_t i = (gSystem->GetPid()+ First)%256;
00197       geant->Do(Form("rndm %i; rndm;",i));
00198     }
00199     if ( Opt.Contains("laser",TString::kIgnoreCase)) {
00200       gSystem->Load("gstar.so");
00201       geant->Do("call gstar");
00202       geant->Do("gkine 1 170   1   1  0   0   0  0    180.00    180.00;");
00203       geant->Do("gprint kine");
00204       geant->Do("gvert 0  54   0");
00205       geant->Do("mode TRAC prin 15");
00206     } else {
00207       Int_t    NTRACK = 100;
00208       Int_t    ID = 5;
00209       Double_t mass = 0.1057;
00210       if      (Opt.Contains("muon",TString::kIgnoreCase))     {ID =  5;                 if (Opt.Contains("muon-",TString::kIgnoreCase)) ID = 6;}
00211       else if (Opt.Contains("electron",TString::kIgnoreCase)) {ID =  2; mass = 0.5110E-03;}
00212       else if (Opt.Contains("positron",TString::kIgnoreCase)) {ID =  3; mass = 0.5110E-03;}
00213       else if (Opt.Contains("pion",TString::kIgnoreCase))     {ID =  8; mass = 0.1396; if (Opt.Contains("pion-",TString::kIgnoreCase)) ID =  9;}
00214       else if (Opt.Contains("kaon",TString::kIgnoreCase))     {ID = 11; mass = 0.4937; if (Opt.Contains("kaon-",TString::kIgnoreCase)) ID = 12;}
00215       else if (Opt.Contains("proton",TString::kIgnoreCase))   {ID = 14; mass = 0.9383; if (Opt.Contains("antiproton",TString::kIgnoreCase)) ID = 15;}
00216       else if (Opt.Contains("deuteron",TString::kIgnoreCase)) {ID = 45; mass = 1.876;}
00217       else if (Opt.Contains("triton",TString::kIgnoreCase))   {ID = 46; mass = 2.80925;}
00218       else if (Opt.Contains("He3",TString::kIgnoreCase))      {ID = 49; mass = 2.80925;}
00219       else if (Opt.Contains("alpha",TString::kIgnoreCase))    {ID = 47; mass = 3.727;}
00220       Double_t bgMin  = 1e-2;
00221       Double_t bgMax  = 1e5;
00222       Double_t pTmin = mass*bgMin; if (pTmin <   0.1) pTmin =   0.1;
00223       Double_t pTmax = mass*bgMax; if (pTmax > 100.0) pTmax = 100.0;
00224       if (Opt.Contains("pionMIP",TString::kIgnoreCase))       {pTmin = 0.2; pTmax = 0.5;}
00225       if (Opt.Contains("50muons1GeV",TString::kIgnoreCase))   {NTRACK =   50; pTmin = pTmax = 1.0;}
00226       if (Opt.Contains("50muons0.5GeV",TString::kIgnoreCase)) {NTRACK =   50; pTmin = pTmax = 0.5;}
00227       if (Opt.Contains("1000muons1GeV",TString::kIgnoreCase)) {NTRACK = 1000; pTmin = pTmax = 1.0;}
00228       if (Opt.Contains("1muon",TString::kIgnoreCase))          NTRACK = 1;
00229       if (Opt.Contains("1alpha1GeV",TString::kIgnoreCase))    {NTRACK = 1;    pTmin = pTmax = 1.0;}
00230       TString Kine(Form("gkine %i %i %f %f -2  2 0 %f -50 50;",NTRACK,ID,pTmin,pTmax,TMath::TwoPi()));
00231       cout << "Set kinematics: " << Kine.Data() << endl;
00232       geant->Do(Kine.Data());
00233     }
00234     if (kuip) {
00235       TString Kuip(kuip);
00236       geant->Do(kuip);
00237     }
00238   }
00239   if (Last > 0)  chain->EventLoop(First,Last);
00240 }
00241 //________________________________________________________________________________
00242 void TpcRS(Int_t Last=100,
00243            const Char_t *Run = "y2009,TpcRS",//trs,fcf", // "TpcRS,fcf",
00244            //      const Char_t *fileIn = "/star/rcf/simu/rcf1207_01_225evts.fzd",
00245            const Char_t *fileIn = 0,
00246            //"/star/rcf/simu/auau200/hijing/b0_20/inverse/year2001/hadronic_on/gstardata/rcf0191_01_380evts.fzd",
00247            const Char_t *opt = "Bichsel", const Char_t *kuip = 0) {
00248   //  /star/data03/daq/2004/093/st_physics_adc_5093007_raw_2050001.daq
00249   //  /star/data03/daq/2004/fisyak/st_physics_adc_5114043_raw_2080001.daq
00250   // nofield /star/data03/daq/2004/076/st_physics_adc_5076061_raw_2060001.daq
00251   //                                   st_physics_adc_5076061_raw_4050001.daq
00252   TpcRS(1,Last,Run,fileIn,opt,kuip);
00253 }