StFmsFilterMaker.cxx

//*-- Author : Jan Balewski
// 
// $Id: StFmsFilterMaker.cxx,v 1.1 2018/05/22 19:55:42 akio Exp $

#include <StThreeVector.hh>
#include <StPrimaryVertex.h>

#include "StFmsDbMaker/StFmsDbMaker.h"
#include "StEvent/StFmsCollection.h"
#include "StEvent/StFmsHit.h"

#include "StEvent.h"
#include "StMessMgr.h"

#include "StFmsFilterMaker.h"

ClassImp(StFmsFilterMaker)

//_____________________________________________________________
//_____________________________________________________________
StFmsFilterMaker::StFmsFilterMaker(const char *name):StMaker(name){
  par_E_thres= 30;
  par_Z0_vert=-9999;
  par_delZ_vert=0.;
  //  myMode=kUnknown;
}


//_____________________________________________________________
//_____________________________________________________________
StFmsFilterMaker::~StFmsFilterMaker(){
  
  //
}


//_____________________________________________________________
//_____________________________________________________________

Int_t StFmsFilterMaker::Init(){
  SetAttr(".Privilege",1); //this maker so BFC 'listens' to it and skip events generated by this maker

  //  mGeomE = new FmsGeomSimple(); //for trackXEndCap
  mFmsDbMaker=static_cast<StFmsDbMaker*>(GetMaker("fmsDb"));
  
  mH0=new TH1F("mH0","Event counter",5,0.5,5.5);
  nInpEve=nRecVert=nZverOK=nAccEve=0;
  LOG_INFO << Form("Init ee-filter cuts: ET>%.2f  Zvert=%.2f +/-%.2f (cm)", par_E_thres, par_Z0_vert, par_delZ_vert)<<endm;
  //  if(myMode==kUseFixedVertex)  LOG_INFO << Form("  FIXED vertex  mode is activated")<<endl;
  //  assert  (myMode!=kUnknown); // probably you forgot to set all params for this maker, Jan
  return StMaker::Init();
}

//_____________________________________________________________
//_____________________________________________________________
Int_t StFmsFilterMaker::FinishRun(int runumber){
  LOG_INFO << Form("Finish cuts: ET>%.2f  Zvert=%.2f +/-%.2f (cm) fixVertexMode=%d", par_E_thres, par_Z0_vert, par_delZ_vert,myMode)<<endm;
  LOG_INFO << Form("%s::Finish run=%d nInp=%d,nRecVer=%d, nZverOK=%d nAcc=%d",GetName(),runumber,nInpEve,nRecVert,nZverOK,nAccEve) << endm;
  return kStOK;
}; 


//_____________________________________________________________
//_____________________________________________________________
Int_t StFmsFilterMaker::Make(){
  LOG_INFO << Form("in::Make inp=%d,nVer=%d, nZverOK=%d nAcc=%d",nInpEve,nRecVert,nZverOK,nAccEve) << endm;
  nInpEve++;
  mH0->Fill(1);
  StEvent* mEvent = (StEvent*)GetInputDS("StEvent");
  assert(mEvent);//fix your chain or open the right event file 
  Float_t vertexPosZ =par_Z0_vert ;
  /*
  if(myMode==kUseRecoVertex) {
    int nV=mEvent->numberOfPrimaryVertices();
    if(nV<=0) return  kStSKIP;
    int iv=0; // pick always the first vertex, may be wrong w/ pileup, fix it,JB
    nRecVert++;
    mH0->Fill(2); 
    StPrimaryVertex*V=mEvent->primaryVertex(iv);
    assert(V);
    StThreeVectorF r=V->position();
    vertexPosZ = r.z();
    if(fabs(vertexPosZ - par_Z0_vert)> par_delZ_vert)  return  kStSKIP;
  LOG_INFO << Form("eveID=%d  nPrimVert=%d  zVert=%.2f\n", mEvent->id(),nV, vertexPosZ);
  }
  */
  nZverOK++;
  mH0->Fill(3); 

  StEvent* event = (StEvent*)GetInputDS("StEvent");
  if(!event) {LOG_ERROR << "StFmsJetMaker::Make did not find StEvent"<<endm; return kStErr;}
  mFmsColl = event->fmsCollection();

  
  Float_t triggerPatchEt;
  Int_t triggerConditionReturn = triggerCondition();


  if(triggerConditionReturn != 1)    {
    //  printf("janKill-%d Me### eve=%d  zVer=%.2f highT_et=%.3f TP_ET=%.3f\n",myMode,nInpEve,vertexPosZ, highTEt,triggerPatchEt);
    return kStSKIP;
  }

  nAccEve++;
  mH0->Fill(4); 
  
 //process event only if et above threshold
  //printf("janKeep-%d Me### eve=%d  zVer=%.2f highT_et=%.3f TP_ET=%.3f\n",myMode,nInpEve,vertexPosZ, highTEt,triggerPatchEt);

  return kStOK;
}
//_____________________________________________________________
//_____________________________________________________________

Int_t StFmsFilterMaker::triggerCondition()
{  
  int nhits=mFmsColl->numberOfHits();

  StSPtrVecFmsHit& hits = mFmsColl->hits();

  Float_t M_pi=TMath::Pi();
  Float_t energyQuadrant[8] = {0,0,0,0,0,0,0,0};
  Int_t binPhi = 8;
  for(int i=0; i<nhits; i++){
    Int_t dId = hits[i]->detectorId();
    Float_t energy = hits[i]->energy();
    if(dId>7&&dId<12 && energy>0.2){
      int ch= hits[i]->channel();
      
      StThreeVectorF v3  =mFmsDbMaker->getStarXYZ(dId,ch);
      Float_t x = v3.x();
      Float_t y = v3.y();
      Int_t id = dId*1000+ch;
      Float_t phi =  v3.phi();
      
      binPhi = 8;
      if((phi< 1.0*(M_pi/4))    && (phi > -1.0*(M_pi/4))) binPhi = 0;
      if((phi<   0*(M_pi/4))    && (phi > -2.0*(M_pi/4))) binPhi = 1;
      if((phi< -1.0*(M_pi/4))   && (phi > -3.0*(M_pi/4))) binPhi = 2;
      if((phi< -2.0*(M_pi/4))   && (phi > -4.0*(M_pi/4))) binPhi = 3;

      if((phi < -3.0*(M_pi/4))  && (phi > 3.0*(M_pi/4))) binPhi = 4;
      if((phi <  4.0*(M_pi/4))  && (phi > 2.0*(M_pi/4))) binPhi = 5;
      if((phi <  3.0*(M_pi/4))  && (phi >  1.*(M_pi/4))) binPhi = 6;
      if((phi <  2.0*(M_pi/4))  && (phi >   0*(M_pi/4))) binPhi = 7;

      if(binPhi<8)   energyQuadrant[binPhi] += energy;


    }
  }
  Int_t iOrd[8];
  TMath::Sort(8,  energyQuadrant, iOrd);
  if( energyQuadrant[iOrd[0]] < par_E_thres ) return 0;

  return 1;
  
}

/*
//_____________________________________________________________
//_____________________________________________________________
Float_t StFmsFilterMaker::transverseNRG(Float_t vertexPosZ, StFmsTower *tower)
{
  //find the transverse energy for a tower when vertex != 0
  //now we do everything with respect to the Smd plane




  TVector3 towVector = mGeomE->getTowerCenter(tower->sector(),tower->subsector(),tower->etabin());
  //get x,y coords at Smd plane
  Float_t xNew = towVector.x() - (288.2 - 279.54)*(towVector.x()/288.2);
  Float_t yNew = towVector.y() - (288.2 - 279.54)*(towVector.y()/288.2);

  TVector3 *eventTowVector = new TVector3(0,0,0);//make vector from vertex to smd plane
  eventTowVector->SetXYZ(xNew,yNew,(279.54 - vertexPosZ));

  //Float_t newEta = eventTowVector->Eta();
  //Float_t newPhi = eventTowVector->Phi();
  Float_t newTheta = eventTowVector->Theta();

  Float_t factor = sin(newTheta);//multiply energy by this to get event ET
  if(factor > 1.0 || factor < 0.0)
    {
      return -1.0;
    }

  Float_t energy = tower->energy();

  Float_t transverseEnergy = energy*factor;

  return transverseEnergy;
  

}
*/
// $Log: StFmsFilterMaker.cxx,v $
// Revision 1.1  2018/05/22 19:55:42  akio
// Initial version of after GEANT filter for FMS by Mriganka
//
// Revision 1.4  2008/06/10 12:59:12  balewski
// fix printout
//
// Revision 1.3  2008/05/11 18:49:17  balewski
// merged 2 makers, now one can us it before and after TPC tracking, run 2 independent copies
//
// Revision 1.2  2008/05/09 22:14:35  balewski
// new there are 2 filters
//
// Revision 1.1  2008/04/21 15:47:09  balewski
// star
//