StMuDSTMaker/COMMON/StMuEvent.cxx

/***************************************************************************
 *
 * $Id: StMuEvent.cxx,v 1.26 2011/08/30 14:41:26 fisyak Exp $
 * Author: Frank Laue, BNL, laue@bnl.gov
 *
 ***************************************************************************/

#include <string.h> 
#include "StEvent/StEvent.h" 
#include "StEvent/StEventTypes.h" 
#include "StEvent/StEventSummary.h" 
#include "StEvent/StEventInfo.h" 
#include "StEvent/StDetectorState.h" 

#include "StEventUtilities/StuRefMult.hh"
#include "StEventUtilities/StuFtpcRefMult.hh"

#include "StarClassLibrary/SystemOfUnits.h"
#include "StarClassLibrary/StTimer.hh"

#include "StMuException.hh"
#include "StMuEvent.h"
#include "StMuTrack.h"
#include "StMuL3EventSummary.h"
#include "StMuCut.h"
#include "StMuDebug.h"
#include "StMuDst.h"
#include "StMuPrimaryVertex.h"
#include "StMuBTofHit.h"

#include "TClonesArray.h"
#include "TObject.h"
#include "TClass.h"

ClassImp(StMuEvent)

//-----------------------------------------------------------------------
//-----------------------------------------------------------------------
//-----------------------------------------------------------------------
  StMuEvent::StMuEvent() : mTriggerData(0), mPrimaryVertexError(-999.,-999.,-999) { 
  DEBUGMESSAGE("");
  int n = (char*)mReactionPlanePtWgt - (char*)&mRefMultPos+sizeof(mReactionPlanePtWgt);
  memset(&mRefMultPos,0,n);
}
//-----------------------------------------------------------------------
//-----------------------------------------------------------------------
//-----------------------------------------------------------------------
StMuEvent::StMuEvent(const StEvent* event) : mPrimaryVertexError(-999.,-999.,-999.) { 
  try {
    fill(event);
  }
  catch (StMuException e) {
    throw e;
  }
}
//-----------------------------------------------------------------------
//-----------------------------------------------------------------------
//-----------------------------------------------------------------------
StMuEvent::~StMuEvent(){
  DEBUGMESSAGE("");
}
//-----------------------------------------------------------------------
//-----------------------------------------------------------------------
//-----------------------------------------------------------------------
void StMuEvent::clear(){
  DEBUGMESSAGE1("");
}
//-----------------------------------------------------------------------
//-----------------------------------------------------------------------
//-----------------------------------------------------------------------
void StMuEvent::fill(const StEvent* event){
  DEBUGMESSAGE("");
  if ( !event ) throw StMuExceptionNullPointer("no StEvent",__PRETTYF__);
  if ( !event->info() ) throw StMuExceptionNullPointer("no event info",__PRETTYF__);
  if ( !event->runInfo() ) throw StMuExceptionNullPointer("no run info",__PRETTYF__);
  if ( !event->summary() ) throw StMuExceptionNullPointer("no event summary",__PRETTYF__);
//  if (event->numberOfPrimaryVertices() != 1 ) throw StMuExceptionBadValue("!= 1 primary vertex");


  mRunInfo = *event->runInfo();
  mEventInfo = *event->info();
  mEventSummary = *event->summary();
  const StPrimaryVertex *p_vtx=event->primaryVertex();
  if (p_vtx) {
    mPrimaryVertexError=p_vtx->positionError();
  }

  if ( !event->triggerDetectorCollection() ) {
    DEBUGVALUE2(event->type());
    DEBUGVALUE2(event->info()->time());
    DEBUGMESSAGE2("no trigger detector collection, creating dummy");
  }
  else {
    mVpdTriggerDetector = event->triggerDetectorCollection()->vpd();
    mMtdTriggerDetector = event->triggerDetectorCollection()->mtd();
    mCtbTriggerDetector = event->triggerDetectorCollection()->ctb();
    mZdcTriggerDetector = event->triggerDetectorCollection()->zdc();
    mBbcTriggerDetector = event->triggerDetectorCollection()->bbc();
    mEmcTriggerDetector = event->triggerDetectorCollection()->emc();
    mFpdTriggerDetector = event->triggerDetectorCollection()->fpd();
    mFmsTriggerDetector = event->triggerDetectorCollection()->fms();
  }

  if (event->fpdCollection())
    mFpdCollection = *event->fpdCollection();
  if (event->l0Trigger())
    mL0Trigger = *event->l0Trigger();
  mL3EventSummary.fill(event);

  mTriggerIdCollection.fill( event->triggerIdCollection() ); // pointer check done inside


  mRefMultPos = uncorrectedNumberOfPositivePrimaries(*event);
  mRefMultNeg = uncorrectedNumberOfNegativePrimaries(*event); 
  mRefMultFtpcEast = uncorrectedNumberOfFtpcEastPrimaries(*event);
  mRefMultFtpcWest = uncorrectedNumberOfFtpcWestPrimaries(*event); 

  if (event->triggerData()) // MC has no triggerData
    mL2Result.Set(event->triggerData()->l2ResultLength(),(const Int_t*) event->triggerData()->l2Result());

  // calibrated vpd for TOF - X.Dong
  mVpdEast = mVpdWest = 0;
  mVpdTstart = mVpdTdiff = 0.;
  if (event->tofCollection()) {
    mVpdEast = event->tofCollection()->vpdEast();
    mVpdWest = event->tofCollection()->vpdWest();
    mVpdTstart = event->tofCollection()->tstart();
    mVpdTdiff = event->tofCollection()->tdiff();
    mVpdVz = event->tofCollection()->vzVpd();
  }
  // trigger data
  mTriggerData = const_cast<StTriggerData*>(event->triggerData());  
  if(mTriggerData!=0) mTriggerData->setDebug(0);
} 

unsigned short StMuEvent::refMultPos(int vtx_id) {
  // Check old format (no StMuPrimaryVertex stored)
  if (StMuDst::numberOfPrimaryVertices()==0 && (vtx_id == 0 || vtx_id == -1))
     return mRefMultPos;
  if (vtx_id == -1)
     vtx_id = StMuDst::currentVertexIndex();
  if (StMuDst::primaryVertex(vtx_id))
     return StMuDst::primaryVertex(vtx_id)->refMultPos();
  return 0;
}

unsigned short StMuEvent::refMultNeg(int vtx_id) {
  // Check old format (no StMuPrimaryVertex stored)
  if (StMuDst::numberOfPrimaryVertices()==0 && (vtx_id == 0 || vtx_id == -1))
     return mRefMultNeg;
  if (vtx_id == -1)
     vtx_id = StMuDst::currentVertexIndex();
  if (StMuDst::primaryVertex(vtx_id))
     return StMuDst::primaryVertex(vtx_id)->refMultNeg();
  return 0;
}

unsigned short StMuEvent::refMult(int vtx_id) {return refMultPos(vtx_id)+refMultNeg(vtx_id);}

unsigned short StMuEvent::refMultFtpcEast(int vtx_id) {
  // Check old format (no StMuPrimaryVertex stored)
  if (StMuDst::numberOfPrimaryVertices()==0 && (vtx_id == 0 || vtx_id == -1))
     return mRefMultFtpcEast;
  if (vtx_id == -1)
     vtx_id = StMuDst::currentVertexIndex();
  if (StMuDst::primaryVertex(vtx_id))
     return StMuDst::primaryVertex(vtx_id)->refMultFtpcEast();
  return 0;
}

unsigned short StMuEvent::refMultFtpcWest(int vtx_id) {
  // Check old format (no StMuPrimaryVertex stored)
  if (StMuDst::numberOfPrimaryVertices()==0 && (vtx_id == 0 || vtx_id == -1))
     return mRefMultFtpcWest;
  if (vtx_id == -1)
     vtx_id = StMuDst::currentVertexIndex();
  if (StMuDst::primaryVertex(vtx_id))
     return StMuDst::primaryVertex(vtx_id)->refMultFtpcWest();
  return 0;
}

unsigned short StMuEvent::refMultFtpc(int vtx_id) {return refMultFtpcEast(vtx_id)+refMultFtpcWest(vtx_id);}

StThreeVectorF StMuEvent::primaryVertexPosition(int vtx_id) { 
        StThreeVectorF vz(-999,-999,-999);
  // Check old format (no StMuPrimaryVertex stored)
  if (StMuDst::numberOfPrimaryVertices()==0 && (vtx_id == 0 || vtx_id == -1)){
         if(fabs(mEventSummary.primaryVertexPosition().x()) < 1.e-5 && fabs(mEventSummary.primaryVertexPosition().y()) < 1.e-5 && fabs(mEventSummary.primaryVertexPosition().z()) < 1.e-5) return vz;   
         else return mEventSummary.primaryVertexPosition();
        }
  if (vtx_id == -1)
     vtx_id = StMuDst::currentVertexIndex();
  if (StMuDst::primaryVertex(vtx_id))
     return StMuDst::primaryVertex(vtx_id)->position();
  return vz;
}
StThreeVectorF StMuEvent::primaryVertexErrors(int vtx_id) {
        StThreeVectorF vz(-999,-999,-999);
  // Check old format (no StMuPrimaryVertex stored)
  if (StMuDst::numberOfPrimaryVertices()==0 && (vtx_id == 0 || vtx_id == -1)){
         if(fabs(mEventSummary.primaryVertexPosition().x()) < 1.e-5 && fabs(mEventSummary.primaryVertexPosition().y()) < 1.e-5 && fabs(mEventSummary.primaryVertexPosition().z()) < 1.e-5) return vz;   
         else return mPrimaryVertexError;
        }
  if (vtx_id == -1)
     vtx_id = StMuDst::currentVertexIndex();
  if (StMuDst::primaryVertex(vtx_id))
     return StMuDst::primaryVertex(vtx_id)->posError();
  return vz;
}

unsigned short StMuEvent::grefmult(int vtx_id){
    unsigned short grefmult = 0;
        StMuTrack *glob;
        //For old MuDsts where there was one vertex per event
        if (StMuDst::numberOfPrimaryVertices()==0 && (vtx_id == 0 || vtx_id == -1)){
                if(!(fabs(mEventSummary.primaryVertexPosition().x()) < 1.e-5 && fabs(mEventSummary.primaryVertexPosition().y()) < 1.e-5 && fabs(mEventSummary.primaryVertexPosition().z()) < 1.e-5)){   
                        for (int i=0;i<StMuDst::globalTracks()->GetEntries();i++){
                                glob = StMuDst::globalTracks(i);
                                if (fabs(glob->eta()) <  0.5 && fabs(glob->dcaGlobal().mag()) < 3 && glob->nHitsFit(kTpcId) >= 10) grefmult++;            
                        }
                        return grefmult;
                }       
                else return 0;
        }       

        if (vtx_id == -1)
                vtx_id = StMuDst::currentVertexIndex();

        if (StMuDst::primaryVertex(vtx_id)){
        for (int i=0;i<StMuDst::globalTracks()->GetEntries();i++){
                        glob = StMuDst::globalTracks(i);
                        if (fabs(glob->eta()) <  0.5 && fabs(glob->dcaGlobal(vtx_id).mag()) < 3 && glob->nHitsFit(kTpcId) >= 10) grefmult++;            
        }
        return grefmult;
    }
        else return 0;
}

unsigned short StMuEvent::btofTrayMultiplicity(){

        unsigned short btofmult = (unsigned short)StMuDst::numberOfBTofHit();
        for(unsigned int i=0;i< StMuDst::numberOfBTofHit();i++) if(StMuDst::btofHit(i)->tray() > 120) btofmult--;
        return btofmult;

}

float StMuEvent::nearestVertexZ(int vtx_id){

        float dz = 999.0;
        //For old MuDsts where there was one vertex per event
        if (StMuDst::numberOfPrimaryVertices()==0) return dz;
        const int Nvert = StMuDst::primaryVertices()->GetEntries();
        if(Nvert < 2) return dz;
        
        if (vtx_id == -1) vtx_id = StMuDst::currentVertexIndex();       
        float z =  primaryVertexPosition(vtx_id).z();
        for(int i=0;i<Nvert;i++){
                if(vtx_id!=i) {
                        if(fabs(z-StMuDst::primaryVertex(i)->position().z()) < dz) dz = fabs(z-StMuDst::primaryVertex(i)->position().z());
                }
        }
        return dz;
}

/***************************************************************************
 *
 * $Log: StMuEvent.cxx,v $
 * Revision 1.26  2011/08/30 14:41:26  fisyak
 * Keep DAQ time for the event untouched
 *
 * Revision 1.25  2010/05/29 16:36:52  tone421
 * Added pointer protection to StTriggerData member
 *
 * Revision 1.24  2010/05/28 19:47:51  tone421
 * Removed a cout needed for test purposes in StMuDstMaker. Made sure StTriggerData objects copied into the MuDst have a debug value of 0..
 *
 * Revision 1.23  2010/05/26 04:25:50  tone421
 * Added StTriggerData arrays in muevent and fixed an issue with PMD arrays being read....
 *
 * Revision 1.22  2010/02/03 17:16:22  tone421
 * Added function StMuEvent::nearestVertexZ(int vtx_id) which returns the z distance of the nearest vertex in relation to vertex vtx_id
 *
 * Revision 1.21  2010/02/03 04:54:45  tone421
 * Added StMuEvent::btofTrayMultiplicity() to return only TOF hits from trays. Should be looked at instead of ctbSum for run 9 and beyond.
 *
 * Revision 1.20  2009/01/09 19:43:47  tone421
 * OAdded gremult in StMuEvent (globals tracks with DCA < 3cm, >= 10 TPC fit hits and |eta| < 0.5)
 *
 * Revision 1.19  2008/11/18 15:34:32  tone421
 * 2 changes. The first ensures StMuEvent::primaryVertexPosition() returns the position of current vertex (set by StMuDst::setVertexIndex(Int_t vtx_id)): previously it returned the position of best ranked vertex. The second insures events with no vertex have a PVx=PYy=PYz=-999 rather than 0.
 *
 * Revision 1.18  2008/06/26 15:46:13  tone421
 *
 * Add getter and setter for the vpd z vertex position
 *
 * Revision 1.17  2008/02/20 09:00:48  mvl
 * Included FMS data (StFMSTriggerDetector) (code by Akio)
 *
 * Revision 1.16  2007/09/21 02:27:12  mvl
 * Added calibrated VPD info from StTofCollection (run-8 prep)
 *
 * Revision 1.15  2007/09/05 23:21:21  mvl
 * Added StMtdTriggerDetector
 *
 * Revision 1.14  2007/04/20 06:25:21  mvl
 * Removed Q-vectors (will implement utility class).
 * Added Vpd info.
 *
 * Revision 1.12  2006/09/20 17:23:39  mvl
 * Added protected for events with StTriggerData (e.g. simulation)
 *
 * Revision 1.11  2006/09/20 01:50:35  mvl
 * Added data member and code for L2Result array (TArrayI).
 *
 * Revision 1.10  2005/08/19 19:46:05  mvl
 * Further updates for multiple vertices. The main changes are:
 * 1) StMudst::primaryTracks() now returns a list (TObjArray*) of tracks
 *    belonging to the 'current' primary vertex. The index number of the
 *    'current' vertex can be set using StMuDst::setCurrentVertex().
 *    This also affects StMuDst::primaryTracks(int i) and
 *    StMuDst::numberOfprimaryTracks().
 * 2) refMult is now stored for all vertices, in StMuPrimaryVertex. The
 *    obvious way to access these numbers is from the StMuprimaryVertex structures,
 *    but for ebakcward compatibility a function is provided in StMuEvent as well
 *    (this is the only function taht works for existing MuDst)
 *
 * As an aside, I've also changes the internals of StMuDst::createStEvent and
 * StMuDst::fixTrackIndices() to be able to deal with a larger range of index numbers for tracks as generated by Ittf.
 *
 * BIG FAT WARNING: StMudst2StEventMaker and StMuDstFilterMaker
 * do not fully support the multiple vertex functionality yet.
 *
 * Revision 1.9  2004/12/02 00:19:52  mvl
 * Added error on primary vertex
 *
 * Revision 1.8  2004/08/04 17:57:13  mvl
 * Added EMC trigger information and fpd trigger (tower) information
 *
 * Revision 1.7  2004/02/17 04:56:36  jeromel
 * Extended help, added crs support, restored __GNUC__ for PRETTY_FUNCTION(checked once
 * more and yes, it is ONLY defined in GCC and so is __FUCTION__),  use of a consistent
 * internal __PRETTYF__, return NULL if no case selected (+message) and protected against
 * NULL mChain.
 *
 * Revision 1.6  2003/10/20 19:50:13  perev
 * workaround added for TClonesArray::Delete + some cleanup of MuEmc
 *
 * Revision 1.5  2003/07/22 19:14:41  laue
 * multiplicities for FTPC added
 *
 * Revision 1.4  2003/02/20 15:29:42  laue
 * StMuTriggerIdCollection added
 *
 * Revision 1.3  2003/02/19 13:51:58  laue
 * added the StTriggerIdCollection
 *
 * Revision 1.2  2002/08/27 19:05:57  laue
 * Minor updates to make the muDst from simulation work
 *
 * Revision 1.1  2002/03/08 17:04:17  laue
 * initial revision
 *
 *
 **************************************************************************/

---