AlgReco.cxx

Go to the documentation of this file.

// $$
//
// AlgReco.cxx
//
// Begin_Html<img src="../../pedestrians.gif" align=center>
// <a href="../source_warning.html">Warning for beginners</a>.<br> 
//
// This is an abstract Algorithm reconstruction object.
// This class should and cannot be instantiated.
//
//
// Author:  R. Lee
//
// Also see <a href="../../root_crib/index.html">The ROOT Crib</a> and 
// <a href="../AlgReco.html"> AlgReco Classes</a> (part of
// <a href="../index.html">The MINOS Class User Guide</a>)End_Html

#include "Algorithm/AlgHandle.h"
#include "Calibrator/Calibrator.h"
#include "Conventions/PlaneView.h"
#include "RecoBase/AlgReco.h"
#include "RecoBase/CandRecoHandle.h"
#include "RecoBase/CandStripHandle.h"
#include "RecoBase/CandTrackHandle.h"
#include "RecoBase/CandShowerHandle.h"
#include "RecoBase/CandEventHandle.h"
#include "MessageService/MsgService.h"
#include "Conventions/Munits.h"

#include <cmath>

ClassImp(AlgReco)

//______________________________________________________________________
CVSID("$Id: AlgReco.cxx,v 1.14 2006/07/14 17:13:52 musser Exp $");

void AlgReco::Calibrate(CandRecoHandle *candreco)
{
  Calibrator& cal = Calibrator::Instance();

  Bool_t isEvent = candreco->InheritsFrom("CandEventHandle");
  const CandEventHandle *candevent = 0;
  if (isEvent) {
      
    candevent = dynamic_cast<const CandEventHandle*>(candreco);    
    for(Int_t i=0;i<=candevent->GetLastShower();i++){
  
      const CandShowerHandle *candshower = candevent->GetShower(i);
      TIter stripItr(candshower->GetDaughterIterator());
      while (CandStripHandle *strip = dynamic_cast<CandStripHandle*>(stripItr())) {
        FloatErr tpos=0.;
        int plane = strip->GetPlane();
        switch (strip->GetPlaneView()) {
        case PlaneView::kX: case PlaneView::kU:
          tpos = candshower->GetV(plane);
          tpos.SetError(1.0 * Munits::meter); // Since this next call doesn't work:
          //tpos.SetError( 0.5 * (candshower->GetMaxV(plane) - candshower->GetMinV(plane))  );
          break;
        case PlaneView::kY: case PlaneView::kV:
          tpos = candshower->GetU(plane);
          tpos.SetError(1.0 * Munits::meter); // Since this next call doesn't work:
          //tpos.SetError( 0.5 * (candshower->GetMaxU(plane) - candshower->GetMinU(plane))  );
          break;
        default:
          continue;
        }
        
        if (fabs(tpos)>4.) {
          MSG("RecoBase",Msg::kDebug) << "tpos = " << tpos << " outside detector" << endl;
          if (tpos>0.) {
            tpos = 4.;
          } else {
            tpos = -4.;
          }
        }
      
        FloatErr sigcorr;
        FloatErr sigmap, sigmip;
        PlexStripEndId seid;
        if (strip->GetNDigit(StripEnd::kNegative)>0) {
          seid    = strip->GetStripEndId(StripEnd::kNegative);
          sigcorr = strip->GetCharge(StripEnd::kNegative,CalDigitType::kSigCorr);
          sigmap =  cal.GetAttenCorrectedTpos(sigcorr, tpos, seid); 
          sigmip =  cal.GetMIP               (sigmap,        seid);
         
          candreco->CalibrateSigMapped(seid.GetEncoded(),sigmap);
          candreco->CalibrateMIP      (seid.GetEncoded(),sigmip);
        }
        if (strip->GetNDigit(StripEnd::kPositive)>0) {
          seid    = strip->GetStripEndId(StripEnd::kPositive);
          sigcorr = strip->GetCharge(StripEnd::kPositive,CalDigitType::kSigCorr);
          sigmap =  cal.GetAttenCorrectedTpos(sigcorr, tpos, seid); 
          sigmip =  cal.GetMIP               (sigmap,        seid);
        
          
          candreco->CalibrateSigMapped(seid.GetEncoded(),sigmap);
          candreco->CalibrateMIP      (seid.GetEncoded(),sigmip);
        }
      }
    }
    for(Int_t i=0;i<=candevent->GetLastTrack();i++){
      const CandTrackHandle *candtrack = candevent->GetTrack(i);
        
      TIter stripItr(candtrack->GetDaughterIterator());
      while (CandStripHandle *strip = dynamic_cast<CandStripHandle*>(stripItr())) {
        FloatErr tpos=0.;
        switch (strip->GetPlaneView()) {
        case PlaneView::kX: case PlaneView::kU:
          tpos = candtrack->GetV(strip->GetPlane());
          tpos.SetError(candtrack->GetTrackPointError(strip->GetPlane())); // Set error on tpos.
          break;
        case PlaneView::kY: case PlaneView::kV:
          tpos = candtrack->GetU(strip->GetPlane());
          tpos.SetError(candtrack->GetTrackPointError(strip->GetPlane())); // Set error on tpos.
          break;
        default:
          continue;
        }

        if (fabs(tpos)>4.) {
          MSG("RecoBase",Msg::kDebug) << "tpos = " << tpos << " outside detector" << endl;
          if (tpos>0.) {
            tpos = 4.;
          } else {
            tpos = -4.;
          }
        }
        
        FloatErr sigcorr;
        FloatErr sigmap, sigmip;
        PlexStripEndId seid;
        if (strip->GetNDigit(StripEnd::kNegative)>0) {
          seid    = strip->GetStripEndId(StripEnd::kNegative);
          sigcorr = strip->GetCharge(StripEnd::kNegative,CalDigitType::kSigCorr);
          sigmap =  cal.GetAttenCorrectedTpos(sigcorr, tpos, seid); 
          sigmip =  cal.GetMIP(               sigmap,        seid);
          
          candreco->CalibrateSigMapped(seid.GetEncoded(),sigmap);
          candreco->CalibrateMIP      (seid.GetEncoded(),sigmip);

        }
        if (strip->GetNDigit(StripEnd::kPositive)>0) {
          seid    = strip->GetStripEndId(StripEnd::kPositive);
          sigcorr = strip->GetCharge(StripEnd::kPositive,CalDigitType::kSigCorr);
          sigmap =  cal.GetAttenCorrectedTpos(sigcorr, tpos, seid); 
          sigmip =  cal.GetMIP               (sigmap,        seid);
          candreco->CalibrateSigMapped(seid.GetEncoded(),sigmap);
          candreco->CalibrateMIP      (seid.GetEncoded(),sigmip);
        }
      }
    }
  }

  else{
    TIter stripItr(candreco->GetDaughterIterator());
    Bool_t isTrack = candreco->InheritsFrom("CandTrackHandle");
    const CandTrackHandle *candtrack = 0;
    if (isTrack) {
      candtrack = dynamic_cast<const CandTrackHandle*>(candreco);
  }
    Bool_t isShower = candreco->InheritsFrom("CandShowerHandle");
    const CandShowerHandle *candshower = 0;
    if (isShower) {
      candshower = dynamic_cast<const CandShowerHandle*>(candreco);
    }
    Double_t vtxuvpos[2] = {candreco->GetVtxU(),candreco->GetVtxV()};
    while (CandStripHandle *strip = dynamic_cast<CandStripHandle*>(stripItr())) {
      Int_t iview=-1;
      switch (strip->GetPlaneView()) {
      case PlaneView::kX: case PlaneView::kU:
        iview = 1;
      break;
      case PlaneView::kY: case PlaneView::kV:
        iview = 0;
        break;
      default:
        continue;
      }
      FloatErr tpos = 0.;
      int plane = strip->GetPlane();
      // get position of hit along strip from orthogonal view info.
      if(candtrack){
        if (!iview) {
          tpos = candtrack->GetU(plane);
          tpos.SetError(candtrack->GetTrackPointError(plane)); // Set error on tpos.
        }
        else {
          tpos = candtrack->GetV(plane);
          tpos.SetError(candtrack->GetTrackPointError(plane)); // Set error on tpos.
        }
      }
      else if (candshower){
        if (!iview) {
          tpos = candshower->GetU(plane);
          tpos.SetError(1.0*Munits::meter); // Since this other call doesn't work:
          //tpos.SetError( 0.5 * (candshower->GetMaxU(plane) - candshower->GetMinU(plane))  );
        }
        else {
          tpos = candshower->GetV(plane);
          tpos.SetError(1.0*Munits::meter); // Since this other call doesn't work:
          //tpos.SetError( 0.5 * (candshower->GetMaxV(plane) - candshower->GetMinV(plane))  );
        }
      }
      else {
        tpos = vtxuvpos[iview];
      }
      
      if (fabs(tpos)>4.) {
        MSG("RecoBase",Msg::kDebug) << "tpos = " << tpos << " outside detector" << endl;
        if (tpos>0.) {
          tpos = 4.;
        } else {
          tpos = -4.;
        }
      }
      
      
      FloatErr sigcorr;
      FloatErr sigmap, sigmip;
      PlexStripEndId seid;
      if (strip->GetNDigit(StripEnd::kNegative)>0) {
        seid    = strip->GetStripEndId(StripEnd::kNegative);
        sigcorr = strip->GetCharge(StripEnd::kNegative,CalDigitType::kSigCorr);
        sigmap =  cal.GetAttenCorrectedTpos(sigcorr, tpos, seid); 
        sigmip =  cal.GetMIP(               sigmap,        seid);
        
        candreco->CalibrateSigMapped(seid.GetEncoded(),sigmap);
        candreco->CalibrateMIP      (seid.GetEncoded(),sigmip);
      }
      if (strip->GetNDigit(StripEnd::kPositive)>0) {
        seid    = strip->GetStripEndId(StripEnd::kPositive);
        sigcorr = strip->GetCharge(StripEnd::kPositive,CalDigitType::kSigCorr);
        sigmap =  cal.GetAttenCorrectedTpos(sigcorr, tpos, seid); 
        sigmip =  cal.GetMIP               (sigmap,        seid);
        
        candreco->CalibrateSigMapped(seid.GetEncoded(),sigmap);
        candreco->CalibrateMIP      (seid.GetEncoded(),sigmip); 
      }
    }
  }
}

---