VtxClusterList.cxx

Go to the documentation of this file.

//                      VtxClusterList.cxx  -  description
//                             -------------------
//    begin                : Thu June 3 2004
//    copyright            : (C) 2004 by Joshua Boehm
//    email                : boehm@physics.harvard.edu
//    Last Modified        : August 10, 2004
//    VtxClusterList is a reconstruction class used by VertexFinder
//    to organize and determine the vertex position in z.

#include "VtxClusterList.h"
#include "MessageService/MsgService.h"

//Inline Functions Defined in VtxClusterList.h
//   Int_t Size() {return list.size();};  //Returns the number of groups
//   Int_t Size(int i) {return list[i].Size();}; Returns size of cluster i
//   void Clear() {list.clear(); };    //Clears the group listing

                                                                                
CVSID("$Id: VtxClusterList.cxx,v 1.3 2006/12/01 20:30:53 rhatcher Exp $");


Float_t VtxClusterList::GetTotalEventEnergy()
{
  Float_t total = 0;
  for(Int_t i = 0; i < Size(); i++) total += list[i].GetTotalEnergy();

  return total;
}

VtxClusterList::VtxClusterList(vector<Float_t> &lenepl, Int_t N_PLANES, Float_t kNoiseCut)
{
  bool inGroup = false;
  Float_t binEnergy = 0;
  Float_t prevBinEnergy = -1;

  for(Int_t plane = 0; plane < N_PLANES; plane++)
  {
     binEnergy = lenepl[plane];
     if(binEnergy > kNoiseCut)
     {
        if(inGroup){
           list[TMath::Max(Size()-1, 0)].AddPlane(plane, binEnergy);
        }else{
            inGroup = true;
            VtxCluster newCluster(plane, binEnergy);
            list.push_back(newCluster);
        }
     }else{
        if(inGroup && prevBinEnergy < kNoiseCut) inGroup = false;
     }
       
     prevBinEnergy = binEnergy;
  }

  if(Size() == 0){
    MSG("VtxClusterList", Msg::kError)
          <<"No Clusters could be formed"<<endl;
  
    for(int i = 0; i < N_PLANES; i++)
    {
       MSG("VtxClusterList", Msg::kError)
           <<"Plane "<<i<<" of "<<N_PLANES<<": "<<lenepl[i]<<endl;
    }
  }

}

VtxClusterList::~VtxClusterList()
{
  Clear();
}

void VtxClusterList::Report()
{
  for(Int_t i = 0; i < Size(); i++){
    cout<<i<<"   "<<list[i].Size()<<"  "<<list[i].GetFirstPlane()
        <<"    "<<list[i].GetLastPlane()<<endl;
    for(Int_t j = 0; j < list[i].Size(); j++)
        cout<<j<<"  "<<list[i].GetPlane(j)<<"  "<<list[i].GetEnergy(j)<<endl;
  }
}

Int_t VtxClusterList::RisingEdgeCandidate()
{
   //First go and check for strong rising edges
   Int_t candidate = -1;
   bool found = false;
 
   for(Int_t i = 0; i < Size() && !found; i++)
   {
      if(list[i].HasRisingEdge())
      {
         if(list[i].Size() > 15 || list[i].GetTotalEnergy() > 200)//i0.5*GetTotalEventEnergy())
         {
           candidate = i;
           found = true;
         }

         Float_t subTotal = 0;
         for(Int_t j = 0; j < i; j++)
              subTotal += list[j].GetTotalEnergy();
 
         if(subTotal/GetTotalEventEnergy() > 0.5){
             candidate = -1;
             found = true;  //What i have found though is this isnt
                            //   meant to be a rising edge system
             continue;
         }
         if(candidate == -1) candidate = i;
         else{
           if(list[i].Size() > 4 && 
               list[i].GetTotalEnergy() > list[candidate].GetTotalEnergy())
              candidate = i;
         }
      }
   }

   return candidate;
}

Int_t VtxClusterList::Process()
{
   if(Size() == 0) return -1;

   bool debug = false;
   Int_t maxGroupSize = 0;
   Int_t candidate = -1;
   Int_t candPlane = -1;
   Int_t startGuess = -1;
   for(Int_t i = 0; i < Size(); i++)
     maxGroupSize = TMath::Max(maxGroupSize, list[i].Size());

   candidate = RisingEdgeCandidate();
   
   if(debug) cout<<candidate<<endl;
     
   if(candidate != -1){
     startGuess = list[candidate].GetEdgeTrigger();
   }else {
     candidate = 0;
     for(Int_t i = 0; i < Size(); i++)
     {
        Float_t recent = list[i].GetEnergy(list[i].GetNCTrigger());
        Float_t old = list[candidate].GetEnergy(list[candidate].GetNCTrigger());
        if(debug) cout<<old<<"    "<<recent<<"    "<<endl;        
   
        Float_t subTotal = 0;
        for(Int_t j = 0; j < i; j++)
           subTotal += list[j].GetTotalEnergy();

        //We cant have  70% of the energy backscattered
        if(recent > 2*old && subTotal/GetTotalEventEnergy() < 0.5)
           candidate = i;
     }
     startGuess = list[candidate].GetNCTrigger();
   }
 
   if(debug) cout<<"start  "<<startGuess<<endl;

   if(startGuess > 0)
   {
     candPlane = list[candidate].GetPlane(startGuess);
    if(list[candidate].Size() > 3) {
     //If there is a blip then the rising edges, ignore that repeat as necessary
     //     rewinding
     Int_t i = startGuess;
     if(i > 0)
     {
        if(i > 0 && list[candidate].GetPlane(i)-list[candidate].GetPlane(i-1) == 1
                 && list[candidate].GetEnergy(i) > list[candidate].GetEnergy(i-1))
                i--;

        candPlane = list[candidate].GetPlane(i);
     }
    }
   }else
      candPlane = list[candidate].GetFirstPlane();
 
   if(debug) cout<<candidate<<"  "<<candPlane<<endl;
   if(candPlane == -1) candPlane = list[candidate].GetFirstPlane();
   if(debug) cout<<candidate<<"  "<<candPlane<<endl;
   return candPlane;
}

---