AlgFarDetEvent Class Reference

#include <AlgFarDetEvent.h>

Inheritance diagram for AlgFarDetEvent:

AlgBase List of all members.

Public Member Functions

 AlgFarDetEvent ()
virtual ~AlgFarDetEvent ()
virtual void RunAlg (AlgConfig &ac, CandHandle &ch, CandContext &cx)
virtual void Trace (const char *c) const

Private Member Functions

void SetTrackProperties (FarDetEventHandle *evt, CandFitTrackHandle *fit, CandTrackHandle *trk)
void SetShowerProperties (FarDetEventHandle *evt, CandShowerHandle *shw)
bool Contained (const CandTrackHandle *track)
double DeDx (double emu)
void GetShower (const CandShowerHandle *shower, const CandTrackHandle *track, vector< double > &TempShw, vector< double > &TempTrk)
void GetPulseHeight (vector< double > &TempShw, vector< double > &TempTrk, double &totph, const double &deweight)
double GetDeWeightFactor (const double &LinearEhad, AlgConfig &ac)
double CalibrateShowerEnergy (const double &totph, const double &deweight, AlgConfig &ac)

Detailed Description

Definition at line 15 of file AlgFarDetEvent.h.


Constructor & Destructor Documentation

AlgFarDetEvent::AlgFarDetEvent (  ) 

Definition at line 32 of file AlgFarDetEvent.cxx.

00033 {
00034 
00035 }   

AlgFarDetEvent::~AlgFarDetEvent (  )  [virtual]

Definition at line 37 of file AlgFarDetEvent.cxx.

00038 {
00039 
00040 }


Member Function Documentation

double AlgFarDetEvent::CalibrateShowerEnergy ( const double &  totph,
const double &  deweight,
AlgConfig ac 
) [private]

Definition at line 597 of file AlgFarDetEvent.cxx.

References Registry::GetDouble(), Registry::KeyExists(), Msg::kVerbose, Msg::kWarning, MSG, and warn().

Referenced by RunAlg().

00598 {
00599 
00600   // Initialize shower energy
00601   double cal_e(0.0);
00602 
00603   if( totph<=0 ) 
00604     {
00605       MSG("FarDetEvent",Msg::kVerbose) << " 0 MIPs -> 0 GeV" << endl;
00606       return cal_e;
00607     }
00608 
00609   // Get Calibration Constants from AlgConfig.
00610   if( deweight>=1.0 )
00611     {
00612       // USE LINEAR CONSTANTS
00613       // ====================
00614       double fCCLinLowScale(0.0); 
00615       double fCCLinLowC0(0.0); 
00616       double fCCLinLowC1(0.0); 
00617       double fCCLinLowC2(0.0); 
00618       double fCCLinHiC0(0.0); 
00619       double fCCLinHiC1(0.0);
00620  
00621       bool warn =false;
00622       TString errstring="AlgFarDetEvent::CalibrateShowerEnergy - Linear - Failed to get AlgConfig info for the following keys: ";
00623       
00624       if( ac.KeyExists("CCLinLowScale") ) fCCLinLowScale = ac.GetDouble("CCLinLowScale"); 
00625       else {errstring+="CCLinLowScale, "; warn =true; }
00626       if( ac.KeyExists("CCLinLowC0") )    fCCLinLowC0 = ac.GetDouble("CCLinLowC0");
00627       else {errstring+="CCLinLowC0, "; warn =true; }
00628       if( ac.KeyExists("CCLinLowC1") )    fCCLinLowC1 = ac.GetDouble("CCLinLowC1");
00629       else {errstring+="CCLinLowC1, "; warn =true; }
00630       if( ac.KeyExists("CCLinLowC2") )    fCCLinLowC2 = ac.GetDouble("CCLinLowC2");
00631       else {errstring+="CCLinLowC2, "; warn =true; }
00632       if( ac.KeyExists("CCLinHiC0") )     fCCLinHiC0 = ac.GetDouble("CCLinHiC0");
00633       else {errstring+="CCLinHiC0, "; warn =true; }
00634       if( ac.KeyExists("CCLinHiC1") )     fCCLinHiC1 = ac.GetDouble("CCLinHiC1"); 
00635       else {errstring+="CCLinHiC1"; warn =true; }
00636 
00637       if(warn) MSG("FarDetEvent",Msg::kWarning) << errstring.Data() <<endl;
00638 
00639       if( totph<fCCLinLowScale )
00640         {
00641           MSG("FarDetEvent",Msg::kVerbose) << "  using linear constants, low scale... " << endl;
00642           cal_e = fCCLinLowC0 + fCCLinLowC1*totph + fCCLinLowC2*totph*totph
00643                    + 0.21*sqrt(totph)*(1.0-totph/25.0)*exp(-totph/25.0);
00644         }
00645       else
00646         {
00647           MSG("FarDetEvent",Msg::kVerbose) << "  using linear constants, high scale... " << endl;
00648           cal_e = fCCLinHiC0 + fCCLinHiC1*totph;
00649         }
00650     }
00651   else
00652     {
00653       // USE DEWEIGHTED CONSTANTS
00654       // ========================
00655       double fCCWtLowScale(0.0);
00656       double fCCWtMidScale(0.0);
00657       double fCCWtLowC0(0.0);
00658       double fCCWtLowC1(0.0);
00659       double fCCWtLowC2(0.0);
00660       double fCCWtLowC3(0.0);
00661       double fCCWtLowC4(0.0);
00662       double fCCWtMidC0(0.0);
00663       double fCCWtMidC1(0.0);
00664       double fCCWtHiC0(0.0);
00665       double fCCWtHiC1(0.0);
00666 
00667       bool warn =false;
00668       TString errstring="AlgFarDetEvent::CalibrateShowerEnergy - Deweight - Failed to get AlgConfig info for the following keys: ";
00669       
00670       if( ac.KeyExists("CCWtLowScale") ) fCCWtLowScale = ac.GetDouble("CCWtLowScale");  
00671       else {errstring+="CCWtLowScale, "; warn =true; }
00672       if( ac.KeyExists("CCWtMidScale") ) fCCWtMidScale = ac.GetDouble("CCWtMidScale");  
00673       else {errstring+="CCWtMidScale, "; warn =true; } 
00674       if( ac.KeyExists("CCWtLowC0") )    fCCWtLowC0 = ac.GetDouble("CCWtLowC0");
00675       else {errstring+="CCWtLowC0, "; warn =true; }
00676       if( ac.KeyExists("CCWtLowC1") )    fCCWtLowC1 = ac.GetDouble("CCWtLowC1");
00677       else {errstring+="CCWtLowC1, "; warn =true; }
00678       if( ac.KeyExists("CCWtLowC2") )    fCCWtLowC2 = ac.GetDouble("CCWtLowC2");
00679       else {errstring+="CCWtLowC2, "; warn =true; }
00680       if( ac.KeyExists("CCWtLowC3") )    fCCWtLowC3 = ac.GetDouble("CCWtLowC3");
00681       else {errstring+="CCWtLowC3, "; warn =true; }
00682       if( ac.KeyExists("CCWtLowC4") )    fCCWtLowC4 = ac.GetDouble("CCWtLowC4");
00683       else {errstring+="CCWtLowC4, "; warn =true; }
00684       if( ac.KeyExists("CCWtMidC0") )    fCCWtMidC0 = ac.GetDouble("CCWtMidC0");
00685       else {errstring+="CCWtMidC0, "; warn =true; }
00686       if( ac.KeyExists("CCWtMidC1") )    fCCWtMidC1 = ac.GetDouble("CCWtMidC1"); 
00687       else {errstring+="CCWtMidC1"; warn =true; }
00688       if( ac.KeyExists("CCWtHiC0") )     fCCWtHiC0 = ac.GetDouble("CCWtHiC0");
00689       else {errstring+="CCWtHiC0, "; warn =true; }
00690       if( ac.KeyExists("CCWtHiC1") )     fCCWtHiC1 = ac.GetDouble("CCWtHiC1"); 
00691       else {errstring+="CCWtHiC1"; warn =true; }
00692 
00693       if(warn) MSG("FarDetEvent",Msg::kWarning) << errstring.Data() <<endl;
00694 
00695       if( totph<fCCWtLowScale )
00696         {
00697           MSG("FarDetEvent",Msg::kVerbose) << "  using deweighted constants, low scale... " << endl;
00698           cal_e = fCCWtLowC0 +
00699                   fCCWtLowC1*totph + 
00700                   fCCWtLowC2*totph*totph +
00701                   fCCWtLowC3*totph*totph*totph +
00702                   fCCWtLowC4*totph*totph*totph*totph
00703                    + 0.18*sqrt(totph)*(1.0-totph/5.0)*exp(-totph/10.0);
00704         }
00705       else if( totph<fCCWtMidScale )
00706         {
00707           MSG("FarDetEvent",Msg::kVerbose) << "  using deweighted constants, middle scale... " << endl;
00708           cal_e = fCCWtMidC0 + fCCWtMidC1*totph;
00709         }
00710       else 
00711         {
00712           MSG("FarDetEvent",Msg::kVerbose) << "  using deweighted constants, high scale... " << endl;
00713           cal_e = fCCWtHiC0 + fCCWtHiC1*totph;
00714         }
00715 
00716     }
00717 
00718   MSG("FarDetEvent",Msg::kVerbose) << "   " << totph << " MIPs -> " << cal_e << " GeV" << endl;
00719 
00720   return cal_e;
00721 }

bool AlgFarDetEvent::Contained ( const CandTrackHandle track  )  [private]

Definition at line 293 of file AlgFarDetEvent.cxx.

References CandRecoHandle::GetEndPlane(), CandRecoHandle::GetEndU(), CandRecoHandle::GetEndV(), CandRecoHandle::GetVtxPlane(), CandRecoHandle::GetVtxU(), and CandRecoHandle::GetVtxV().

00294 {
00295   // vertex containment
00296   double VtxU = track->GetVtxU();
00297   double VtxV = track->GetVtxV();
00298   double VtxR = sqrt(VtxU*VtxU+VtxV*VtxV);
00299   bool VtxContained = (VtxR<3.75 && VtxR>0.4) && ((track->GetVtxPlane()>5 && track->GetVtxPlane()<243) || (track->GetVtxPlane()>254 && track->GetVtxPlane()<480));
00300 
00301   // end containment
00302   double EndU = track->GetEndU();
00303   double EndV = track->GetEndV();
00304   double EndR = sqrt(EndU*EndU+EndV*EndV);
00305   bool EndContained = (EndR<3.75 && EndR>0.4) && ((track->GetEndPlane()>5 && track->GetEndPlane()<243) || (track->GetEndPlane()>254 && track->GetEndPlane()<480));
00306 
00307   return (VtxContained && EndContained);
00308 }

double AlgFarDetEvent::DeDx ( double  emu  )  [private]

Definition at line 310 of file AlgFarDetEvent.cxx.

References MuELoss::a, MuELoss::a_2, Munits::m, Munits::mm, MuELoss::Na, and Mphysical::pi.

Referenced by GetShower().

00311 {
00312   double dedx = 0.;
00313   double mm = 105.658389;
00314   double mm_2 = mm*mm;
00315   
00316   if(emu*emu-mm_2<0.){return 0.;}
00317   double a_2     = 1./(137.036*137.036);
00318   double pi = 3.141;
00319   double Na = 6.023;
00320   double lamda_2 = 3.8616*3.8616;
00321   double Z_A = 0.5377;
00322   double me  = 0.51099906;
00323   double me_2 = me*me; 
00324   double beta = sqrt(emu*emu-mm_2)/emu;
00325   double beta_2 = beta*beta;
00326   double gamma = emu/105.658389;
00327   double gamma_2 = gamma*gamma;
00328   double I = 68.7e-6;
00329   double I_2 = I*I;
00330   double p = sqrt(emu*emu-mm_2);
00331   double p_2 = p*p;
00332   double Em = 2 * me * p_2 / ( me_2 + mm_2 + 2*me*emu );
00333   double Em_2= Em*Em;
00334   double emu_2 = emu*emu;
00335   double X0 = 0.165;
00336   double X1 = 2.503;
00337   double X = log10(beta*gamma);
00338   double a = 0.165;
00339   double C = -3.3;
00340   double m = 3.222;
00341   double d = 0;
00342   if(X0 < X && X < X1){d = 4.6052 * X + a * pow(X1-X,m) + C;}
00343   if(X > X1){d = 4.6052 * X + C;}
00344 
00345   dedx =  a_2 * 2*pi * Na * lamda_2 * Z_A * (me / beta_2) *( log( 2*me*beta_2*gamma_2*Em/I_2 ) - 2*beta_2 + 0.25*(Em_2/emu_2) - d );
00346 
00347   return 10*dedx;
00348 }

double AlgFarDetEvent::GetDeWeightFactor ( const double &  LinearEhad,
AlgConfig ac 
) [private]

Definition at line 509 of file AlgFarDetEvent.cxx.

References Registry::GetDouble(), Registry::KeyExists(), Msg::kVerbose, Msg::kWarning, MSG, and warn().

Referenced by RunAlg().

00510 {
00511   double deweightfactorCC(1.0);
00512 
00513   double fdeweightLowScale(0.0);
00514   double fdeweightC0(0.0);
00515   double fdeweightC1(0.0);
00516   double fdeweightC2(0.0);
00517   double fdeweightC3(0.0);
00518 
00519   bool warn =false;
00520   TString errstring="AlgFarDetEvent::GetDeWeightFactor - Failed to get AlgConfig info for the following keys: ";
00521   
00522   if( ac.KeyExists("deweightLowScale") ) fdeweightLowScale = ac.GetDouble("deweightLowScale");
00523   else {errstring+="deweightLowScale, "; warn =true; }
00524   if( ac.KeyExists("deweightC0") )fdeweightC0 = ac.GetDouble("deweightC0");
00525   else {errstring+="deweightC0, "; warn =true; }
00526   if( ac.KeyExists("deweightC1") )fdeweightC1 = ac.GetDouble("deweightC1");
00527   else {errstring+="deweightC1, "; warn =true; }
00528   if( ac.KeyExists("deweightC2") )fdeweightC2 = ac.GetDouble("deweightC2");
00529   else {errstring+="deweightC2, "; warn =true; }
00530   if( ac.KeyExists("deweightC3") )fdeweightC3 = ac.GetDouble("deweightC3");
00531   else {errstring+="deweightC3"; warn =true; }
00532   if(warn) MSG("FarDetEvent",Msg::kWarning) << errstring.Data() <<endl;
00533 
00534   if(LinearEhad < fdeweightLowScale)
00535     {
00536       deweightfactorCC = fdeweightC0 + 
00537                          fdeweightC1*LinearEhad +
00538                          fdeweightC2*LinearEhad*LinearEhad +
00539                          fdeweightC3*LinearEhad*LinearEhad*LinearEhad;
00540     }
00541   
00542   MSG("FarDetEvent",Msg::kVerbose) << " LinearEnergy=" << LinearEhad << " DeweightFactor=" << deweightfactorCC << endl;
00543 
00544   return deweightfactorCC;
00545 }

void AlgFarDetEvent::GetPulseHeight ( vector< double > &  TempShw,
vector< double > &  TempTrk,
double &  totph,
const double &  deweight 
) [private]

Definition at line 551 of file AlgFarDetEvent.cxx.

References Msg::kVerbose, and MSG.

Referenced by RunAlg().

00552 {
00553 
00554   // zero pulse height
00555   totph=0.0;
00556   
00557   double ph;
00558   double pulseheight_check_1(0);
00559   double pulseheight_check_2(0);
00560   double pulseheight_check_3(0);
00561 
00562   unsigned int nstrips = TempShw.size();
00563   if( nstrips<1 ) return;
00564 
00565   if( deweight>=1.0 )
00566     {
00567       for(unsigned int i=0; i<nstrips; ++i)
00568         {
00569           ph = TempShw[i]-TempTrk[i];
00570           pulseheight_check_1 += TempShw[i];
00571           pulseheight_check_2 += TempTrk[i];
00572           totph += ph;
00573         }
00574       MSG("FarDetEvent",Msg::kVerbose) << "    ... linear check: ShwPH=" << pulseheight_check_1 << ", TrkPH=" << pulseheight_check_2 << endl;
00575     }
00576   else
00577     {
00578       for(unsigned int i=0; i<nstrips; ++i)
00579         {
00580           if( TempShw[i]-1.4*TempTrk[i]>0.0 )
00581             {
00582               ph = TempShw[i]-1.4*TempTrk[i];
00583               pulseheight_check_3 += ph;
00584               totph += pow(ph,deweight);
00585             }
00586         }
00587         MSG("FarDetEvent",Msg::kVerbose) << "    ... deweight check: ShwPH-TrkPH=" << pulseheight_check_3 << endl;   
00588     }
00589 
00590   return;
00591 }

void AlgFarDetEvent::GetShower ( const CandShowerHandle shower,
const CandTrackHandle track,
vector< double > &  TempShw,
vector< double > &  TempTrk 
) [private]

Definition at line 357 of file AlgFarDetEvent.cxx.

References DeDx(), Calibrator::GetAttenCorrected(), CandStripHandle::GetCharge(), CandHandle::GetDaughterIterator(), Calibrator::GetMIP(), CandTrackHandle::GetMomentum(), CandFitTrackHandle::GetMomentumCurve(), CandFitTrackHandle::GetMomentumRange(), CandStripHandle::GetPlane(), CandStripHandle::GetPlaneView(), CandStripHandle::GetStrip(), CandStripHandle::GetStripEndId(), CandShowerHandle::GetU(), CandShowerHandle::GetV(), CandRecoHandle::GetVtxDirCosZ(), Calibrator::Instance(), StripEnd::kNegative, StripEnd::kPositive, CalDigitType::kSigCorr, PlaneView::kU, PlaneView::kV, Msg::kVerbose, and MSG.

Referenced by RunAlg().

00361 {
00362 
00363   Calibrator& cal = Calibrator::Instance();
00364   int indx(-1),pln(-1),strp(-1);
00365 
00366 
00367   // make map of shower strips
00368   // ========================
00369   map<int,int> isShowerStrip;
00370   if( shower ) 
00371     {
00372       TIter shwitr(shower->GetDaughterIterator());
00373       while(CandStripHandle* strip = (CandStripHandle*)(shwitr()))
00374         { 
00375           if( strip )
00376             {
00377               pln = strip->GetPlane();
00378               strp = strip->GetStrip();
00379               indx = strp + 192*pln;
00380               isShowerStrip[indx] = 1;
00381             }
00382         }
00383     }
00384 
00385 
00386   // make map of track strips
00387   // ========================
00388   map<int,int> isTrackStrip;
00389   
00390   double sharedplanes(0.0);
00391   if( track ) 
00392     {
00393       int minpln(-999),maxpln(-999);
00394       TIter trkitr(track->GetDaughterIterator());
00395       while(CandStripHandle* strip = (CandStripHandle*)(trkitr()))
00396         { 
00397           if( strip )
00398             {
00399               pln = strip->GetPlane();
00400               strp = strip->GetStrip();
00401               indx = strp + 192*pln;
00402               isTrackStrip[indx] = 1;
00403 
00404               if( isShowerStrip[indx]==1 )
00405                 {
00406                   if(minpln<0 || pln<minpln) minpln=pln;
00407                   if(maxpln<0 || pln>maxpln) maxpln=pln;
00408                 }
00409             }
00410         }
00411     
00412       if( maxpln>=0 && minpln>=0 && maxpln-minpln>=0 )
00413         {
00414           sharedplanes = 1+maxpln-minpln;
00415         }
00416     }
00417 
00418   MSG("FarDetEvent",Msg::kVerbose) << "    ... shared planes = " << sharedplanes << endl;    
00419 
00420 
00421   // Calculate track correction (MIPs per strip) to remove from shower
00422   // =================================================================
00423   double trackcorrection = 0.0;
00424   if(track)
00425     { 
00426        double pmu = track->GetMomentum();
00427        double trackdircosZ = fabs(1.0/(track->GetVtxDirCosZ())); 
00428        MSG("FarDetEvent",Msg::kVerbose) << "    ... track direction = " << trackdircosZ << endl;
00429         
00430        const CandFitTrackHandle* fit = dynamic_cast<const CandFitTrackHandle*>(track);
00431        if( fit )
00432          {
00433            if( this->Contained(track) ) pmu = fit->GetMomentumRange(); 
00434            else pmu = fit->GetMomentumCurve();
00435          }
00436 
00437        // (Note: alternative dEdX calculation)
00438        // Material  poly(MinosMaterial::ePolystyreneMinos);
00439        // BetheBlochModel bbmodel(poly);
00440        // double dedx = bbmodel.dE_dx(emu*1000.)/1.985; 
00441 
00442        double dedx = DeDx(1000.0*pmu)/1.985;
00443        double dedx_end = dedx;
00444        double pmuend = pmu - sharedplanes/(30.*trackdircosZ);
00445        MSG("FarDetEvent",Msg::kVerbose) << "    ... muon momentum: vtx = " << pmu << " end = " << pmuend << endl;
00446        if( pmuend>0.3 ) // beta*gamma~3 => p~300MeV
00447          {
00448            dedx_end = DeDx(1000.0*pmuend)/1.985;
00449          }
00450        
00451        double avg_dedx = 0.5*(dedx+dedx_end);
00452        MSG("FarDetEvent",Msg::kVerbose) << "    ... average dedx per strip (MIPs) = " << avg_dedx << endl;
00453 
00454        double temp_trackdircosZ = fabs(trackdircosZ);
00455        if( temp_trackdircosZ<0.5 ) temp_trackdircosZ = 0.5;
00456        trackcorrection = avg_dedx/temp_trackdircosZ;
00457     }
00458   
00459   MSG("FarDetEvent",Msg::kVerbose) << "    ... track correction = " << trackcorrection << endl;
00460 
00461   // Calculate the MIP values for each strip in the shower
00462   // =====================================================
00463   
00464   if(shower)
00465     {
00466       double opos(0.0);
00467       double sigmapN(0.0),sigmapP(0.0);
00468       double shwMIP(0.0),trkMIP(0.0); 
00469       TIter shwitr(shower->GetDaughterIterator());
00470       while(CandStripHandle* strip = (CandStripHandle*)(shwitr()))
00471         { 
00472           if( strip )
00473             {
00474               pln = strip->GetPlane();
00475               strp = strip->GetStrip();
00476               indx = strp + 192*pln;
00477 
00478               shwMIP = 0.0; 
00479               trkMIP = 0.0;
00480               opos = 0.0;
00481 
00482               if(strip->GetPlaneView()==PlaneView::kU) opos = -1.0*shower->GetV(pln); 
00483               if(strip->GetPlaneView()==PlaneView::kV) opos =  1.0*shower->GetU(pln); 
00484               if( fabs(opos)<999.9 )
00485                 {
00486                   sigmapN=cal.GetAttenCorrected(strip->GetCharge(StripEnd::kNegative,CalDigitType::kSigCorr), opos, strip->GetStripEndId(StripEnd::kNegative) );
00487                   shwMIP+=cal.GetMIP(sigmapN, strip->GetStripEndId(StripEnd::kNegative) );
00488                   sigmapP=cal.GetAttenCorrected(strip->GetCharge(StripEnd::kPositive,CalDigitType::kSigCorr), opos, strip->GetStripEndId(StripEnd::kPositive) );
00489                   shwMIP+=cal.GetMIP(sigmapP, strip->GetStripEndId(StripEnd::kPositive) );
00490                 }
00491               
00492               if( isTrackStrip[indx]==1 )
00493                 {
00494                   trkMIP = trackcorrection;
00495                 }
00496               TempShw.push_back(shwMIP);
00497               TempTrk.push_back(trkMIP);
00498             }   
00499         }  
00500     }
00501   
00502   return;
00503 }

void AlgFarDetEvent::RunAlg ( AlgConfig ac,
CandHandle ch,
CandContext cx 
) [virtual]

Implements AlgBase.

Definition at line 54 of file AlgFarDetEvent.cxx.

References CalibrateShowerEnergy(), CandContext::GetCandRecord(), CandContext::GetDataIn(), GetDeWeightFactor(), VHS::GetPlane(), CandEventHandle::GetPrimaryShower(), CandEventHandle::GetPrimaryTrack(), GetPulseHeight(), GetShower(), RecMinos::GetVldContext(), Calibrator::Instance(), Msg::kDebug, Msg::kWarning, MSG, CandEventHandle::SetEnergy(), FarDetEventHandle::SetMaxPlane(), FarDetEventHandle::SetMaxZ(), FarDetEventHandle::SetMinPlane(), FarDetEventHandle::SetMinZ(), SetShowerProperties(), FarDetEventHandle::SetShwEnergyGeV(), FarDetEventHandle::SetShwEnergyGeVDeweighted(), FarDetEventHandle::SetShwEnergyGeVLinear(), FarDetEventHandle::SetShwMipsDeweighted(), FarDetEventHandle::SetShwMipsLinear(), and SetTrackProperties().

00055 {
00056   MSG("FarDetEvent",Msg::kDebug) << " AlgFarDetEvent::RunAlg(....) " << endl;
00057 
00058   FarDetEventHandle& myevent = dynamic_cast<FarDetEventHandle&>(ch);
00059 
00060   Calibrator& cal = Calibrator::Instance();
00061   cal.ReInitialise(*(cx.GetCandRecord()->GetVldContext()));
00062 
00063   // Unpack CandContext
00064   // ==================
00065   const TObjArray* eventArray = dynamic_cast<const TObjArray*>(cx.GetDataIn());
00066 
00067   CandFitTrackHandle* fit = (CandFitTrackHandle*)(eventArray->At(0));
00068   CandTrackHandle* trk = (CandTrackHandle*)(eventArray->At(1));
00069   CandShowerHandle* shw = (CandShowerHandle*)(eventArray->At(2));
00070   // TObjArray* HadronicShowers = (TObjArray*)(eventArray->At(3));
00071   // TObjArray* BremsAndDeltas = (TObjArray*)(eventArray->At(4));
00072   // TObjArray* AssociatedStrips = (TObjArray*)(eventArray->At(5));
00073   TObjArray* AllStrips = (TObjArray*)(eventArray->At(6));
00074   
00075   // Set primary track and shower properties
00076   // =======================================
00077   MSG("FarDetEvent",Msg::kDebug) << " Setting track properties " << endl;
00078   this->SetTrackProperties(&myevent, fit, trk);
00079 
00080   MSG("FarDetEvent",Msg::kDebug) << " Setting shower properties " << endl;
00081   this->SetShowerProperties(&myevent, shw);
00082 
00083   if( myevent.GetPrimaryTrack()==0 && myevent.GetPrimaryShower()==0 )
00084     {
00085       MSG("FarDetEvent",Msg::kWarning) << " AlgFarDetEvent: Event has no primary track or shower!" << endl; 
00086     }
00087 
00088   // Calculate other event properties
00089   // ================================
00090   MSG("FarDetEvent",Msg::kDebug) << " Setting event properties " << endl;
00091 
00092   int minplane(-999), maxplane(-999);
00093   double minz(-999.9), maxz(-999.9);
00094 
00095   for(Int_t i=0; i<AllStrips->GetEntries(); i++)
00096     {
00097       CandStripHandle* strip = (CandStripHandle*)(AllStrips->At(i));
00098     
00099       if( minplane<0 || strip->GetPlane()<minplane ){
00100         minplane = strip->GetPlane(); 
00101         minz = strip->GetZPos();
00102       }
00103  
00104       if( maxplane<0 || strip->GetPlane()>maxplane ){
00105         maxplane = strip->GetPlane();  
00106         maxz = strip->GetZPos();
00107       }
00108 
00109     }
00110 
00111   myevent.SetMinPlane(minplane); 
00112   myevent.SetMaxPlane(maxplane);
00113   myevent.SetMinZ(minz); 
00114   myevent.SetMaxZ(maxz);
00115 
00116   // Calculate event kinematics
00117   // ==========================
00118   MSG("FarDetEvent",Msg::kDebug) << " Calculating event kinematics " << endl;
00119 
00120   double EhadLinear(0.), EhadDeweight(0.);
00121 
00122   // (1) get MIP values for each shower strip
00123   // ========================================
00124   // input: primary shower, primary track
00125   // output: list of MIPs, track correction, shared planes 
00126 
00127   MSG("FarDetEvent",Msg::kDebug) << " (1) calculate MIP values for each strip " << endl;
00128 
00129   vector<double> TempShw;
00130   vector<double> TempTrk;
00131 
00132   this->GetShower( myevent.GetPrimaryShower(), myevent.GetPrimaryTrack(), TempShw, TempTrk );
00133 
00134   MSG("FarDetEvent",Msg::kDebug) << "    ... shower strips = " << TempShw.size() << endl;
00135 
00136   // (2) calculate pulse height and energy
00137   // =====================================
00138   // input: list of MIPs, track correction, shared planes, deweight factor
00139   // output: corrected pulse height, then shower energy
00140 
00141   MSG("FarDetEvent",Msg::kDebug) << " (2) calculate pulse height and energy " << endl;
00142 
00143   // linear pulse height and energy
00144   double linearfactor(1.0), linearph(0.0);
00145 
00146   this->GetPulseHeight( TempShw, TempTrk, linearph, linearfactor );
00147 
00148   EhadLinear = this->CalibrateShowerEnergy( linearph, linearfactor, ac );
00149 
00150   MSG("FarDetEvent",Msg::kDebug) << "    ... linear pulse height = " << linearph << endl;
00151   MSG("FarDetEvent",Msg::kDebug) << "        linear energy = " << EhadLinear << endl;
00152 
00153   // deweighted pulse height and energy
00154   double deweightfactor(1.0), deweightph(0.0);
00155 
00156   deweightfactor = this->GetDeWeightFactor( EhadLinear, ac );
00157 
00158   this->GetPulseHeight( TempShw, TempTrk, deweightph, deweightfactor );
00159 
00160   EhadDeweight = this->CalibrateShowerEnergy( deweightph, deweightfactor, ac );
00161 
00162   MSG("FarDetEvent",Msg::kDebug) << "    ... deweighted pulse height = " << deweightph << endl;
00163   MSG("FarDetEvent",Msg::kDebug) << "        deweighted energy = " << EhadDeweight << endl;
00164 
00165   // set reconstructed shower energy
00166   myevent.SetShwMipsLinear(linearph);            
00167   myevent.SetShwMipsDeweighted(deweightph);        
00168   myevent.SetShwEnergyGeVLinear(EhadLinear); 
00169   myevent.SetShwEnergyGeVDeweighted(EhadDeweight); 
00170   myevent.SetShwEnergyGeV(EhadLinear); 
00171 
00172   // (3) calculate overall energy
00173   // ============================
00174   double Enu(0.0);
00175   
00176   myevent.SetEnergy(Enu);
00177 
00178   return;
00179 }

void AlgFarDetEvent::SetShowerProperties ( FarDetEventHandle evt,
CandShowerHandle shw 
) [private]

Definition at line 250 of file AlgFarDetEvent.cxx.

References CandEventHandle::AddShower(), CandRecoHandle::GetDirCosU(), CandRecoHandle::GetDirCosV(), CandRecoHandle::GetDirCosZ(), FarDetEventHandle::GetMuReco(), FarDetEventHandle::GetShwDirCosU(), FarDetEventHandle::GetShwDirCosV(), FarDetEventHandle::GetShwDirCosZ(), FarDetEventHandle::GetShwReco(), FarDetEventHandle::GetVtxPlane(), CandRecoHandle::GetVtxPlane(), CandRecoHandle::GetVtxT(), FarDetEventHandle::GetVtxTime(), FarDetEventHandle::GetVtxU(), CandRecoHandle::GetVtxU(), FarDetEventHandle::GetVtxV(), CandRecoHandle::GetVtxV(), FarDetEventHandle::GetVtxZ(), CandRecoHandle::GetVtxZ(), Msg::kVerbose, MSG, CandEventHandle::SetPrimaryShower(), FarDetEventHandle::SetShwDirCosU(), FarDetEventHandle::SetShwDirCosV(), FarDetEventHandle::SetShwDirCosZ(), FarDetEventHandle::SetShwReco(), FarDetEventHandle::SetVtxPlane(), FarDetEventHandle::SetVtxTime(), FarDetEventHandle::SetVtxU(), FarDetEventHandle::SetVtxV(), and FarDetEventHandle::SetVtxZ().

Referenced by RunAlg().

00251 {
00252   MSG("FarDetEvent",Msg::kVerbose) << " AlgFarDetEvent::SetShowerProperties(....) " << endl;
00253 
00254   if( shw )
00255     {
00256       // set primary shower properties
00257       MSG("FarDetEvent",Msg::kVerbose) << "   (found primary shower) " << endl;
00258       evt->SetShwReco(1);
00259       evt->AddShower(shw);
00260       evt->SetPrimaryShower(shw);
00261       evt->SetShwDirCosU(shw->GetDirCosU()); 
00262       evt->SetShwDirCosV(shw->GetDirCosV()); 
00263       evt->SetShwDirCosZ(shw->GetDirCosZ());
00264       if( evt->GetMuReco()==0 )
00265         {
00266           evt->SetVtxTime(shw->GetVtxT());
00267           evt->SetVtxU(shw->GetVtxU()); 
00268           evt->SetVtxV(shw->GetVtxV()); 
00269           evt->SetVtxZ(shw->GetVtxZ());
00270           evt->SetVtxPlane(shw->GetVtxPlane()); 
00271         }
00272     }
00273 
00274   MSG("FarDetEvent",Msg::kVerbose) << endl
00275     << " Shower properties: " << endl
00276     << "  ShwReco = " << evt->GetShwReco() << endl
00277     << "  ShwDirCosU = " << evt->GetShwDirCosU() << endl
00278     << "  ShwDirCosV = " << evt->GetShwDirCosV() << endl
00279     << "  ShwDirCosZ = " << evt->GetShwDirCosZ() << endl   
00280     << "  (VtxTime = " << evt->GetVtxTime() << ")" << endl
00281     << "  (VtxU = " << evt->GetVtxU() << ")" << endl
00282     << "  (VtxV = " << evt->GetVtxV() << ")" << endl
00283     << "  (VtxZ = " << evt->GetVtxZ() << ")" << endl
00284     << "  (VtxPlane = " << evt->GetVtxPlane() << ")" << endl;
00285 
00286   return;
00287 }

void AlgFarDetEvent::SetTrackProperties ( FarDetEventHandle evt,
CandFitTrackHandle fit,
CandTrackHandle trk 
) [private]

Definition at line 184 of file AlgFarDetEvent.cxx.

References CandEventHandle::AddTrack(), CandRecoHandle::GetDirCosU(), CandRecoHandle::GetDirCosV(), CandRecoHandle::GetDirCosZ(), FarDetEventHandle::GetEMCharge(), CandFitTrackHandle::GetEMCharge(), CandTrackHandle::GetMomentum(), CandFitTrackHandle::GetMomentumCurve(), CandFitTrackHandle::GetMomentumRange(), FarDetEventHandle::GetMuDirCosU(), FarDetEventHandle::GetMuDirCosV(), FarDetEventHandle::GetMuDirCosZ(), FarDetEventHandle::GetMuMomentumCurve(), FarDetEventHandle::GetMuMomentumRange(), FarDetEventHandle::GetMuReco(), FarDetEventHandle::GetVtxPlane(), CandRecoHandle::GetVtxPlane(), CandRecoHandle::GetVtxT(), FarDetEventHandle::GetVtxTime(), FarDetEventHandle::GetVtxU(), CandRecoHandle::GetVtxU(), FarDetEventHandle::GetVtxV(), CandRecoHandle::GetVtxV(), FarDetEventHandle::GetVtxZ(), CandRecoHandle::GetVtxZ(), Msg::kVerbose, MSG, FarDetEventHandle::SetEMCharge(), FarDetEventHandle::SetMuDirCosU(), FarDetEventHandle::SetMuDirCosV(), FarDetEventHandle::SetMuDirCosZ(), FarDetEventHandle::SetMuMomentumCurve(), FarDetEventHandle::SetMuMomentumRange(), FarDetEventHandle::SetMuReco(), CandEventHandle::SetPrimaryTrack(), FarDetEventHandle::SetVtxPlane(), FarDetEventHandle::SetVtxTime(), FarDetEventHandle::SetVtxU(), FarDetEventHandle::SetVtxV(), and FarDetEventHandle::SetVtxZ().

Referenced by RunAlg().

00185 {
00186   MSG("FarDetEvent",Msg::kVerbose) << " AlgFarDetEvent::SetTrackProperties(....) " << endl;
00187 
00188   CandTrackHandle* track = 0; 
00189 
00190   double erange,ecurve,emcharge;
00191   // Use the fitted track if we have one...
00192   if(fit) 
00193     {
00194       track = fit;
00195       erange = fit->GetMomentumRange();
00196       ecurve = fabs(fit->GetMomentumCurve());
00197       emcharge = fit->GetEMCharge();
00198     }  
00199   // ...or else fall back to found track
00200   else if( trk ) 
00201     {
00202       track = trk;
00203       erange = trk->GetMomentum();
00204       ecurve = 0.0;
00205       emcharge = 0.0;
00206     }
00207   
00208   if( track )
00209     {
00210       // set primary track properties
00211       MSG("FarDetEvent",Msg::kVerbose) << "   (found primary track) " << endl;
00212       evt->AddTrack(track);
00213       evt->SetPrimaryTrack(track);
00214       evt->SetMuReco(1);
00215       evt->SetMuMomentumRange(erange);
00216       evt->SetMuMomentumCurve(ecurve);
00217       evt->SetEMCharge(emcharge);
00218       evt->SetMuDirCosU(track->GetDirCosU()); 
00219       evt->SetMuDirCosV(track->GetDirCosV()); 
00220       evt->SetMuDirCosZ(track->GetDirCosZ());
00221       evt->SetVtxTime(track->GetVtxT());
00222       evt->SetVtxU(track->GetVtxU()); 
00223       evt->SetVtxV(track->GetVtxV()); 
00224       evt->SetVtxZ(track->GetVtxZ());
00225       evt->SetVtxPlane(track->GetVtxPlane()); 
00226     }
00227 
00228   MSG("FarDetEvent",Msg::kVerbose) << endl
00229     << " Track Properties: " << endl
00230     << "  MuReco = " << evt->GetMuReco() << endl
00231     << "  MuMomentumRange = " << evt->GetMuMomentumRange() << endl
00232     << "  MuMomentumCurve = " << evt->GetMuMomentumCurve() << endl
00233     << "  EMCharge = " << evt->GetEMCharge() << endl
00234     << "  MuDirCosU = " << evt->GetMuDirCosU() << endl
00235     << "  MuDirCosV = " << evt->GetMuDirCosV() << endl
00236     << "  MuDirCosZ = " << evt->GetMuDirCosZ() << endl   
00237     << "  (VtxTime = " << evt->GetVtxTime() << ")" << endl
00238     << "  (VtxU = " << evt->GetVtxU() << ")" << endl
00239     << "  (VtxV = " << evt->GetVtxV() << ")" << endl
00240     << "  (VtxZ = " << evt->GetVtxZ() << ")" << endl
00241     << "  (VtxPlane = " << evt->GetVtxPlane() << ")" << endl;
00242 
00243   return;
00244 }

void AlgFarDetEvent::Trace ( const char *  c  )  const [virtual]

Reimplemented from AlgBase.

Definition at line 723 of file AlgFarDetEvent.cxx.

00724 {
00725 
00726 }


The documentation for this class was generated from the following files:
Generated on Thu Jul 10 22:52:16 2014 for loon by  doxygen 1.4.7