00001
00002
00003 #include "RecoTreeModule.h"
00004
00005 #include "MessageService/MsgService.h"
00006 #include "JobControl/JobCModuleRegistry.h"
00007 #include "JobControl/JobCommand.h"
00008 #include "MinosObjectMap/MomNavigator.h"
00009
00010 #include "RerootExodus/RerootExodus.h"
00011 #include "REROOT_Classes/REROOT_NeuKin.h"
00012 #include "REROOT_Classes/REROOT_FLSHit.h"
00013 #include "Record/SimSnarlRecord.h"
00014 #include "Record/SimSnarlHeader.h"
00015 #include "Digitization/DigiSignal.h"
00016 #include "Digitization/DigiScintHit.h"
00017
00018 #include "RawData/RawRecord.h"
00019 #include "RawData/RawSnarlHeaderBlock.h"
00020 #include "RawData/RawDigit.h"
00021 #include "RawData/RawDigitDataBlock.h"
00022 #include "RawData/RawDaqSnarlHeader.h"
00023
00024 #include "CandData/CandRecord.h"
00025 #include "CandData/CandHeader.h"
00026 #include "CandDigit/CandDeMuxDigitHandle.h"
00027 #include "CandDigit/CandDeMuxDigitListHandle.h"
00028 #include "RecoBase/CandStripHandle.h"
00029 #include "RecoBase/CandStripListHandle.h"
00030 #include "RecoBase/CandTrackHandle.h"
00031 #include "RecoBase/CandTrackListHandle.h"
00032 #include "RecoBase/CandShowerHandle.h"
00033 #include "RecoBase/CandShowerListHandle.h"
00034 #include "RecoBase/CandFitTrackHandle.h"
00035 #include "RecoBase/CandFitTrackListHandle.h"
00036 #include "RecoBase/Vertex.h"
00037
00038 #include "UgliGeometry/UgliGeomHandle.h"
00039 #include "UgliGeometry/UgliPlnNode.h"
00040 #include "Validity/VldContext.h"
00041 #include "Plex/PlexStripEndId.h"
00042 #include "Plex/PlexPlaneId.h"
00043 #include "Plex/PlexSEIdAltL.h"
00044 #include "Plex/PlexSEIdAltLItem.h"
00045
00046 #include "TDirectory.h"
00047 #include "TClonesArray.h"
00048 #include "TObjArray.h"
00049 #include "TParticle.h"
00050
00051 CVSID("$Id: RecoTreeModule.cxx,v 1.6 2006/05/22 16:44:42 rhatcher Exp $");
00052
00053 JOBMODULE(RecoTreeModule,"RecoTreeModule","Testing CandTracks and CandShowers");
00054
00055 RecoTreeModule::RecoTreeModule() :
00056 fRecoFile(0),
00057 fRecoTree(0),
00058 fListInTrk("CandTrackListHandle"),
00059 fListInShw("CandShowerListHandle"),
00060 fListInFit("CandFitTrackListHandle")
00061 {
00062
00063 }
00064
00065 RecoTreeModule::~RecoTreeModule()
00066 {
00067
00068 }
00069
00070 void RecoTreeModule::BeginJob()
00071 {
00072 Int_t k;
00073
00074 Int_t tmpStrpExtr[8]={ 0, 28, 56, 76, 96, 116, 136, 164 };
00075 Int_t tmpStrpCell[28]={ 0, 1, 2, 3, 4, 5, 6,
00076 7, 8, 9, 10, 11, 12, 13,
00077 14, 15, 16, 17, 18, 19, 20,
00078 21, 22, 23, 24, 25, 26, 27 };
00079
00080 for(k=0;k<8;k++){
00081 StrpExtr[k]=tmpStrpExtr[k];
00082 }
00083
00084 for(k=0;k<28;k++){
00085 StrpCell[k]=tmpStrpCell[k];
00086 }
00087 }
00088
00089 JobCResult RecoTreeModule::Ana(const MomNavigator* mom)
00090 {
00091
00092 MSG("RecoTreeModule",Msg::kInfo) << " *** RecoTreeModule::Ana() *** " << endl;
00093 JobCResult result(JobCResult::kPassed);
00094 Int_t i,j,k;
00095
00096
00097
00098
00099
00100 run = -1;
00101 subrun = -1;
00102 snarl = -1;
00103 runtype = -1;
00104 trigsrc = -1;
00105 date = -1;
00106 time = -1;
00107 timeframe = -1;
00108 trigtimeraw = -999.9;
00109 trigtimecorr = -999.9;
00110
00111 Ndigits = 0;
00112 Nstrps = 0;
00113 Ustrps = 0;
00114 Vstrps = 0;
00115 Nplns = 0;
00116 Uplns = 0;
00117 Vplns = 0;
00118 begpln = -1;
00119 endpln = -1;
00120 Nmuplns = 0;
00121
00122 totPH = 0.0;
00123
00124 IDnu = 0;
00125 IDact = 0;
00126 IDboson = 0;
00127 IDres = 0;
00128 IDtarget = 0;
00129 x = -1.0;
00130 y = -1.0;
00131 Q2 = -1.0;
00132 W2 = -1.0;
00133 Xsection = -1.0;
00134 EMfrac = -1.0;
00135 PnuX = 0.0;
00136 PnuY = 0.0;
00137 PnuZ = 0.0;
00138 Enu = 0.0;
00139 PmuX = 0.0;
00140 PmuY = 0.0;
00141 PmuZ = 0.0;
00142 Emu = 0.0;
00143 PelX = 0.0;
00144 PelY = 0.0;
00145 PelZ = 0.0;
00146 Eel = 0.0;
00147 PhadX = 0.0;
00148 PhadY = 0.0;
00149 PhadZ = 0.0;
00150 Ehad = 0.0;
00151 PtargX = 0.0;
00152 PtargY = 0.0;
00153 PtargZ = 0.0;
00154 Etarg = 0.0;
00155 vtxX = 0.0;
00156 vtxY = 0.0;
00157 vtxZ = 0.0;
00158 muvtxX = 0.0;
00159 muvtxY = 0.0;
00160 muvtxZ = 0.0;
00161 muvtxR = 0.0;
00162 muvtxpln=-1;
00163 muendvtxX = 0.0;
00164 muendvtxY = 0.0;
00165 muendvtxZ = 0.0;
00166 muendvtxR = 0.0;
00167 muendvtxpln=-1;
00168
00169 Ntrks = 0;
00170 TRKlast = 0;
00171 for(j=0;j<3;j++){
00172 TRKplns[j] = 0;
00173 TRKstrps[j] = 0;
00174 TRKdigits[j] = 0;
00175 TRKvtxpln[j] = -1;
00176 TRKbegpln[j] = -1;
00177 TRKendpln[j] = -1;
00178 TRKpbegX[j] = 0.0;
00179 TRKpbegY[j] = 0.0;
00180 TRKpbegZ[j] = 0.0;
00181 TRKpendX[j] = 0.0;
00182 TRKpendY[j] = 0.0;
00183 TRKpendZ[j] = 0.0;
00184 TRKbegvtxX[j] = 0.0;
00185 TRKbegvtxY[j] = 0.0;
00186 TRKbegvtxZ[j] = 0.0;
00187 TRKbegtrace[j] = -999.9;
00188 TRKbegtraceZ[j] = 0.0;
00189 TRKendvtxX[j] = 0.0;
00190 TRKendvtxY[j] = 0.0;
00191 TRKendvtxZ[j] = 0.0;
00192 TRKendtrace[j] = -999.9;
00193 TRKendtraceZ[j] = -999.9;
00194 TRKdirTrue[j] = 0.0;
00195 TRKdirTimeSlope[j] = 0.0;
00196 TRKdirTimeOffset[j] = 0.0;
00197 TRKrange[j] = 0.0;
00198 TRKmomentumRange[j] = 0.0;
00199 TRKscore[j] = -1.0;
00200 TRKscoreMuEff[j] = -1.0;
00201 TRKscoreMuPur[j] = -1.0;
00202 TRKscorePiEff[j] = -1.0;
00203 TRKscorePiPur[j] = -1.0;
00204 }
00205
00206 Nfits = 0;
00207 FITlast = 0;
00208 for(j=0;j<3;j++){
00209 FITvtxpln[j] = -1;
00210 FITbegpln[j] = -1;
00211 FITendpln[j] = -1;
00212 FITcharge[j] = 0.0;
00213 FITchisq[j] = -1.0;
00214 FITmomentumCurve[j] = 0.0;
00215 FITpass[j] = 0;
00216 }
00217
00218 Nshws = 0;
00219 SHWlast = 0;
00220 for(j=0;j<3;j++){
00221 SHWplns[j] = 0;
00222 SHWstrps[j] = 0;
00223 SHWdigits[j] = 0;
00224 SHWvtxpln[j] = -1;
00225 SHWbegpln[j] = -1;
00226 SHWendpln[j] = -1;
00227 SHWtotPH[j] = 0.0;
00228 SHWpshwX[j] = 0.0;
00229 SHWpshwY[j] = 0.0;
00230 SHWpshwZ[j] = 0.0;
00231 SHWvtxX[j] = 0.0;
00232 SHWvtxY[j] = 0.0;
00233 SHWvtxZ[j] = 0.0;
00234 SHWdirTrue[j] = 0.0;
00235 SHWdirTimeSlope[j] = 0.0;
00236 SHWdirTimeOffset[j] = 0.0;
00237 SHWenergyPH[j] = 0.0;
00238 }
00239
00240
00241
00242
00243
00244
00245 MSG("RecoTreeModule",Msg::kDebug) << " *** SIMSNARL RECORD *** " << endl;
00246
00247 SimSnarlRecord* simrec = dynamic_cast<SimSnarlRecord *>(mom->GetFragment("SimSnarlRecord"));
00248 if(simrec){
00249
00250 VldContext *vldc = (VldContext*)(simrec->GetVldContext());
00251 UgliGeomHandle ugh(*vldc);
00252
00253 if(run<0) run=RerootExodus::GetRunNo();
00254 if(snarl<0) snarl=RerootExodus::GetEventNo();
00255
00256 const SimSnarlHeader* hdr = dynamic_cast<const SimSnarlHeader*>(simrec->GetSimSnarlHeader());
00257 if(hdr){
00258 if(run<0) run = hdr->GetRun();
00259 if(subrun<0) subrun = hdr->GetSubRun();
00260 if(snarl<0) snarl = hdr->GetSnarl();
00261 if(runtype<0) runtype = hdr->GetRunType();
00262 }
00263
00264 TObjArray arr(simrec->GetComponents());
00265 TIter iter(arr.MakeIterator());
00266 while(TObject* tob = (TObject*)(iter())){
00267 MSG("RecoTreeModule",Msg::kDebug) << tob->GetName() << endl;
00268
00269 if(tob->GetName()==TString("REROOT_NeuKin")){
00270 REROOT_NeuKin* nukin = dynamic_cast<REROOT_NeuKin*>(tob);
00271 MSG("RecoTreeModule",Msg::kDebug) << " ... found REROOT_NeuKin " << endl;
00272 IDnu = nukin->INu();
00273 IDact = nukin->IAction();
00274 IDtarget = nukin->ITg();
00275 IDboson = nukin->IBoson();
00276 IDres = nukin->IResonance();
00277 Xsection = nukin->Sigma();
00278 EMfrac = nukin->EMFrac();
00279 W2 = nukin->W2();
00280 Q2 = nukin->Q2();
00281 x = nukin->X();
00282 y = nukin->Y();
00283 PnuX = nukin->P4Neu()[0];
00284 PnuY = nukin->P4Neu()[1];
00285 PnuZ = nukin->P4Neu()[2];
00286 Enu = nukin->P4Neu()[3];
00287 PmuX = nukin->P4Mu1()[0];
00288 PmuY = nukin->P4Mu1()[1];
00289 PmuZ = nukin->P4Mu1()[2];
00290 Emu = nukin->P4Mu1()[3];
00291 PelX = nukin->P4El1()[0];
00292 PelY = nukin->P4El1()[1];
00293 PelZ = nukin->P4El1()[2];
00294 Eel = nukin->P4El1()[3];
00295 PhadX = nukin->P4Shw()[0];
00296 PhadY = nukin->P4Shw()[1];
00297 PhadZ = nukin->P4Shw()[2];
00298 Ehad = nukin->P4Shw()[3];
00299 PtargX = nukin->P4Tgt()[0];
00300 PtargY = nukin->P4Tgt()[1];
00301 PtargZ = nukin->P4Tgt()[2];
00302 Etarg = nukin->P4Tgt()[3];
00303 }
00304
00305 if(tob->GetName()==TString("StdHep")){
00306 TClonesArray* tpart = (TClonesArray*)(tob);
00307 MSG("RecoTreeModule",Msg::kDebug) << " found StdHep " << endl;
00308 TParticle* apart = dynamic_cast<TParticle*>(tpart->At(0));
00309 vtxX=apart->Vx(); vtxY=apart->Vy(); vtxZ=apart->Vz();
00310 }
00311 }
00312
00313 MSG("RecoTreeModule",Msg::kDebug) << " *** FLSHIT LIST *** " << endl;
00314 Int_t packedPEC,trueplane,trueext,truecell,truestrip,trueview=0,ipdg,itrack;
00315 Int_t itrkneg=-999,itrkpos=999;
00316 const TClonesArray* FLShits = (TClonesArray*)(RerootExodus::GetFLSHitList());
00317 for(i=0;i<1+FLShits->GetLast();i++){
00318 REROOT_FLSHit* FLS_hit = (REROOT_FLSHit*)FLShits->At(i);
00319 packedPEC=FLS_hit->IPackedPEC();
00320 trueplane=(Int_t)(packedPEC/65536.0);
00321 fMCStrpList[trueplane].Add(FLS_hit);
00322 if(FLS_hit->ITrack()>0 && FLS_hit->ITrack()<itrkpos) itrkpos=FLS_hit->ITrack();
00323 if(FLS_hit->ITrack()<0 && FLS_hit->ITrack()>itrkneg) itrkneg=FLS_hit->ITrack();
00324 }
00325
00326 Double_t ptot,eloss;
00327 Double_t tpos=0.,opos=0.,xpos=0.,ypos=0.,upos=0.,vpos=0.,zpos=0.,rpos=0.;
00328 Double_t um,up,vm,vp,xm,xp,ym,yp;
00329 Double_t minmuvtxX=0.,minmuvtxY=0.,minmuvtxZ=0.,minmuvtxR=0.,minptot=-1.0;
00330 Double_t maxmuvtxX=0.,maxmuvtxY=0.,maxmuvtxZ=0.,maxmuvtxR=0.,maxptot=-1.0;
00331 Int_t minmuvtxpln=0,maxmuvtxpln=0;
00332 for(i=0;i<500;i++){
00333 for(j=0;j<1+fMCStrpList[i].GetLast();j++){
00334 REROOT_FLSHit* FLS_hit = (REROOT_FLSHit*)(fMCStrpList[i].At(j));
00335 ipdg = FLS_hit->IPDG();
00336 itrack = FLS_hit->ITrack();
00337 ptot = FLS_hit->Ptot();
00338 eloss = FLS_hit->ELoss();
00339 if(ptot<0) ptot=-ptot;
00340
00341 packedPEC=FLS_hit->IPackedPEC();
00342 trueplane=(Int_t)(packedPEC/65536.0);
00343 trueext=(Int_t)((packedPEC-trueplane*65536)/256.0);
00344 truecell=packedPEC-trueplane*65536-trueext*256;
00345 truestrip=StrpExtr[trueext-1]+StrpCell[truecell-1];
00346
00347 PlexPlaneId trueplnid(vldc->GetDetector(),trueplane,0);
00348 UgliScintPlnHandle trueplnhandle = ugh.GetScintPlnHandle(trueplnid);
00349 if(trueplnhandle.GetPlaneView()==PlaneView::kU) trueview = 0;
00350 if(trueplnhandle.GetPlaneView()==PlaneView::kV) trueview = 1;
00351 PlexStripEndId truestrpid(trueplnid,truestrip);
00352 UgliStripHandle truestrphandle = ugh.GetStripHandle(truestrpid);
00353 opos=0.005*(FLS_hit->XBegin()+FLS_hit->XEnd());
00354 tpos=truestrphandle.GetTPos();
00355 if(trueview==0){ upos=tpos; vpos=-opos; }
00356 if(trueview==1){ vpos=tpos; upos=opos; }
00357 xpos=0.7071*(upos-vpos); ypos=0.7071*(upos+vpos);
00358 zpos=trueplnhandle.GetZ0();
00359
00360 rpos=4.0;
00361 up=4.0-upos; if(up<rpos) rpos=up;
00362 um=4.0+upos; if(um<rpos) rpos=um;
00363 vp=4.0-vpos; if(vp<rpos) rpos=vp;
00364 vm=4.0+vpos; if(vm<rpos) rpos=vm;
00365 xp=4.0-xpos; if(xp<rpos) rpos=xp;
00366 xm=4.0+xpos; if(xm<rpos) rpos=xm;
00367 yp=4.0-ypos; if(yp<rpos) rpos=yp;
00368 ym=4.0+ypos; if(ym<rpos) rpos=ym;
00369
00370 if( itrack>0
00371 && ( (IDnu==0&&(ipdg==13||ipdg==-13))
00372 ||(IDnu==14&&ipdg==13)||(IDnu==-14&&ipdg==-13) ) ){
00373 if(minptot<0||ptot<minptot){
00374 minmuvtxX=xpos; minmuvtxY=ypos; minmuvtxZ=zpos;
00375 minmuvtxR=rpos; minmuvtxpln=trueplane; minptot=ptot;
00376 }
00377 if(maxptot<0||ptot>maxptot){
00378 maxmuvtxX=xpos; maxmuvtxY=ypos; maxmuvtxZ=zpos;
00379 maxmuvtxR=rpos; maxmuvtxpln=trueplane; maxptot=ptot;
00380 }
00381 }
00382 }
00383 }
00384
00385 if(minptot>=0.0 && maxptot>=0.0){
00386 muvtxX = maxmuvtxX;
00387 muvtxY = maxmuvtxY;
00388 muvtxZ = maxmuvtxZ;
00389 muvtxR = maxmuvtxR;
00390 muvtxpln = maxmuvtxpln;
00391 muendvtxX = minmuvtxX;
00392 muendvtxY = minmuvtxY;
00393 muendvtxZ = minmuvtxZ;
00394 muendvtxR = minmuvtxR;
00395 muendvtxpln = minmuvtxpln;
00396 }
00397 }
00398
00399 MSG("RecoTreeModule",Msg::kDebug) << " *** END OF SIMRECORD *** " << endl;
00400
00401
00402
00403
00404
00405
00406 MSG("RecoTreeModule",Msg::kDebug) << " *** RAW RECORD *** " << endl;
00407
00408 RawRecord *rawrec = dynamic_cast<RawRecord *>(mom->GetFragment("RawRecord"));
00409 if(rawrec){
00410
00411 Double_t timestamp;
00412 const RawDaqSnarlHeader* hdr = dynamic_cast<const RawDaqSnarlHeader*>(rawrec->GetRawHeader());
00413 if(hdr){
00414 if(run<0) run = hdr->GetRun();
00415 if(subrun<0) subrun = hdr->GetSubRun();
00416 if(snarl<0) snarl = hdr->GetSnarl();
00417 if(runtype<0) runtype = hdr->GetRunType();
00418 if(trigsrc<0) trigsrc = hdr->GetTrigSrc();
00419 if(timeframe<0) timeframe = hdr->GetTimeFrameNum();
00420 }
00421
00422 TIter rdbit(rawrec->GetRawBlockIter());
00423 while(TObject* tob = (TObject*)(rdbit())){
00424 MSG("RecoTreeModule",Msg::kDebug) << tob->GetName() << endl;
00425
00426 if(tob->InheritsFrom("RawSnarlHeaderBlock")){
00427 MSG("RecoTreeModule",Msg::kDebug) << " ... found RawSnarlHeaderBlock " << endl;
00428 RawSnarlHeaderBlock* rdb = (RawSnarlHeaderBlock*)(tob);
00429 if(run<0) run = rdb->GetRun();
00430 if(subrun<0) subrun = rdb->GetSubRun();
00431 if(snarl<0) snarl = rdb->GetSnarl();
00432 if(runtype<0) runtype = rdb->GetRunType();
00433 if(trigsrc<0) trigsrc = rdb->GetTriggerSource();
00434 if(timeframe<0) timeframe = rdb->GetTimeFrameNum();
00435 timestamp = ((VldTimeStamp)(rdb->GetTriggerTime())).GetNanoSec();
00436 date = (((VldTimeStamp)(rdb->GetTriggerTime())).GetSec()-1059696000)/(3600*24);
00437 time = (((VldTimeStamp)(rdb->GetTriggerTime())).GetSec()-1059696000)%(3600*24);
00438 if(trigtimeraw<0.0) trigtimeraw=0.001*timestamp;
00439 }
00440 }
00441 }
00442
00443 MSG("RecoTreeModule",Msg::kDebug) << " *** END OF RAWRECORD *** " << endl;
00444
00445
00446
00447
00448
00449
00450 MSG("RecoTreeModule",Msg::kDebug) << " *** CANDRECORD *** " << endl;
00451
00452 CandRecord* candrec = dynamic_cast<CandRecord*>(mom->GetFragment("CandRecord","PrimaryCandidateRecord"));
00453 if(candrec){
00454
00455 const CandHeader* hdr = candrec->GetCandHeader();
00456 if(hdr){
00457 if(run<0) run = hdr->GetRun();
00458 if(snarl<0) snarl = hdr->GetSnarl();
00459 }
00460
00461 VldContext *vldc = (VldContext*)(candrec->GetVldContext());
00462 UgliGeomHandle ugh(*vldc);
00463
00464 Int_t SHW=0,TRK=0,FIT=0;
00465 TIter candbit(candrec->GetCandHandleIter());
00466 while(TObject* tob = (TObject*)(candbit())){
00467 MSG("RecoTreeModule",Msg::kDebug) << tob->GetName() << endl;
00468
00469
00470 if(tob->InheritsFrom("CandDeMuxDigitListHandle")){
00471 CandDeMuxDigitListHandle* digitlist = (CandDeMuxDigitListHandle*)(tob);
00472 MSG("RecoTreeModule",Msg::kDebug) << " ... found CandDeMuxDigitListHandle " << endl;
00473
00474 Int_t pln,xtalk;
00475 if(trigtimecorr<0.0) trigtimecorr=1.0e6*digitlist->GetAbsTime();
00476 TIter digitr(digitlist->GetDaughterIterator());
00477 while(CandDeMuxDigitHandle* digit = (CandDeMuxDigitHandle*)(digitr())){
00478 xtalk=0;
00479 if( (digit->GetDeMuxDigitFlagWord()<8)
00480 && (digit->GetDeMuxDigitFlagWord() & CandDeMuxDigit::kXTalk)==(CandDeMuxDigit::kXTalk) ){
00481 xtalk=1;
00482 }
00483 pln=digit->GetPlexSEIdAltL().GetBestSEId().GetPlane();
00484 if( pln>0 && pln<500 ){
00485 totPH+=digit->GetCharge(CalDigitType::kPE);
00486 }
00487 }
00488 }
00489
00490
00491 if(tob->InheritsFrom("CandStripListHandle")){
00492 CandStripListHandle* strplist = (CandStripListHandle*)(tob);
00493 MSG("RecoTreeModule",Msg::kDebug) << " ... found CandStripListHandle " << endl;
00494
00495 Int_t pln,vuw;
00496 TIter stritr(strplist->GetDaughterIterator());
00497 while(CandStripHandle* strip = dynamic_cast<CandStripHandle*>(stritr())){
00498 if(strip){
00499 pln=strip->GetPlane();
00500 if( pln>0 && pln<500 && strip->GetCharge()>2.0){
00501 if(begpln<0||pln<begpln) begpln=pln;
00502 if(endpln<0||pln>endpln) endpln=pln;
00503 fTotStrpList[pln].Add(strip);
00504 }
00505 }
00506 }
00507 for(i=0;i<500;i++){
00508 if(1+fTotStrpList[i].GetLast()>0){
00509 vuw=-1;
00510 for(j=0;j<1+fTotStrpList[i].GetLast();j++){
00511 CandStripHandle* strip = (CandStripHandle*)(fTotStrpList[i].At(j));
00512
00513 if(strip->GetPlaneView()==PlaneView::kU){
00514 vuw=0; Ustrps++; Ndigits+=strip->GetNDaughters();
00515 }
00516 if(strip->GetPlaneView()==PlaneView::kV){
00517 vuw=1; Vstrps++; Ndigits+=strip->GetNDaughters();
00518 }
00519 }
00520 if(vuw==0) Uplns++; if(vuw==1) Vplns++;
00521 }
00522 }
00523 Nplns=Uplns+Vplns; Nstrps=Ustrps+Vstrps;
00524 }
00525
00526 Int_t ipdg,muflag=0;
00527 if(simrec){
00528 for(i=0;i<500;i++){
00529 if(1+fMCStrpList[i].GetLast()>0 && 1+fTotStrpList[i].GetLast()>0){
00530 for(j=0;j<1+fMCStrpList[i].GetLast();j++){
00531 REROOT_FLSHit* FLS_hit = (REROOT_FLSHit*)(fMCStrpList[i].At(j));
00532 ipdg = FLS_hit->IPDG();
00533 if( (IDnu==0&&(ipdg==13||ipdg==-13))
00534 ||(IDnu==14&&ipdg==13)||(IDnu==-14&&ipdg==-13) ){
00535 muflag=1;
00536 }
00537 }
00538 if(muflag) Nmuplns++;
00539 }
00540 }
00541 }
00542
00543
00544 if(tob->InheritsFrom(fListInTrk.Data())&&!TRK){
00545 CandTrackListHandle* trklist = (CandTrackListHandle*)(tob);
00546 MSG("RecoTreeModule",Msg::kDebug) << " ... found CandTrackListHandle (" << fListInTrk.Data() << ")" << endl;
00547 TRK=1;
00548
00549 Int_t pln,flag;
00550 Double_t dir,scr,maxscr,tmpscr;
00551
00552 Ntrks=0; TRKlast=0;
00553 CandTrackHandle* trk=0;
00554 flag=1; maxscr=9999;
00555 while(trk||flag){
00556 trk=0;
00557 flag=0; tmpscr=0;
00558 TIter trkitr(trklist->GetDaughterIterator());
00559 while(CandTrackHandle* trackhandle = dynamic_cast<CandTrackHandle*>(trkitr())){
00560 if(trackhandle){
00561 if(trackhandle->GetEndPlane()>trackhandle->GetBegPlane()) scr=trackhandle->GetEndPlane()-trackhandle->GetBegPlane(); else scr=trackhandle->GetBegPlane()-trackhandle->GetEndPlane();
00562 if(scr>0&&scr<maxscr&&scr>tmpscr){
00563 tmpscr=scr; trk=trackhandle;
00564 }
00565 }
00566 }
00567 if(trk){
00568 Ntrks++; maxscr=tmpscr;
00569 if(TRKlast<3){
00570 MSG("RecoTreeModule",Msg::kDebug) << " track: " << TRKlast << endl;
00571 dir=3.0e8*trk->GetTimeSlope();
00572
00573 TRKvtxpln[TRKlast]=trk->GetVtxPlane();
00574 TRKdirTimeSlope[TRKlast]=trk->GetTimeSlope();
00575 TRKdirTimeOffset[TRKlast]=trk->GetTimeOffset();
00576 TRKmomentumRange[TRKlast]=trk->GetMomentum();
00577
00578 if(dir>=0){
00579 TRKbegpln[TRKlast]=trk->GetBegPlane();
00580 TRKbegvtxX[TRKlast]=0.7071*(trk->GetVtxU()-trk->GetVtxV());
00581 TRKbegvtxY[TRKlast]=0.7071*(trk->GetVtxU()+trk->GetVtxV());
00582 TRKbegvtxZ[TRKlast]=trk->GetVtxZ();
00583 TRKpbegX[TRKlast]=0.7071*(trk->GetDirCosU()-trk->GetDirCosV());
00584 TRKpbegY[TRKlast]=0.7071*(trk->GetDirCosU()+trk->GetDirCosV());
00585 TRKpbegZ[TRKlast]=trk->GetDirCosZ();
00586 TRKbegtrace[TRKlast]=trk->GetVtxTrace();
00587 TRKbegtraceZ[TRKlast]=trk->GetVtxTraceZ();
00588 TRKendpln[TRKlast]=trk->GetEndPlane();
00589 TRKendvtxX[TRKlast]=0.7071*(trk->GetEndU()-trk->GetEndV());
00590 TRKendvtxY[TRKlast]=0.7071*(trk->GetEndU()+trk->GetEndV());
00591 TRKendvtxZ[TRKlast]=trk->GetEndZ();
00592 TRKpendX[TRKlast]=0.7071*(trk->GetEndDirCosU()-trk->GetEndDirCosV());
00593 TRKpendY[TRKlast]=0.7071*(trk->GetEndDirCosU()+trk->GetEndDirCosV());
00594 TRKpendZ[TRKlast]=trk->GetEndDirCosZ();
00595 TRKendtrace[TRKlast]=trk->GetEndTrace();
00596 TRKendtraceZ[TRKlast]=trk->GetEndTraceZ();
00597 }
00598 else{
00599 TRKbegpln[TRKlast]=trk->GetEndPlane();
00600 TRKbegvtxX[TRKlast]=0.7071*(trk->GetEndU()-trk->GetEndV());
00601 TRKbegvtxY[TRKlast]=0.7071*(trk->GetEndU()+trk->GetEndV());
00602 TRKbegvtxZ[TRKlast]=trk->GetEndZ();
00603 TRKpbegX[TRKlast]=-0.7071*(trk->GetEndDirCosU()-trk->GetEndDirCosV());
00604 TRKpbegY[TRKlast]=-0.7071*(trk->GetEndDirCosU()+trk->GetEndDirCosV());
00605 TRKpbegZ[TRKlast]=-trk->GetEndDirCosZ();
00606 TRKbegtrace[TRKlast]=trk->GetEndTrace();
00607 TRKbegtraceZ[TRKlast]=trk->GetEndTraceZ();
00608 TRKendpln[TRKlast]=trk->GetBegPlane();
00609 TRKendvtxX[TRKlast]=0.7071*(trk->GetVtxU()-trk->GetVtxV());
00610 TRKendvtxY[TRKlast]=0.7071*(trk->GetVtxU()+trk->GetVtxV());
00611 TRKendvtxZ[TRKlast]=trk->GetVtxZ();
00612 TRKpendX[TRKlast]=-0.7071*(trk->GetDirCosU()-trk->GetDirCosV());
00613 TRKpendY[TRKlast]=-0.7071*(trk->GetDirCosU()+trk->GetDirCosV());
00614 TRKpendZ[TRKlast]=-trk->GetDirCosZ();
00615 TRKendtrace[TRKlast]=trk->GetVtxTrace();
00616 TRKendtraceZ[TRKlast]=trk->GetVtxTraceZ();
00617 }
00618
00619 TIter stritr(trk->GetDaughterIterator());
00620 while(CandStripHandle* strp = dynamic_cast<CandStripHandle*>(stritr())){
00621 if(strp){
00622 pln=strp->GetPlane();
00623 if(pln>0 && pln<500){
00624 fTrkStrpList[pln].Add(strp);
00625 }
00626 }
00627 }
00628
00629 Double_t range=0.0,ds=0.0;
00630 for(i=0;i<500;i++){
00631 if(1+fTrkStrpList[i].GetLast()>0){
00632 if(trk->IsTPosValid(i)){
00633 ds=trk->GetdS(i);
00634 if(ds>range) range=ds;
00635 }
00636 TRKplns[TRKlast]++;
00637 for(j=0;j<1+fTrkStrpList[i].GetLast();j++){
00638 CandStripHandle* strip = (CandStripHandle*)(fTrkStrpList[i].At(j));
00639 TRKstrps[TRKlast]++; TRKdigits[TRKlast]+=strip->GetNDaughters();
00640 }
00641 }
00642 }
00643 TRKrange[TRKlast]=range;
00644
00645 if(simrec){
00646 Int_t muflag,muflagtrk,piflag,piflagtrk;
00647 Int_t trueview,packedPEC,trueplane,truecell,trueext,truestrip,ipdg;
00648 Double_t trk_pln=0,trk_mu=0.0,trk_pi=0.0,notrk_mu=0.0,notrk_pi=0.0;
00649 for(i=0;i<500;i++){
00650 if(1+fMCStrpList[i].GetLast()>0 && 1+fTotStrpList[i].GetLast()>0){
00651 muflag=0; muflagtrk=0; piflag=0; piflagtrk=0;
00652 for(j=0;j<1+fMCStrpList[i].GetLast();j++){
00653 REROOT_FLSHit* FLS_hit = (REROOT_FLSHit*)(fMCStrpList[i].At(j));
00654 packedPEC=FLS_hit->IPackedPEC();
00655 trueplane=(Int_t)(packedPEC/65536.0);
00656 trueext=(Int_t)((packedPEC-trueplane*65536)/256.0);
00657 truecell=packedPEC-trueplane*65536-trueext*256;
00658 truestrip=StrpExtr[trueext-1]+StrpCell[truecell-1];
00659 ipdg = FLS_hit->IPDG();
00660 PlexPlaneId trueplnid(vldc->GetDetector(),trueplane,0);
00661 UgliScintPlnHandle trueplnhandle = ugh.GetScintPlnHandle(trueplnid);
00662 if(trueplnhandle.GetPlaneView()==PlaneView::kU) trueview = 0;
00663 if(trueplnhandle.GetPlaneView()==PlaneView::kV) trueview = 1;
00664 if( (IDnu==0&&(ipdg==13||ipdg==-13))
00665 ||(IDnu==14&&ipdg==13)||(IDnu==-14&&ipdg==-13) ){
00666 muflag=1;
00667 for(k=0;k<1+fTrkStrpList[i].GetLast();k++){
00668 CandStripHandle* strp = (CandStripHandle*)(fTrkStrpList[i].At(k));
00669 if(strp->GetStrip()-truestrip>-2 && strp->GetStrip()-truestrip<2){
00670 muflagtrk=1;
00671 }
00672 }
00673 }
00674 if(ipdg==211||ipdg==-211){
00675 piflag=1;
00676 for(k=0;k<1+fTrkStrpList[i].GetLast();k++){
00677 CandStripHandle* strp = (CandStripHandle*)(fTrkStrpList[i].At(k));
00678 if(strp->GetStrip()-truestrip>-2 && strp->GetStrip()-truestrip<2){
00679 piflagtrk=1;
00680 }
00681 }
00682 }
00683 }
00684 if(muflag){ if(muflagtrk) trk_mu+=1.0; else notrk_mu+=1.0; }
00685 if(piflag){ if(piflagtrk) trk_pi+=1.0; else notrk_pi+=1.0; }
00686 if(1+fTrkStrpList[i].GetLast()>0) trk_pln+=1.0;
00687 }
00688 }
00689 if(trk_pln>0){
00690 TRKscore[TRKlast] = (trk_mu)/(trk_pln+notrk_mu);
00691 if(trk_mu+notrk_mu>0.0){
00692 TRKscoreMuEff[TRKlast] = (trk_mu)/(trk_mu+notrk_mu);
00693 TRKscoreMuPur[TRKlast]= (trk_mu)/(trk_pln);
00694 }
00695 if(trk_pi+notrk_pi>0.0){
00696 TRKscorePiEff[TRKlast] = (trk_pi)/(trk_pi+notrk_pi);
00697 TRKscorePiPur[TRKlast] = (trk_pi)/(trk_pln);
00698 }
00699 }
00700
00701 Double_t begZ,endZ;
00702 begZ=TRKbegvtxZ[TRKlast]; endZ=TRKendvtxZ[TRKlast];
00703 if(endZ>begZ){
00704 TRKdirTrue[TRKlast]=2.0*(0.5*(endZ+begZ)-vtxZ)/(endZ-begZ);
00705 }
00706 }
00707
00708 MSG("RecoTreeModule",Msg::kDebug)
00709 << " RESULTS FROM TRACKS " << endl
00710 << " ------------------- " << endl
00711 << " TRKplns = " << TRKplns[TRKlast] << endl
00712 << " TRKstrps = " << TRKstrps[TRKlast] << endl
00713 << " TRKdigits = " << TRKdigits[TRKlast] << endl
00714 << " TRKvtxpln = " << TRKvtxpln[TRKlast] << endl
00715 << " TRKbegpln = " << TRKbegpln[TRKlast] << endl
00716 << " TRKendpln = " << TRKendpln[TRKlast] << endl
00717 << " TRKpbegX = " << TRKpbegX[TRKlast] << endl
00718 << " TRKpbegY = " << TRKpbegY[TRKlast] << endl
00719 << " TRKpbegZ = " << TRKpbegZ[TRKlast] << endl
00720 << " TRKpendX = " << TRKpendX[TRKlast] << endl
00721 << " TRKpendY = " << TRKpendY[TRKlast] << endl
00722 << " TRKpendZ = " << TRKpendZ[TRKlast] << endl
00723 << " TRKbegvtxX = " << TRKbegvtxX[TRKlast] << endl
00724 << " TRKbegvtxY = " << TRKbegvtxY[TRKlast] << endl
00725 << " TRKbegvtxZ = " << TRKbegvtxZ[TRKlast] << endl
00726 << " TRKbegtrace = " << TRKbegtrace[TRKlast] << endl
00727 << " TRKbegtraceZ = " << TRKbegtraceZ[TRKlast] << endl
00728 << " TRKendvtxX = " << TRKendvtxX[TRKlast] << endl
00729 << " TRKendvtxY = " << TRKendvtxY[TRKlast] << endl
00730 << " TRKendvtxZ = " << TRKendvtxZ[TRKlast] << endl
00731 << " TRKendtrace = " << TRKendtrace[TRKlast] << endl
00732 << " TRKendtraceZ = " << TRKendtraceZ[TRKlast] << endl
00733 << " TRKdirTrue = " << TRKdirTrue[TRKlast] << endl
00734 << " TRKdirTimeSlope = " << TRKdirTimeSlope[TRKlast] << endl
00735 << " TRKdirTimeOffset = " << TRKdirTimeOffset[TRKlast] << endl
00736 << " TRKrange = " << TRKrange[TRKlast] << endl
00737 << " TRKmomentumRange = " << TRKmomentumRange[TRKlast] << endl
00738 << " TRKscore = " << TRKscore[TRKlast] << endl
00739 << " TRKscoreMuEff = " << TRKscoreMuEff[TRKlast] << endl
00740 << " TRKscoreMuPur = " << TRKscoreMuPur[TRKlast] << endl
00741 << " TRKscorePiEff = " << TRKscorePiEff[TRKlast] << endl
00742 << " TRKscorePiPur = " << TRKscorePiPur[TRKlast] << endl;
00743
00744 TRKlast++;
00745 }
00746 }
00747 }
00748 }
00749
00750
00751
00752 if(tob->InheritsFrom(fListInShw.Data())&&!SHW){
00753 CandShowerListHandle* shwlist = (CandShowerListHandle*)(tob);
00754 MSG("RecoTreeModule",Msg::kDebug) << " ... found CandShowerListHandle (" << fListInShw.Data() << ")" << endl;
00755 SHW=1;
00756
00757 Int_t pln,flag;
00758 Double_t dir,scr,tmpscr,maxscr;
00759
00760 Nshws=0; SHWlast=0;
00761 CandShowerHandle* shw=0;
00762 flag=1; maxscr=9999;
00763 while(shw||flag){
00764 shw=0;
00765 flag=0; tmpscr=0;
00766 TIter shwitr(shwlist->GetDaughterIterator());
00767 while(CandShowerHandle* showerhandle = dynamic_cast<CandShowerHandle*>(shwitr())){
00768 if(showerhandle){
00769 scr=0;
00770 TIter stritr(showerhandle->GetDaughterIterator());
00771 while(CandStripHandle* strp = dynamic_cast<CandStripHandle*>(stritr())){
00772 if(strp) scr++;
00773 }
00774 if(scr>0&&scr<maxscr&&scr>tmpscr){
00775 tmpscr=scr; shw=showerhandle;
00776 }
00777 }
00778 }
00779 if(shw){
00780 Nshws++; maxscr=tmpscr;
00781 if(SHWlast<3){
00782 MSG("RecoTreeModule",Msg::kDebug) << " shower: " << SHWlast << endl;
00783 dir=3.0e8*shw->GetTimeSlope();
00784
00785 SHWvtxpln[SHWlast]=shw->GetVtxPlane();
00786 SHWvtxX[SHWlast]=0.7071*(shw->GetVtxU()-shw->GetVtxV());
00787 SHWvtxY[SHWlast]=0.7071*(shw->GetVtxU()+shw->GetVtxV());
00788 SHWvtxZ[SHWlast]=shw->GetVtxZ();
00789 SHWdirTimeSlope[SHWlast]=shw->GetTimeSlope();
00790 SHWdirTimeOffset[SHWlast]=shw->GetTimeOffset();
00791 SHWenergyPH[SHWlast]=shw->GetEnergy();
00792
00793 if(dir>=0){
00794 SHWbegpln[SHWlast]=shw->GetBegPlane();
00795 SHWendpln[SHWlast]=shw->GetEndPlane();
00796 SHWpshwX[SHWlast]=0.7071*(shw->GetDirCosU()-shw->GetDirCosV());
00797 SHWpshwY[SHWlast]=0.7071*(shw->GetDirCosU()+shw->GetDirCosV());
00798 SHWpshwZ[SHWlast]=shw->GetDirCosZ();
00799 }
00800 else{
00801 SHWbegpln[SHWlast]=shw->GetEndPlane();
00802 SHWendpln[SHWlast]=shw->GetBegPlane();
00803 SHWpshwX[SHWlast]=-0.7071*(shw->GetDirCosU()-shw->GetDirCosV());
00804 SHWpshwY[SHWlast]=-0.7071*(shw->GetDirCosU()+shw->GetDirCosV());
00805 SHWpshwZ[SHWlast]=-shw->GetDirCosZ();
00806 }
00807
00808 TIter stritr(shw->GetDaughterIterator());
00809 while(CandStripHandle* strp = dynamic_cast<CandStripHandle*>(stritr())){
00810 if(strp){
00811 pln=strp->GetPlane();
00812 if(pln>0 && pln<500){
00813 fShwStrpList[pln].Add(strp);
00814 }
00815 }
00816 }
00817
00818 for(i=0;i<500;i++){
00819 if(1+fShwStrpList[i].GetLast()>0){
00820 SHWplns[SHWlast]++;
00821 for(j=0;j<1+fShwStrpList[i].GetLast();j++){
00822 CandStripHandle* strip = (CandStripHandle*)(fShwStrpList[i].At(j));
00823 SHWtotPH[SHWlast]+=strip->GetCharge();
00824 SHWstrps[SHWlast]++; SHWdigits[SHWlast]+=strip->GetNDaughters();
00825 }
00826 }
00827 }
00828
00829 if(simrec){
00830 if(PhadX*PhadX+PhadY*PhadY+PhadZ*PhadZ>0.0){
00831 SHWdirTrue[SHWlast]=PhadZ/sqrt(PhadX*PhadX+PhadY*PhadY+PhadZ*PhadZ);
00832 }
00833 }
00834
00835 MSG("RecoTreeModule",Msg::kDebug)
00836 << " RESULTS FROM SHOWER " << endl
00837 << " ------------------- " << endl
00838 << " SHWplns = " << SHWplns[SHWlast] << endl
00839 << " SHWstrps = " << SHWstrps[SHWlast] << endl
00840 << " SHWdigits = " << SHWdigits[SHWlast] << endl
00841 << " SHWvtxpln = " << SHWvtxpln[SHWlast] << endl
00842 << " SHWbegpln = " << SHWbegpln[SHWlast] << endl
00843 << " SHWendpln = " << SHWendpln[SHWlast] << endl
00844 << " SHWtotPH = " << SHWtotPH[SHWlast] << endl
00845 << " SHWpshwX = " << SHWpshwX[SHWlast] << endl
00846 << " SHWpshwY = " << SHWpshwY[SHWlast] << endl
00847 << " SHWpshwZ = " << SHWpshwZ[SHWlast] << endl
00848 << " SHWvtxX = " << SHWvtxX[SHWlast] << endl
00849 << " SHWvtxY = " << SHWvtxY[SHWlast] << endl
00850 << " SHWvtxZ = " << SHWvtxZ[SHWlast] << endl
00851 << " SHWdirTrue = " << SHWdirTrue[SHWlast] << endl
00852 << " SHWdirTimeSlope = " << SHWdirTimeSlope[SHWlast] << endl
00853 << " SHWdirTimeOffset = " << SHWdirTimeOffset[SHWlast] << endl
00854 << " SHWenergyPH = " << SHWenergyPH[SHWlast] << endl;
00855
00856 SHWlast++;
00857 }
00858 }
00859 }
00860 }
00861
00862
00863 if(tob->InheritsFrom(fListInFit.Data())&&!FIT){
00864 CandFitTrackListHandle* fitlist = (CandFitTrackListHandle*)(tob);
00865 MSG("RecoTreeModule",Msg::kDebug) << " ... found CandFitTrackListHandle (" << fListInFit.Data() << ")" << endl;
00866 FIT=1;
00867
00868 Int_t flag;
00869 Double_t dir,scr,tmpscr,maxscr;
00870
00871 Nfits=0; FITlast=0;
00872 CandFitTrackHandle* fit=0;
00873 flag=1; maxscr=9999;
00874 while(fit||flag){
00875 fit=0;
00876 flag=0; tmpscr=0;
00877 TIter fititr(fitlist->GetDaughterIterator());
00878 while(CandFitTrackHandle* fithandle = dynamic_cast<CandFitTrackHandle*>(fititr())){
00879 if(fithandle){
00880 if(fithandle->GetEndPlane()>fithandle->GetBegPlane()) scr=fithandle->GetEndPlane()-fithandle->GetBegPlane(); else scr=fithandle->GetBegPlane()-fithandle->GetEndPlane();
00881 if(scr>0&&scr<maxscr&&scr>tmpscr){
00882 tmpscr=scr; fit=fithandle;
00883 }
00884 }
00885 }
00886 if(fit){
00887 Nfits++; maxscr=tmpscr;
00888 if(FITlast<3){
00889 MSG("RecoTreeModule",Msg::kDebug) << " fitted track: " << FITlast << endl;
00890 dir=3.0e8*fit->GetTimeSlope();
00891
00892 FITcharge[FITlast]=fit->GetEMCharge();
00893 FITchisq[FITlast]=fit->GetChi2();
00894 FITmomentumCurve[FITlast]=fit->GetMomentumCurve();
00895 FITpass[FITlast]=fit->GetPass();
00896
00897 if(dir>=0.0){
00898 FITvtxpln[FITlast]=fit->GetVtxPlane();
00899 FITbegpln[FITlast]=fit->GetBegPlane();
00900 FITendpln[FITlast]=fit->GetEndPlane();
00901 }
00902 else{
00903 FITvtxpln[FITlast]=fit->GetVtxPlane();
00904 FITbegpln[FITlast]=fit->GetEndPlane();
00905 FITendpln[FITlast]=fit->GetBegPlane();
00906 }
00907
00908 MSG("RecoTreeModule",Msg::kDebug)
00909 << " RESULTS FROM FITTED TRACK " << endl
00910 << " ------------------------- " << endl
00911 << " FITvtxpln = " << FITvtxpln[FITlast] << endl
00912 << " FITbegpln = " << FITbegpln[FITlast] << endl
00913 << " FITendpln = " << FITendpln[FITlast] << endl
00914 << " FITcharge = " << FITcharge[FITlast] << endl
00915 << " FITchisq = " << FITchisq[FITlast] << endl
00916 << " FITmomentumCurve = " << FITmomentumCurve[FITlast] << endl
00917 << " FITpass = " << FITpass[FITlast] << endl;
00918
00919 FITlast++;
00920 }
00921 }
00922 }
00923 }
00924
00925 MSG("RecoTreeModule",Msg::kDebug) << " ... next object?" << endl;
00926 }
00927 }
00928
00929 MSG("RecoTreeModule",Msg::kDebug) << " *** END OF CANDRECORD *** " << endl;
00930
00931 for(i=0;i<500;i++){
00932 fTotStrpList[i].Clear();
00933 fTrkStrpList[i].Clear();
00934 fShwStrpList[i].Clear();
00935 fMCStrpList[i].Clear();
00936 }
00937
00938
00939
00940
00941
00942
00943 MSG("RecoTreeModule",Msg::kDebug) << " *** SAVE TO FILE *** " << endl;
00944 if(!fRecoFile && run>0){
00945 TDirectory* tmpd = gDirectory;
00946
00947 TString filename("./results/recotree");
00948 if(run>=0){ filename.Append("."); filename+=run; }
00949 if(subrun>=0){ filename.Append("."); filename+=subrun; }
00950 if(fListInTrk==TString("CandTrackAtNuListHandle")) filename.Append(".AtNu");
00951 if(fListInTrk==TString("CandTrackSRListHandle")) filename.Append(".SR");
00952 if(fListInFit==TString("CandFitTrackAtNuListHandle")) filename.Append(".AtNu");
00953 if(fListInFit==TString("CandFitTrackSRListHandle")) filename.Append(".SR");
00954 filename.Append(".root");
00955 MSG("RecoTreeModule",Msg::kDebug) << " ..." << filename.Data() << endl;
00956
00957 fRecoFile = new TFile(filename.Data(),"RECREATE");
00958 fRecoTree = new TTree("RecoTree","RecoTree");
00959
00960 fRecoTree->SetAutoSave(100);
00961
00962 fRecoTree->Branch("run",&run,"run/I");
00963 fRecoTree->Branch("subrun",&subrun,"subrun/I");
00964 fRecoTree->Branch("snarl",&snarl,"snarl/I");
00965 fRecoTree->Branch("runtype",&runtype,"runtype/I");
00966 fRecoTree->Branch("trigsrc",&trigsrc,"trigsrc/I");
00967 fRecoTree->Branch("date",&date,"date/I");
00968 fRecoTree->Branch("time",&time,"time/I");
00969 fRecoTree->Branch("timeframe",&timeframe,"timeframe/I");
00970 fRecoTree->Branch("trigtimeraw",&trigtimeraw,"trigtimeraw/D");
00971 fRecoTree->Branch("trigtimecorr",&trigtimecorr,"trigtimecorr/D");
00972
00973 fRecoTree->Branch("Ndigits",&Ndigits,"Ndigits/I");
00974 fRecoTree->Branch("Nstrps",&Nstrps,"Nstrps/I");
00975 fRecoTree->Branch("Nplns",&Nplns,"Nplns/I");
00976 fRecoTree->Branch("Ustrps",&Ustrps,"Ustrps/I");
00977 fRecoTree->Branch("Uplns",&Uplns,"Uplns/I");
00978 fRecoTree->Branch("Vstrps",&Vstrps,"Vstrps/I");
00979 fRecoTree->Branch("Vplns",&Vplns,"Vplns/I");
00980 fRecoTree->Branch("begpln",&begpln,"begpln/I");
00981 fRecoTree->Branch("endpln",&endpln,"endpln/I");
00982 fRecoTree->Branch("Nmuplns",&Nmuplns,"Nmuplns/I");
00983
00984 fRecoTree->Branch("totPH",&totPH,"totPH/D");
00985
00986 fRecoTree->Branch("IDnu",&IDnu,"IDnu/I");
00987 fRecoTree->Branch("IDact",&IDact,"IDact/I");
00988 fRecoTree->Branch("IDboson",&IDboson,"IDboson/I");
00989 fRecoTree->Branch("IDtarget",&IDtarget,"IDtarget/I");
00990 fRecoTree->Branch("IDres",&IDres,"IDres/I");
00991 fRecoTree->Branch("x",&x,"x/D");
00992 fRecoTree->Branch("y",&y,"y/D");
00993 fRecoTree->Branch("Q2",&Q2,"Q2/D");
00994 fRecoTree->Branch("W2",&W2,"W2/D");
00995 fRecoTree->Branch("EMfrac",&EMfrac,"EMfrac/D");
00996 fRecoTree->Branch("Xsection",&Xsection,"Xsection/D");
00997 fRecoTree->Branch("PnuX",&PnuX,"PnuX/D");
00998 fRecoTree->Branch("PnuY",&PnuY,"PnuY/D");
00999 fRecoTree->Branch("PnuZ",&PnuZ,"PnuZ/D");
01000 fRecoTree->Branch("Enu",&Enu,"Enu/D");
01001 fRecoTree->Branch("PmuX",&PmuX,"PmuX/D");
01002 fRecoTree->Branch("PmuY",&PmuY,"PmuY/D");
01003 fRecoTree->Branch("PmuZ",&PmuZ,"PmuZ/D");
01004 fRecoTree->Branch("Emu",&Emu,"Emu/D");
01005 fRecoTree->Branch("PelX",&PelX,"PelX/D");
01006 fRecoTree->Branch("PelY",&PelY,"PelY/D");
01007 fRecoTree->Branch("PelZ",&PelZ,"PelZ/D");
01008 fRecoTree->Branch("Eel",&Eel,"Eel/D");
01009 fRecoTree->Branch("PhadX",&PhadX,"PhadX/D");
01010 fRecoTree->Branch("PhadY",&PhadY,"PhadY/D");
01011 fRecoTree->Branch("PhadZ",&PhadZ,"PhadZ/D");
01012 fRecoTree->Branch("Ehad",&Ehad,"Ehad/D");
01013 fRecoTree->Branch("PtargX",&PtargX,"PtargX/D");
01014 fRecoTree->Branch("PtargY",&PtargY,"PtargY/D");
01015 fRecoTree->Branch("PtargZ",&PtargZ,"PtargZ/D");
01016 fRecoTree->Branch("Etarg",&Etarg,"Etarg/D");
01017 fRecoTree->Branch("vtxX",&vtxX,"vtxX/D");
01018 fRecoTree->Branch("vtxY",&vtxY,"vtxY/D");
01019 fRecoTree->Branch("vtxZ",&vtxZ,"vtxZ/D");
01020 fRecoTree->Branch("muvtxX",&muvtxX,"muvtxX/D");
01021 fRecoTree->Branch("muvtxY",&muvtxY,"muvtxY/D");
01022 fRecoTree->Branch("muvtxZ",&muvtxZ,"muvtxZ/D");
01023 fRecoTree->Branch("muvtxR",&muvtxR,"muvtxR/D");
01024 fRecoTree->Branch("muvtxpln",&muvtxpln,"muvtxpln/I");
01025 fRecoTree->Branch("muendvtxX",&muendvtxX,"muendvtxX/D");
01026 fRecoTree->Branch("muendvtxY",&muendvtxY,"muendvtxY/D");
01027 fRecoTree->Branch("muendvtxZ",&muendvtxZ,"muendvtxZ/D");
01028 fRecoTree->Branch("muendvtxR",&muendvtxR,"muendvtxR/D");
01029 fRecoTree->Branch("muendvtxpln",&muendvtxpln,"muendvtxpln/I");
01030
01031 fRecoTree->Branch("Ntrks",&Ntrks,"Ntrks/I");
01032 fRecoTree->Branch("TRKlast",&TRKlast,"TRKlast/I");
01033 fRecoTree->Branch("TRKplns",TRKplns,"TRKplns[3]/I");
01034 fRecoTree->Branch("TRKstrps",TRKstrps,"TRKstrps[3]/I");
01035 fRecoTree->Branch("TRKdigits",TRKdigits,"TRKdigits[3]/I");
01036 fRecoTree->Branch("TRKvtxpln",TRKvtxpln,"TRKvtxpln[3]/I");
01037 fRecoTree->Branch("TRKbegpln",TRKbegpln,"TRKbegpln[3]/I");
01038 fRecoTree->Branch("TRKendpln", TRKendpln,"TRKendpln[3]/I");
01039 fRecoTree->Branch("TRKpbegX",TRKpbegX,"TRKpbegX[3]/D");
01040 fRecoTree->Branch("TRKpbegY",TRKpbegY,"TRKpbegY[3]/D");
01041 fRecoTree->Branch("TRKpbegZ",TRKpbegZ,"TRKpbegZ[3]/D");
01042 fRecoTree->Branch("TRKpendX",TRKpendX,"TRKpendX[3]/D");
01043 fRecoTree->Branch("TRKpendY",TRKpendY,"TRKpendY[3]/D");
01044 fRecoTree->Branch("TRKpendZ",TRKpendZ,"TRKpendZ[3]/D");
01045 fRecoTree->Branch("TRKbegvtxX",TRKbegvtxX,"TRKbegvtxX[3]/D");
01046 fRecoTree->Branch("TRKbegvtxY",TRKbegvtxY,"TRKbegvtxY[3]/D");
01047 fRecoTree->Branch("TRKbegvtxZ",TRKbegvtxZ,"TRKbegvtxZ[3]/D");
01048 fRecoTree->Branch("TRKbegtrace",TRKbegtrace,"TRKbegtrace[3]/D");
01049 fRecoTree->Branch("TRKbegtraceZ",TRKbegtraceZ,"TRKbegtraceZ[3]/D");
01050 fRecoTree->Branch("TRKendvtxX",TRKendvtxX,"TRKendvtxX[3]/D");
01051 fRecoTree->Branch("TRKendvtxY",TRKendvtxY,"TRKendvtxY[3]/D");
01052 fRecoTree->Branch("TRKendvtxZ",TRKendvtxZ,"TRKendvtxZ[3]/D");
01053 fRecoTree->Branch("TRKendtrace",TRKendtrace,"TRKendtrace[3]/D");
01054 fRecoTree->Branch("TRKendtraceZ",TRKendtraceZ,"TRKendtraceZ[3]/D");
01055 fRecoTree->Branch("TRKdirTrue", TRKdirTrue,"TRKdirTrue[3]/D");
01056 fRecoTree->Branch("TRKdirTimeSlope", TRKdirTimeSlope,"TRKdirTimeSlope[3]/D");
01057 fRecoTree->Branch("TRKdirTimeOffset", TRKdirTimeOffset,"TRKdirTimeOffset[3]/D");
01058 fRecoTree->Branch("TRKrange",TRKrange,"TRKrange[3]/D");
01059 fRecoTree->Branch("TRKmomentumRange", TRKmomentumRange,"TRKmomentumRange[3]/D");
01060 fRecoTree->Branch("TRKscore", TRKscore, "TRKscore[3]/D");
01061 fRecoTree->Branch("TRKscoreMuEff", TRKscoreMuEff, "TRKscoreMuEff[3]/D");
01062 fRecoTree->Branch("TRKscoreMuPur", TRKscoreMuPur, "TRKscoreMuPur[3]/D");
01063 fRecoTree->Branch("TRKscorePiEff", TRKscorePiEff, "TRKscorePiEff[3]/D");
01064 fRecoTree->Branch("TRKscorePiPur", TRKscorePiPur, "TRKscorePiPur[3]/D");
01065
01066 fRecoTree->Branch("Nfits",&Nfits,"Nfits/I");
01067 fRecoTree->Branch("FITlast",&FITlast,"FITlast/I");
01068 fRecoTree->Branch("FITvtxpln",FITvtxpln,"FITvtxpln[3]/I");
01069 fRecoTree->Branch("FITbegpln",FITbegpln,"FITbegpln[3]/I");
01070 fRecoTree->Branch("FITendpln", FITendpln,"FITendpln[3]/I");
01071 fRecoTree->Branch("FITcharge",FITcharge,"FITcharge[3]/D");
01072 fRecoTree->Branch("FITchisq", FITchisq, "FITchisq[3]/D");
01073 fRecoTree->Branch("FITmomentumCurve", FITmomentumCurve,"FITmomentumCurve[3]/D");
01074 fRecoTree->Branch("FITpass", FITpass, "FITpass[3]/I");
01075
01076 fRecoTree->Branch("Nshws",&Nshws,"Nshws/I");
01077 fRecoTree->Branch("SHWlast",&SHWlast,"SHWlast/I");
01078 fRecoTree->Branch("SHWplns",SHWplns,"SHWplns[3]/I");
01079 fRecoTree->Branch("SHWstrps",SHWstrps,"SHWstrps[3]/I");
01080 fRecoTree->Branch("SHWdigits",SHWdigits,"SHWdigits[3]/I");
01081 fRecoTree->Branch("SHWvtxpln",SHWvtxpln,"SHWvtxpln[3]/I");
01082 fRecoTree->Branch("SHWbegpln",SHWbegpln,"SHWbegpln[3]/I");
01083 fRecoTree->Branch("SHWendpln",SHWendpln,"SHWendpln[3]/I");
01084 fRecoTree->Branch("SHWtotPH",SHWtotPH,"SHWtotPH[3]/D");
01085 fRecoTree->Branch("SHWpshwX",SHWpshwX,"SHWpshwX[3]/D");
01086 fRecoTree->Branch("SHWpshwY",SHWpshwY,"SHWpshwY[3]/D");
01087 fRecoTree->Branch("SHWpshwZ",SHWpshwZ,"SHWpshwZ[3]/D");
01088 fRecoTree->Branch("SHWvtxX",SHWvtxX,"SHWvtxX[3]/D");
01089 fRecoTree->Branch("SHWvtxY",SHWvtxY,"SHWvtxY[3]/D");
01090 fRecoTree->Branch("SHWvtxZ",SHWvtxZ,"SHWvtxZ[3]/D");
01091 fRecoTree->Branch("SHWdirTrue",SHWdirTrue,"SHWdirTrue[3]/D");
01092 fRecoTree->Branch("SHWdirTimeSlope",SHWdirTimeSlope,"SHWdirTimeSlope[3]/D");
01093 fRecoTree->Branch("SHWdirTimeOffset",SHWdirTimeOffset,"SHWdirTimeOffset[3]/D");
01094 fRecoTree->Branch("SHWenergyPH",SHWenergyPH,"SHWenergyPH[3]/D");
01095
01096 tmpd->cd();
01097 }
01098
01099 if(fRecoFile){
01100 MSG("RecoTreeModule",Msg::kDebug) << " ***** FILL TREE ***** " << endl;
01101 TDirectory* tmpd = gDirectory;
01102 fRecoFile->cd();
01103 fRecoTree->Fill();
01104 tmpd->cd();
01105 MSG("RecoTreeModule",Msg::kDebug) << " ***** TREE FILLED ***** " << endl;
01106 }
01107
01108 MSG("RecoTreeModule",Msg::kInfo) << " *** RecoTreeModule::Ana() FINISHED *** " << endl;
01109 return result;
01110 }
01111
01112 const Registry& RecoTreeModule::DefaultConfig() const
01113 {
01114 MSG("RecoTreeModule",Msg::kInfo) << " *** RecoTreeModule::DefaultConfig() *** " << endl;
01115 static Registry r;
01116 r.SetName("RecoTreeModule.config.default");
01117 r.UnLockValues();
01118 r.Set("ListInTrk",fListInTrk.Data());
01119 r.Set("ListInShw",fListInShw.Data());
01120 r.Set("ListInFit",fListInFit.Data());
01121 r.LockValues();
01122 return r;
01123 }
01124
01125 void RecoTreeModule::Config(const Registry &r)
01126 {
01127 MSG("RecoTreeModule",Msg::kInfo) << " *** RecoTreeModule::Config() *** " << endl;
01128 const char* tmpchar = 0;
01129 if(r.Get("ListInTrk",tmpchar)) fListInTrk = tmpchar;
01130 if(r.Get("ListInShw",tmpchar)) fListInShw = tmpchar;
01131 if(r.Get("ListInFit",tmpchar)) fListInFit = tmpchar;
01132 MSG("RecoTreeModule",Msg::kDebug) << " configuration: " << endl
01133 << " ListInTrk=" << fListInTrk.Data() << endl
01134 << " ListInShw=" << fListInShw.Data() << endl
01135 << " ListInFit=" << fListInFit.Data() << endl;
01136 }
01137
01138 void RecoTreeModule::HandleCommand(JobCommand *command)
01139 {
01140 TString cmd = command->PopCmd();
01141 if(cmd=="Set"){
01142 TString opt = command->PopOpt();
01143 }
01144 }
01145
01146 void RecoTreeModule::EndJob()
01147 {
01148 if(fRecoFile){
01149 MSG("RecoTreeModule",Msg::kInfo) << " ***** WRITING TO FILE " << fRecoFile->GetName() << " ***** " << endl;
01150 TDirectory* tmpd = gDirectory;
01151 fRecoFile->cd();
01152 fRecoTree->Write();
01153 fRecoFile->Close();
01154 tmpd->cd();
01155 MSG("RecoTreeModule",Msg::kInfo) << " ***** FILE CLOSED ***** " << endl;
01156 }
01157 }
