FGDHybridTrackEffSystematics.cxx

#include "FGDHybridTrackEffSystematics.hxx"
#include "ND280AnalysisUtils.hxx"
#include "EventBoxTracker.hxx"
#include "SystematicUtils.hxx"
#include "ToyBoxTracker.hxx"
#include "SystId.hxx"
#include "EventBoxId.hxx"


//********************************************************************
FGDHybridTrackEffSystematics::FGDHybridTrackEffSystematics(bool comp):EventWeightBase(1){
//********************************************************************

  _computecounters=comp;
  _fgd1 = new BinnedParams("FGD1HybridTrackEff_p6B",BinnedParams::k2D_EFF_ASSYMMETRIC);
  _fgd2 = new BinnedParams("FGD2HybridTrackEff_p6B",BinnedParams::k2D_EFF_ASSYMMETRIC);
  SetNParameters(2*(_fgd1->GetNBins()+_fgd2->GetNBins()));

  if(_computecounters){
    _fgd1->InitializeEfficiencyCounter();
    _fgd2->InitializeEfficiencyCounter();
  }
}

//********************************************************************
Weight_h FGDHybridTrackEffSystematics::ComputeWeight(const ToyExperiment& toy, const AnaEventC& event, const ToyBoxB& boxB){
//********************************************************************

  Weight_h eventWeight=1; 
  BinnedParamsParams params;

  // Cast the ToyBox to the appropriate type
  const ToyBoxTracker& box = *static_cast<const ToyBoxTracker*>(&boxB); 

  // Get the SystBox for this event, and the appropriate selection and branch
  SystBoxB* SystBox = GetSystBox(event,box.SelectionEnabledIndex,box.SuccessfulBranch);

  if(_computecounters){
    _fgd1->InitializeEfficiencyCounter();
    _fgd2->InitializeEfficiencyCounter();
  }
  // Get the direction of the main track
  AnaTrackB *maintrack = box.MainTrack; 
  if(!maintrack)  return eventWeight;

  for (Int_t itrue=0;itrue< SystBox->nRelevantTrueObjects; itrue++){      
    AnaTrueParticleB* trueTrack = static_cast<AnaTrueParticleB*>(SystBox->RelevantTrueObjects[itrue]);            

    // TODO by now the reco track associated to the true is the first one found, we should use instead a criteria to decide.
    // by now since it's very rare it's like that.
    // retrieve the reconstructed track associated
    AnaTrackB* recoTrack = static_cast<AnaTrackB*>(SystBox->RelevantTrueObjectsReco[itrue]);

    // Get the cos (theta), where theta is the angle between the main track and the FGD track
    Float_t costheta = 0;
    if (recoTrack)
      costheta = anaUtils::ScalarProduct(maintrack->DirectionStart, recoTrack->DirectionStart);
    else
      costheta = anaUtils::ScalarProduct(maintrack->DirectionStart, trueTrack->Direction);

    Float_t pdg=(Float_t)trueTrack->PDG;
    bool found = (recoTrack);
    
    // Get the FGD tracking efficiency
    Int_t index;
    if(box.DetectorFV==SubDetId::kFGD1){
      if(!_fgd1->GetBinValues(fabs(pdg),costheta, params, index)) continue;
      // save the index of the varied parameters into the systematics box
      if(_computecounters)
        _fgd1->UpdateEfficiencyCounter(index,found);    
    }
    else if(box.DetectorFV==SubDetId::kFGD2){
      if(!_fgd2->GetBinValues(fabs(pdg),costheta, params, index)) continue;
      if(_computecounters)
        _fgd2->UpdateEfficiencyCounter(index,found);    
      
    }
    else return eventWeight;
    
    eventWeight *= systUtils::ComputeEffLikeWeight(true, toy, GetIndex(), index, params);
  }

  return eventWeight;

}

//********************************************************************
bool FGDHybridTrackEffSystematics::IsRelevantTrueObject(const AnaEventC& event, const AnaTrueObjectC& trueTrack) const{
//********************************************************************

  EventBoxB* EventBox = event.EventBoxes[EventBoxId::kEventBoxTracker];

  //we should not take into account the tracks that are matched to the tpc,
  //this is done in tpc-fgd matching efficiency  
  for (Int_t i=0;i<EventBox->nRecObjectsInGroup[EventBoxTracker::kTracksWithTPC];i++){
    AnaRecObjectC* track = EventBox->RecObjectsInGroup[EventBoxTracker::kTracksWithTPC][i];
    if (CheckTrueRecoAssociation(trueTrack, *track))
      return false;
  }

  return true;
}

//********************************************************************
Int_t FGDHybridTrackEffSystematics::GetRelevantRecObjectGroups(const SelectionBase& sel, Int_t ibranch, Int_t* IDs) const{
//********************************************************************

  SubDetId_h det = sel.GetDetectorFV(ibranch);

  if      (det == SubDetId::kFGD1){
    IDs[0] = EventBoxTracker::kTracksWithFGD1AndNoTPC; 
    return 1;
  }
  else if (det == SubDetId::kFGD2){
    IDs[0] = EventBoxTracker::kTracksWithFGD2AndNoTPC;
    return 1; 
  }
  else if (det == SubDetId::kFGD){
    IDs[0] = EventBoxTracker::kTracksWithFGD1AndNoTPC;
    IDs[1] = EventBoxTracker::kTracksWithFGD2AndNoTPC;
    return 2; 
  }

  return 0;
}

//********************************************************************
Int_t FGDHybridTrackEffSystematics::GetRelevantTrueObjectGroups(const SelectionBase& sel, Int_t ibranch, Int_t* IDs) const{
//********************************************************************

  SubDetId_h det = sel.GetDetectorFV(ibranch);

  if      (det == SubDetId::kFGD1){
    IDs[0] = EventBoxTracker::kTrueParticlesChargedInFGD1AndNoTPCInBunch; 
    return 1; 
  }
  else if (det == SubDetId::kFGD2){
    IDs[0] = EventBoxTracker::kTrueParticlesChargedInFGD2AndNoTPCInBunch;
    return 1; 
  }
  else if (det == SubDetId::kFGD){
    IDs[0] = EventBoxTracker::kTrueParticlesChargedInFGD1AndNoTPCInBunch;
    IDs[1] = EventBoxTracker::kTrueParticlesChargedInFGD2AndNoTPCInBunch;
    return 2; 
  }

  return 0;
}