////////////////////////////////////////////////////////////////////////
// $Id: AlgThruMuon.cxx,v 1.24 2007/02/04 06:30:13 rhatcher Exp $
//
// AlgThruMuon
//
// Begin_Html<img src="../../pedestrians.gif" align=center>
// <a href="../source_warning.html">Warning for beginners</a>.<br> 
//
// This is an Algorithm class to set the two CandStraightClusters
// associated with the CandThruMuon and fit a straight line to the digit
// pairs.
//
// Author:  P.S. Miyagawa 9/2000
//
// Also see <a href="../../root_crib/index.html">The ROOT Crib</a> and 
// <a href="../ALG_Classes.html"> Algorithm Classes</a> (part of
// <a href="../index.html">The MINOS Class User Guide</a>)End_Html
////////////////////////////////////////////////////////////////////////

#include <cassert>

#include "Algorithm/AlgFactory.h"
#include "Algorithm/AlgHandle.h"
#include "Candidate/CandContext.h"
#include "MessageService/MsgService.h"
#include "UgliGeometry/UgliGeomHandle.h"
#include "BubbleSpeak/AlgThruMuon.h"
#include "BubbleSpeak/CandDigiPairHandle.h"
#include "BubbleSpeak/CandStraightCluster.h"
#include "BubbleSpeak/CandStraightClusterHandle.h"
#include "BubbleSpeak/CandThruMuonHandle.h"
#include "TMath.h"

ClassImp(AlgThruMuon)

//......................................................................

CVSID("$Id: AlgThruMuon.cxx,v 1.24 2007/02/04 06:30:13 rhatcher Exp $");

//......................................................................

AlgThruMuon::AlgThruMuon()
{
//
//  Purpose:    Default constructor.
//
//  Arguments:  n/a
//
//  Return:     n/a
//
}

//......................................................................

AlgThruMuon::~AlgThruMuon()
{
//
//  Purpose:    Default destructor.
//
//  Arguments:  n/a
//
//  Return:     n/a
//
}

//......................................................................

void AlgThruMuon::RunAlg(AlgConfig & /* ac */, CandHandle &ch, CandContext &cx)
{
//
//  Purpose:  Fills the CandThruMuon with the two matched
//            CandStraightClusters and fits a straight line to the digit
//            pairs.
//
//  Argument:
//    ac        in    AlgConfig (not used).
//    ch        in    Handle to the CandThruMuon to fill.
//    cx        in    CandContext containing a TObjArray containing the
//                    matched pair of CandStraightClusterHandles.
//
//  Return:  n/a
//

   MSG("BubAlg", Msg::kVerbose)
      << "Starting AlgThruMuon::RunAlg()" << endl;

   assert(ch.InheritsFrom("CandThruMuonHandle"));
   CandThruMuonHandle &cmh = dynamic_cast<CandThruMuonHandle &>(ch);

// Check for TObjArray input.
   assert(cx.GetDataIn());
   assert(cx.GetDataIn()->InheritsFrom("TObjArray"));
   const TObjArray *mpair =
                    dynamic_cast<const TObjArray *>(cx.GetDataIn());

// Loop over elements of pair to retrieve/create U- and V-clusters.
   CandStraightClusterHandle cch[2];   // 0=U, 1=V
   for (Int_t i=0; i<2; i++) {

// Check whether already have CandStraightCluster.
      TObject *tobj = mpair->At(i);
      if (tobj->InheritsFrom("CandStraightClusterHandle")) {
         MSG("BubAlg", Msg::kVerbose)
            << "Adding StraightCluster." << endl;
         cch[i] = *(dynamic_cast<CandStraightClusterHandle *>(tobj));
      }

// Have to create CandStraightCluster.
      else if (tobj->InheritsFrom("CandMSTClusterHandle")
               || tobj->InheritsFrom("TObjArray")) {
         MSG("BubAlg", Msg::kVerbose)
            << "Creating StraightCluster." << endl;
         AlgFactory &af = AlgFactory::GetInstance();
         AlgHandle ach =af.GetAlgHandle("AlgStraightCluster","default");
         CandContext cxx(this, cx.GetMom());
         cxx.SetDataIn(tobj);
         cxx.SetCandRecord(cx.GetCandRecord());
         
         cch[i] = CandStraightCluster::MakeCandidate(ach, cxx);
      }

// Wrong input.
      else assert(0);
   }
   MSG("BubAlg", Msg::kVerbose)
      << "Created all StraightClusters" << endl;

// Set U- and V-clusters.
   MSG("BubAlg", Msg::kVerbose) << "Saving StraightClusters." << endl;
   cch[0].SetCandRecord(cx.GetCandRecord());
   cmh.SetClusterU(cch[0]);
   cch[1].SetCandRecord(cx.GetCandRecord());
   cmh.SetClusterV(cch[1]);

// Save digits in each cluster.
   MSG("BubAlg", Msg::kVerbose) << "Saving DigiPairs." << endl;
   for (int i=0; i<2; i++) {
      TIter chhItr(cch[i].GetDaughterIterator());
      while (CandDigiPairHandle *chh =
                dynamic_cast<CandDigiPairHandle *>(chhItr())) {
         cmh.AddDaughterLink(*chh);
      }
   }
   MSG("BubAlg", Msg::kVerbose)
      << "Saved clusters and digits in CandTrack" << endl;

// Select vertex point to be highest point.
   assert(cch[0].GetVldContext());
   UgliGeomHandle ugh(*cch[0].GetVldContext());

   // Set begin point.
   MSG("BubAlg", Msg::kVerbose) << "Setting begin point" << endl;
   Float_t uval = cch[0].GetZBeg();
   Float_t vval = cch[1].GetZBeg();
   Float_t zbeg = (uval < vval ? uval : vval);
   Float_t ubeg = cch[0].GetFitTPos(zbeg);
   Float_t vbeg = cch[1].GetFitTPos(zbeg);
   Float_t xbeg, ybeg;
   ugh.uv2xy(ubeg, vbeg, xbeg, ybeg);
   MSG("BubAlg", Msg::kVerbose) << "Set begin point" << endl;

   // Set end point.
   MSG("BubAlg", Msg::kVerbose) << "Setting end point" << endl;
   uval = cch[0].GetZEnd();
   vval = cch[1].GetZEnd();
   Float_t zend = (uval > vval ? uval : vval);
   Float_t uend = cch[0].GetFitTPos(zend);
   Float_t vend = cch[1].GetFitTPos(zend);
   Float_t xend, yend;
   ugh.uv2xy(uend, vend, xend, yend);
   MSG("BubAlg", Msg::kVerbose) << "Set end point" << endl;

   // Assume muon is down-going.
   Float_t dirf;
   if ((ybeg > yend) || (cmh.GetBegPlane() <= 2)) {
      cmh.SetVtxU(ubeg);
      cmh.SetVtxV(vbeg);
      cmh.SetVtxZ(zbeg);
      dirf = 1;
   }
   else {
      cmh.SetVtxU(uend);
      cmh.SetVtxV(vend);
      cmh.SetVtxZ(zend);
      dirf = -1;
   }

   // Set vertex plane.
   cmh.SetVtxPlane(dirf > 0 ? cmh.GetBegPlane() : cmh.GetEndPlane());

   // Set vertex time.
   cmh.SetVtxT(cmh.GetT(cmh.GetVtxPlane()));
   MSG("BubAlg", Msg::kVerbose) << "Set vertex info" << endl;

// Set direction cosines.
   Float_t du = dirf;
   Float_t uslope = cch[0].GetFitSlope();
   if (cch[0].GetFitMode()) du /= (uslope!=0. ? uslope : 1e-12);
   else du *= uslope;
   Float_t dv = dirf;
   Float_t vslope = cch[1].GetFitSlope();
   if (cch[1].GetFitMode()) dv /= (vslope!=0. ? vslope : 1e-12);
   else dv *= vslope;
   Float_t dz = dirf;
   Float_t ds = TMath::Sqrt(du*du + dv*dv + dz*dz);
   cmh.SetDirCosU(du/ds);
   cmh.SetDirCosV(dv/ds);
   cmh.SetDirCosZ(dz/ds);
   MSG("BubAlg", Msg::kVerbose) << "Set direction info" << endl;

// Retrieve cosmic counts, if any.
   TObject *tobj = mpair->At(2);
   if (tobj && tobj->InheritsFrom("TObjArray")) {
      TIter chhItr(dynamic_cast<TObjArray *>(tobj));
      while (CandDigiPairHandle *chh =
                dynamic_cast<CandDigiPairHandle *>(chhItr())) {
         cmh.AddDaughterLink(*chh);
      }
   }

// Set CandSlice.
   cmh.SetCandSlice(cch[0].GetCandSlice());

// Set u-, v-, z-coordinates for each plane.
   this->SetUVZT(&cmh);

// Calibrate.
   this->Calibrate(&cmh);

   MSG("BubAlg", Msg::kVerbose)
      << "Completed AlgThruMuon::RunAlg()" << endl;
}

//......................................................................

void AlgThruMuon::Trace(const char *c) const
{
//
//  Purpose:  Trace the AlgThruMuon.
//
//  Arguments:
//    c          in    String tag for the trace.
//
//  Return:   n/a
//

  MSG("BubCand", Msg::kDebug)
     << "**********Begin AlgThruMuon::Trace(\"" << c << "\")" << endl
     << "**********End AlgThruMuon::Trace(\"" << c << "\")" << endl;
}