[u/mrichter/AliRoot.git] / HLT / global / AliHLTGlobalVertexerComponent.cxx

//**************************************************************************
//* This file is property of and copyright by the ALICE HLT Project        * 
//* ALICE Experiment at CERN, All rights reserved.                         *
//*                                                                        *
//* Primary Authors: S.Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de>    *
//*                  for The ALICE HLT Project.                            *
//*                                                                        *
//* Permission to use, copy, modify and distribute this software and its   *
//* documentation strictly for non-commercial purposes is hereby granted   *
//* without fee, provided that the above copyright notice appears in all   *
//* copies and that both the copyright notice and this permission notice   *
//* appear in the supporting documentation. The authors make no claims     *
//* about the suitability of this software for any purpose. It is          *
//* provided "as is" without express or implied warranty.                  *
//**************************************************************************

/** @file   AliHLTGlobalVertexerComponent.cxx
    @author Sergey Gorbunov
    @brief  Component for reconstruct primary vertex and V0's
*/

#if __GNUC__>= 3
using namespace std;
#endif

#include "AliHLTGlobalVertexerComponent.h"
#include "AliCDBEntry.h"
#include "AliCDBManager.h"
#include <TFile.h>
#include <TString.h>
#include "TObjString.h"
#include "TObjArray.h"
#include "TH1F.h"
#include "TH2F.h"
#include "AliESDEvent.h"
#include "AliESDtrack.h"
#include "AliESDVertex.h"
#include "AliESDv0.h"
#include "AliHLTMessage.h"
#include "TMath.h"
#include "AliKFParticle.h"
#include "AliKFVertex.h"


/** ROOT macro for the implementation of ROOT specific class methods */
ClassImp(AliHLTGlobalVertexerComponent)

AliHLTGlobalVertexerComponent::AliHLTGlobalVertexerComponent()
:
  fESD(0),
  fTrackInfos(0),
  fPrimaryVtx(),  
  fHistPrimVertexXY(0),
  fHistPrimVertexZX(0),
  fHistPrimVertexZY(0),
  fNEvents(0),
  fPlotHistograms(1),
  fFitTracksToVertex(1),
  fConstrainedTrackDeviation(4.),
  fV0DaughterPrimDeviation( 2.5 ),
  fV0PrimDeviation( 3.5 ),
  fV0Chi(3.5),
  fV0DecayLengthInSigmas(3.)
{
  // see header file for class documentation
  // or
  // refer to README to build package
  // or
  // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt

}

AliHLTGlobalVertexerComponent::~AliHLTGlobalVertexerComponent()
{
  // see header file for class documentation

  delete[] fTrackInfos;
  delete   fHistPrimVertexXY;
  delete   fHistPrimVertexZX;
  delete   fHistPrimVertexZY;
}

// Public functions to implement AliHLTComponent's interface.
// These functions are required for the registration process

const char* AliHLTGlobalVertexerComponent::GetComponentID()
{
  // see header file for class documentation
  
  return "GlobalVertexer";
}

void AliHLTGlobalVertexerComponent::GetInputDataTypes(AliHLTComponentDataTypeList& list)
{
  // see header file for class documentation
  list.clear();
  list.push_back( kAliHLTDataTypeESDObject|kAliHLTDataOriginOut );
}

AliHLTComponentDataType AliHLTGlobalVertexerComponent::GetOutputDataType()
{
  // see header file for class documentation
  return kAliHLTMultipleDataType;
}

int AliHLTGlobalVertexerComponent::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList)

{
  // see header file for class documentation
  tgtList.clear();
  tgtList.push_back(kAliHLTDataTypeESDObject|kAliHLTDataOriginOut);
  tgtList.push_back(kAliHLTDataTypeHistogram);
  return tgtList.size();
}

void AliHLTGlobalVertexerComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
{
  // see header file for class documentation
  // XXX TODO: Find more realistic values.
  constBase = 80000;
  inputMultiplier = 2.;
}

AliHLTComponent* AliHLTGlobalVertexerComponent::Spawn()
{
  // see header file for class documentation
  return new AliHLTGlobalVertexerComponent;
}

int AliHLTGlobalVertexerComponent::DoInit( int argc, const char** argv )
{
  // init
  
  fHistPrimVertexXY = new TH2F("primVertexXY","HLT:  Primary vertex distribution in XY",60,-2.,2.,60,-2.,2.);
  fHistPrimVertexXY->SetMarkerStyle(8);
  fHistPrimVertexXY->SetMarkerSize(0.4);
  fHistPrimVertexXY->SetYTitle("Y");
  fHistPrimVertexXY->SetXTitle("X");
  fHistPrimVertexXY->SetStats(0);
  //fHistPrimVertexXY->SetOption("COLZ");

  fHistPrimVertexZX = new TH2F("primVertexZX","HLT:  Primary vertex distribution in ZX",60,-15.,15.,60,-2.,2.);
  fHistPrimVertexZX->SetMarkerStyle(8);
  fHistPrimVertexZX->SetMarkerSize(0.4);
  fHistPrimVertexZX->SetYTitle("X");
  fHistPrimVertexZX->SetXTitle("Z");
  fHistPrimVertexZX->SetStats(0);
  //fHistPrimVertexZX->SetOption("COLZ");

  fHistPrimVertexZY = new TH2F("primVertexZY","HLT:  Primary vertex distribution in ZY",60,-10.,10.,60,-2.,2.);
  fHistPrimVertexZY->SetMarkerStyle(8);
  fHistPrimVertexZY->SetMarkerSize(0.4);
  fHistPrimVertexZY->SetYTitle("Y");
  fHistPrimVertexZY->SetXTitle("Z");
  fHistPrimVertexZY->SetStats(0);
  //fHistPrimVertexZY->SetOption("COLZ");

  fNEvents =0;

  int iResult=0;
  TString configuration="";
  TString argument="";
  for (int i=0; i<argc && iResult>=0; i++) {
    argument=argv[i];
    if (!configuration.IsNull()) configuration+=" ";
    configuration+=argument;
  }
  
  if (!configuration.IsNull()) {
    iResult=Configure(configuration.Data());
  }  

  return iResult; 
}
  
int AliHLTGlobalVertexerComponent::DoDeinit()
{
  // see header file for class documentation

  delete[] fTrackInfos;
  delete fHistPrimVertexXY;
  delete fHistPrimVertexZX;
  delete fHistPrimVertexZY;

  fTrackInfos = 0;
  fHistPrimVertexXY = 0;
  fHistPrimVertexZX = 0;
  fHistPrimVertexZY = 0;

  fPlotHistograms = 1;
  fFitTracksToVertex = 1;
  fConstrainedTrackDeviation = 4.;
  fV0DaughterPrimDeviation = 2.5 ;
  fV0PrimDeviation =3.5;
  fV0Chi = 3.5;
  fV0DecayLengthInSigmas = 3.;
  fNEvents = 0;
  return 0;
}

int AliHLTGlobalVertexerComponent::DoEvent(const AliHLTComponentEventData& /*evtData*/, AliHLTComponentTriggerData& /*trigData*/)
{
  
  if ( GetFirstInputBlock( kAliHLTDataTypeSOR ) || GetFirstInputBlock( kAliHLTDataTypeEOR ) )
    return 0;

  fNEvents++;

  for ( const TObject *iter = GetFirstInputObject(kAliHLTDataTypeESDObject); iter != NULL; iter = GetNextInputObject() ) {

    AliESDEvent *event = dynamic_cast<AliESDEvent*>(const_cast<TObject*>( iter ) );
    event->GetStdContent();

    // primary vertex & V0's 
          
    SetESD( event );
    FindPrimaryVertex();
    FindV0s();
    const AliESDVertex *vPrim = event->GetPrimaryVertexTracks();

    // plot histograms

    if( fPlotHistograms ){
      if( vPrim && vPrim->GetNContributors()>=3 ){
	fHistPrimVertexXY->Fill( vPrim->GetX(),vPrim->GetY()  );
	fHistPrimVertexZX->Fill( vPrim->GetZ(),vPrim->GetX()  );
	fHistPrimVertexZY->Fill( vPrim->GetZ(),vPrim->GetY()  );    
      }  
      if( fHistPrimVertexXY ) PushBack( (TObject*) fHistPrimVertexXY, kAliHLTDataTypeHistogram,0);
      if( fHistPrimVertexZX ) PushBack( (TObject*) fHistPrimVertexZX, kAliHLTDataTypeHistogram,0);
      if( fHistPrimVertexZY ) PushBack( (TObject*) fHistPrimVertexZY, kAliHLTDataTypeHistogram,0);
    }
    
    PushBack( (TObject*) event, kAliHLTDataTypeESDObject|kAliHLTDataOriginOut, 0);
  }
    
  return 0;
}

int AliHLTGlobalVertexerComponent::Configure(const char* arguments)
{
  
  int iResult=0;
  if (!arguments) return iResult;
  
  TString allArgs=arguments;
  TString argument;
  
  TObjArray* pTokens=allArgs.Tokenize(" ");
  int bMissingParam=0;

  if (pTokens) {
    for (int i=0; i<pTokens->GetEntries() && iResult>=0; i++) {
      argument=((TObjString*)pTokens->At(i))->GetString();
      if (argument.IsNull()) continue;
      
      if (argument.CompareTo("-fitTracksToVertex")==0) {
	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	HLTInfo("fitTracksToVertex is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	fFitTracksToVertex=((TObjString*)pTokens->At(i))->GetString().Atoi();
	continue;
      }
      else if (argument.CompareTo("-plotHistograms")==0) {
	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	HLTInfo("plotHistograms is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	fPlotHistograms=((TObjString*)pTokens->At(i))->GetString().Atoi();
	continue;
      }
      else if (argument.CompareTo("-constrainedTrackDeviation")==0) {
	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	HLTInfo("constrainedTrackDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	fConstrainedTrackDeviation=((TObjString*)pTokens->At(i))->GetString().Atof();
	continue;
      }
      else if (argument.CompareTo("-v0DaughterPrimDeviation")==0) {
	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	HLTInfo("v0DaughterPrimDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	fV0DaughterPrimDeviation=((TObjString*)pTokens->At(i))->GetString().Atof();
	continue;
      }
      else if (argument.CompareTo("-v0PrimDeviation")==0) {
	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	HLTInfo("v0PrimDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	fV0PrimDeviation=((TObjString*)pTokens->At(i))->GetString().Atof();
	continue;
      }
      else if (argument.CompareTo("-v0Chi")==0) {
	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	HLTInfo("v0Chi is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	fV0Chi=((TObjString*)pTokens->At(i))->GetString().Atof();
	continue;
      }
      else if (argument.CompareTo("-v0DecayLengthInSigmas")==0) {
	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	HLTInfo("v0DecayLengthInSigmas is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	fV0DecayLengthInSigmas=((TObjString*)pTokens->At(i))->GetString().Atof();
	continue;
      }
      else {
	HLTError("unknown argument %s", argument.Data());
	iResult=-EINVAL;
	break;
      }
    }
    delete pTokens;
  }
  if (bMissingParam) {
    HLTError("missing parameter for argument %s", argument.Data());
    iResult=-EINVAL;
  }  
  
  return iResult;
}

int AliHLTGlobalVertexerComponent::Reconfigure(const char* cdbEntry, const char* chainId)
{
  // see header file for class documentation
  int iResult=0;
  const char* path="HLT/ConfigTPC/KryptonHistoComponent";
  const char* defaultNotify="";
  if (cdbEntry) {
    path=cdbEntry;
    defaultNotify=" (default)";
  }
  if (path) {
    HLTInfo("reconfigure from entry %s%s, chain id %s", path, defaultNotify,(chainId!=NULL && chainId[0]!=0)?chainId:"<none>");
    AliCDBEntry *pEntry = AliCDBManager::Instance()->Get(path/*,GetRunNo()*/);
    if (pEntry) {
      TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject());
      if (pString) {
	HLTInfo("received configuration object string: \'%s\'", pString->GetString().Data());
	iResult=Configure(pString->GetString().Data());
      } else {
	HLTError("configuration object \"%s\" has wrong type, required TObjString", path);
      }
    } else {
      HLTError("can not fetch object \"%s\" from CDB", path);
    }
  }

  return iResult;
}


void AliHLTGlobalVertexerComponent::SetESD( AliESDEvent *event )
{
  //* Fill fTrackInfo array

  delete[] fTrackInfos;
  fESD = event;

  AliKFParticle::SetField( fESD->GetMagneticField() );

  Int_t nESDTracks=event->GetNumberOfTracks(); 
  fTrackInfos = new AliESDTrackInfo[ nESDTracks ];

  for (Int_t iTr=0; iTr<nESDTracks; iTr++){ 
  
    AliESDTrackInfo &info = fTrackInfos[iTr];
    info.fOK = 0;
    info.fPrimUsedFlag = 0;
    
    //* track quality check

    AliESDtrack *pTrack = event->GetTrack(iTr);    
    if( !pTrack  ) continue;
    if (pTrack->GetKinkIndex(0)>0) continue;
    if ( !( pTrack->GetStatus()&AliESDtrack::kTPCin ) ) continue;
    
    //* Construct KFParticle for the track

    info.fParticle = AliKFParticle( *pTrack->GetInnerParam(), 211 );    
    info.fOK = 1;
  }
}


void AliHLTGlobalVertexerComponent::FindPrimaryVertex(  )
{
  //* Find event primary vertex

  int nTracks = fESD->GetNumberOfTracks();

  const AliKFParticle **vSelected = new const AliKFParticle*[nTracks]; //* Selected particles for vertex fit
  Int_t *vIndex = new int [nTracks];                    //* Indices of selected particles
  Bool_t *vFlag = new bool [nTracks];                    //* Flags returned by the vertex finder

  fPrimaryVtx.Initialize();
  fPrimaryVtx.SetBeamConstraint(fESD->GetDiamondX(),fESD->GetDiamondY(),0,
				TMath::Sqrt(fESD->GetSigma2DiamondX()),TMath::Sqrt(fESD->GetSigma2DiamondY()),5.3);
  
  Int_t nSelected = 0;
  for( Int_t i = 0; i<nTracks; i++){ 
    if(!fTrackInfos[i].fOK ) continue;
    //if( fESD->GetTrack(i)->GetTPCNcls()<60  ) continue;
    const AliKFParticle &p = fTrackInfos[i].fParticle;
    Double_t chi = p.GetDeviationFromVertex( fPrimaryVtx );      
    if( chi > fConstrainedTrackDeviation ) continue;
    vSelected[nSelected] = &(fTrackInfos[i].fParticle);
    vIndex[nSelected] = i;
    nSelected++;  
  }
  fPrimaryVtx.ConstructPrimaryVertex( vSelected, nSelected, vFlag, fConstrainedTrackDeviation );
  for( Int_t i = 0; i<nSelected; i++){ 
    if( vFlag[i] ) fTrackInfos[vIndex[i]].fPrimUsedFlag = 1;
  }

  for( Int_t i = 0; i<nTracks; i++ ){
    AliESDTrackInfo &info = fTrackInfos[i];
    info.fPrimDeviation = info.fParticle.GetDeviationFromVertex( fPrimaryVtx );   
  }
  
  if( fPrimaryVtx.GetNContributors()>=3 ){
    AliESDVertex vESD( fPrimaryVtx.Parameters(), fPrimaryVtx.CovarianceMatrix(), fPrimaryVtx.GetChi2(), fPrimaryVtx.GetNContributors() );
    fESD->SetPrimaryVertexTracks( &vESD );

    // relate the tracks to vertex

    if( fFitTracksToVertex ){      
      for( Int_t i = 0; i<nTracks; i++ ){
	if( !fTrackInfos[i].fPrimUsedFlag ) continue;	  
	if( fTrackInfos[i].fPrimDeviation > fConstrainedTrackDeviation ) continue;
	fESD->GetTrack(i)->RelateToVertex( &vESD, fESD->GetMagneticField(),100. );
      }
    }

  } else {
    for( Int_t i = 0; i<nTracks; i++)
      fTrackInfos[i].fPrimUsedFlag = 0;
  }


  delete[] vSelected;
  delete[] vIndex;
  delete[] vFlag;
}


void AliHLTGlobalVertexerComponent::FindV0s(  )
{
  //* V0 finder

  int nTracks = fESD->GetNumberOfTracks();
  //AliKFVertex primVtx( *fESD->GetPrimaryVertexTracks() );
  AliKFVertex &primVtx = fPrimaryVtx;
  if( primVtx.GetNContributors()<3 ) return;

  bool *constrainedV0   = new bool[nTracks];
  for( Int_t iTr = 0; iTr<nTracks; iTr++ ){ 
    constrainedV0[iTr] = 0;
  }
  

  for( Int_t iTr = 0; iTr<nTracks; iTr++ ){ //* first daughter

    AliESDTrackInfo &info = fTrackInfos[iTr];
    if( !info.fOK ) continue;    
    if( info.fParticle.GetQ() >0 ) continue;    
    if( info.fPrimDeviation < fV0DaughterPrimDeviation ) continue;

    for( Int_t jTr = 0; jTr<nTracks; jTr++ ){  //* second daughter

      AliESDTrackInfo &jnfo = fTrackInfos[jTr];
      if( !jnfo.fOK ) continue;
      if( jnfo.fParticle.GetQ() < 0 ) continue;
      if( jnfo.fPrimDeviation < fV0DaughterPrimDeviation ) continue;
   
      //* construct V0 mother

      AliKFParticle v0( info.fParticle, jnfo.fParticle );     

      //* check V0 Chi^2
      
      if( v0.GetChi2()<0 || v0.GetChi2() > fV0Chi*fV0Chi*v0.GetNDF() ) continue;

      //* subtruct daughters from primary vertex 

      AliKFVertex primVtxCopy = primVtx;    
       
      if( info.fPrimUsedFlag ){	
	if( primVtxCopy.GetNContributors()<=2 ) continue;
	primVtxCopy -= info.fParticle;
      }
      if( jnfo.fPrimUsedFlag ){
	if( primVtxCopy.GetNContributors()<=2 ) continue;
	primVtxCopy -= jnfo.fParticle;
      }
      //* Check v0 Chi^2 deviation from primary vertex 

      if( v0.GetDeviationFromVertex( primVtxCopy ) > fV0PrimDeviation ) continue;

      //* Add V0 to primary vertex to improve the primary vertex resolution

      primVtxCopy += v0;      

      //* Set production vertex for V0

      v0.SetProductionVertex( primVtxCopy );

      //* Get V0 decay length with estimated error

      Double_t length, sigmaLength;
      if( v0.GetDecayLength( length, sigmaLength ) ) continue;

      //* Reject V0 if it decays too close[sigma] to the primary vertex

      if( length  < fV0DecayLengthInSigmas*sigmaLength ) continue;

      //* add ESD v0 
      
      AliESDv0 v0ESD( *fESD->GetTrack( iTr ), iTr, *fESD->GetTrack( jTr ), jTr );  
      fESD->AddV0( &v0ESD );

      // relate the tracks to vertex

      if( fFitTracksToVertex ){
	if( constrainedV0[iTr] || constrainedV0[jTr]
	    || info.fPrimDeviation < fConstrainedTrackDeviation || jnfo.fPrimDeviation < fConstrainedTrackDeviation ) continue;
	AliESDVertex vESD(v0.Parameters(), v0.CovarianceMatrix(), v0.GetChi2(), 2);
	fESD->GetTrack(iTr)->RelateToVertex( &vESD, fESD->GetMagneticField(),100. );
	fESD->GetTrack(jTr)->RelateToVertex( &vESD, fESD->GetMagneticField(),100. );
	constrainedV0[iTr] = 1;
	constrainedV0[jTr] = 1;	
      }
    }
  }
  delete[] constrainedV0;
}
Commit	Line	Data
4d5ee3db	1	//**************************************************************************
	2	//* This file is property of and copyright by the ALICE HLT Project *
	3	//* ALICE Experiment at CERN, All rights reserved. *
	4	//* *
	5	//* Primary Authors: S.Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de> *
	6	//* for The ALICE HLT Project. *
	7	//* *
	8	//* Permission to use, copy, modify and distribute this software and its *
	9	//* documentation strictly for non-commercial purposes is hereby granted *
	10	//* without fee, provided that the above copyright notice appears in all *
	11	//* copies and that both the copyright notice and this permission notice *
	12	//* appear in the supporting documentation. The authors make no claims *
	13	//* about the suitability of this software for any purpose. It is *
	14	//* provided "as is" without express or implied warranty. *
	15	//**************************************************************************
	16
	17	/** @file AliHLTGlobalVertexerComponent.cxx
	18	@author Sergey Gorbunov
	19	@brief Component for reconstruct primary vertex and V0's
	20	*/
	21
	22	#if __GNUC__>= 3
	23	using namespace std;
	24	#endif
	25
	26	#include "AliHLTGlobalVertexerComponent.h"
	27	#include "AliCDBEntry.h"
	28	#include "AliCDBManager.h"
	29	#include <TFile.h>
	30	#include <TString.h>
	31	#include "TObjString.h"
	32	#include "TObjArray.h"
	33	#include "TH1F.h"
	34	#include "TH2F.h"
	35	#include "AliESDEvent.h"
	36	#include "AliESDtrack.h"
	37	#include "AliESDVertex.h"
	38	#include "AliESDv0.h"
	39	#include "AliHLTMessage.h"
4d5ee3db	40	#include "TMath.h"
	41	#include "AliKFParticle.h"
	42	#include "AliKFVertex.h"
	43
	44
	45	/** ROOT macro for the implementation of ROOT specific class methods */
	46	ClassImp(AliHLTGlobalVertexerComponent)
	47
	48	AliHLTGlobalVertexerComponent::AliHLTGlobalVertexerComponent()
	49	:
	50	fESD(0),
	51	fTrackInfos(0),
	52	fPrimaryVtx(),
	53	fHistPrimVertexXY(0),
	54	fHistPrimVertexZX(0),
	55	fHistPrimVertexZY(0),
	56	fNEvents(0),
	57	fPlotHistograms(1),
608cfbda	58	fFitTracksToVertex(1),
4d5ee3db	59	fConstrainedTrackDeviation(4.),
	60	fV0DaughterPrimDeviation( 2.5 ),
	61	fV0PrimDeviation( 3.5 ),
	62	fV0Chi(3.5),
	63	fV0DecayLengthInSigmas(3.)
	64	{
	65	// see header file for class documentation
	66	// or
	67	// refer to README to build package
	68	// or
	69	// visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
	70
	71	}
	72
	73	AliHLTGlobalVertexerComponent::~AliHLTGlobalVertexerComponent()
	74	{
	75	// see header file for class documentation
	76
	77	delete[] fTrackInfos;
	78	delete fHistPrimVertexXY;
	79	delete fHistPrimVertexZX;
	80	delete fHistPrimVertexZY;
	81	}
	82
	83	// Public functions to implement AliHLTComponent's interface.
	84	// These functions are required for the registration process
	85
	86	const char* AliHLTGlobalVertexerComponent::GetComponentID()
	87	{
	88	// see header file for class documentation
	89
	90	return "GlobalVertexer";
	91	}
	92
	93	void AliHLTGlobalVertexerComponent::GetInputDataTypes(AliHLTComponentDataTypeList& list)
	94	{
	95	// see header file for class documentation
	96	list.clear();
	97	list.push_back( kAliHLTDataTypeESDObject\|kAliHLTDataOriginOut );
	98	}
	99
	100	AliHLTComponentDataType AliHLTGlobalVertexerComponent::GetOutputDataType()
	101	{
	102	// see header file for class documentation
	103	return kAliHLTMultipleDataType;
	104	}
	105
	106	int AliHLTGlobalVertexerComponent::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList)
	107
	108	{
	109	// see header file for class documentation
	110	tgtList.clear();
	111	tgtList.push_back(kAliHLTDataTypeESDObject\|kAliHLTDataOriginOut);
	112	tgtList.push_back(kAliHLTDataTypeHistogram);
	113	return tgtList.size();
	114	}
	115
	116	void AliHLTGlobalVertexerComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
	117	{
	118	// see header file for class documentation
	119	// XXX TODO: Find more realistic values.
	120	constBase = 80000;
	121	inputMultiplier = 2.;
	122	}
123
124	AliHLTComponent* AliHLTGlobalVertexerComponent::Spawn()
125	{
126	// see header file for class documentation
127	return new AliHLTGlobalVertexerComponent;
128	}
129
130	int AliHLTGlobalVertexerComponent::DoInit( int argc, const char** argv )
131	{
132	// init
133
608cfbda	134	fHistPrimVertexXY = new TH2F("primVertexXY","HLT: Primary vertex distribution in XY",60,-2.,2.,60,-2.,2.);
4d5ee3db	135	fHistPrimVertexXY->SetMarkerStyle(8);
	136	fHistPrimVertexXY->SetMarkerSize(0.4);
	137	fHistPrimVertexXY->SetYTitle("Y");
	138	fHistPrimVertexXY->SetXTitle("X");
	139	fHistPrimVertexXY->SetStats(0);
608cfbda	140	//fHistPrimVertexXY->SetOption("COLZ");
4d5ee3db	141
608cfbda	142	fHistPrimVertexZX = new TH2F("primVertexZX","HLT: Primary vertex distribution in ZX",60,-15.,15.,60,-2.,2.);
4d5ee3db	143	fHistPrimVertexZX->SetMarkerStyle(8);
	144	fHistPrimVertexZX->SetMarkerSize(0.4);
	145	fHistPrimVertexZX->SetYTitle("X");
	146	fHistPrimVertexZX->SetXTitle("Z");
	147	fHistPrimVertexZX->SetStats(0);
608cfbda	148	//fHistPrimVertexZX->SetOption("COLZ");
4d5ee3db	149
608cfbda	150	fHistPrimVertexZY = new TH2F("primVertexZY","HLT: Primary vertex distribution in ZY",60,-10.,10.,60,-2.,2.);
4d5ee3db	151	fHistPrimVertexZY->SetMarkerStyle(8);
	152	fHistPrimVertexZY->SetMarkerSize(0.4);
	153	fHistPrimVertexZY->SetYTitle("Y");
	154	fHistPrimVertexZY->SetXTitle("Z");
	155	fHistPrimVertexZY->SetStats(0);
608cfbda	156	//fHistPrimVertexZY->SetOption("COLZ");
4d5ee3db	157
	158	fNEvents =0;
	159
	160	int iResult=0;
	161	TString configuration="";
	162	TString argument="";
	163	for (int i=0; i<argc && iResult>=0; i++) {
	164	argument=argv[i];
	165	if (!configuration.IsNull()) configuration+=" ";
	166	configuration+=argument;
	167	}
	168
	169	if (!configuration.IsNull()) {
	170	iResult=Configure(configuration.Data());
	171	}
	172
	173	return iResult;
	174	}
	175
	176	int AliHLTGlobalVertexerComponent::DoDeinit()
	177	{
	178	// see header file for class documentation
	179
	180	delete[] fTrackInfos;
	181	delete fHistPrimVertexXY;
	182	delete fHistPrimVertexZX;
	183	delete fHistPrimVertexZY;
	184
	185	fTrackInfos = 0;
	186	fHistPrimVertexXY = 0;
	187	fHistPrimVertexZX = 0;
	188	fHistPrimVertexZY = 0;
	189
	190	fPlotHistograms = 1;
608cfbda	191	fFitTracksToVertex = 1;
4d5ee3db	192	fConstrainedTrackDeviation = 4.;
	193	fV0DaughterPrimDeviation = 2.5 ;
	194	fV0PrimDeviation =3.5;
	195	fV0Chi = 3.5;
	196	fV0DecayLengthInSigmas = 3.;
	197	fNEvents = 0;
	198	return 0;
	199	}
	200
	201	int AliHLTGlobalVertexerComponent::DoEvent(const AliHLTComponentEventData& /evtData/, AliHLTComponentTriggerData& /trigData/)
	202	{
	203
	204	if ( GetFirstInputBlock( kAliHLTDataTypeSOR ) \|\| GetFirstInputBlock( kAliHLTDataTypeEOR ) )
	205	return 0;
	206
	207	fNEvents++;
	208
	209	for ( const TObject *iter = GetFirstInputObject(kAliHLTDataTypeESDObject); iter != NULL; iter = GetNextInputObject() ) {
	210
	211	AliESDEvent event = dynamic_cast<AliESDEvent>(const_cast<TObject*>( iter ) );
	212	event->GetStdContent();
	213
	214	// primary vertex & V0's
608cfbda	215
	216	SetESD( event );
	217	FindPrimaryVertex();
	218	FindV0s();
4d5ee3db	219	const AliESDVertex *vPrim = event->GetPrimaryVertexTracks();
	220
	221	// plot histograms
	222
	223	if( fPlotHistograms ){
	224	if( vPrim && vPrim->GetNContributors()>=3 ){
	225	fHistPrimVertexXY->Fill( vPrim->GetX(),vPrim->GetY() );
	226	fHistPrimVertexZX->Fill( vPrim->GetZ(),vPrim->GetX() );
	227	fHistPrimVertexZY->Fill( vPrim->GetZ(),vPrim->GetY() );
	228	}
	229	if( fHistPrimVertexXY ) PushBack( (TObject*) fHistPrimVertexXY, kAliHLTDataTypeHistogram,0);
	230	if( fHistPrimVertexZX ) PushBack( (TObject*) fHistPrimVertexZX, kAliHLTDataTypeHistogram,0);
	231	if( fHistPrimVertexZY ) PushBack( (TObject*) fHistPrimVertexZY, kAliHLTDataTypeHistogram,0);
	232	}
	233
f716b815	234	PushBack( (TObject*) event, kAliHLTDataTypeESDObject\|kAliHLTDataOriginOut, 0);
4d5ee3db	235	}
	236
	237	return 0;
	238	}
	239
	240	int AliHLTGlobalVertexerComponent::Configure(const char* arguments)
	241	{
	242
	243	int iResult=0;
	244	if (!arguments) return iResult;
	245
	246	TString allArgs=arguments;
	247	TString argument;
	248
	249	TObjArray* pTokens=allArgs.Tokenize(" ");
	250	int bMissingParam=0;
	251
	252	if (pTokens) {
	253	for (int i=0; i<pTokens->GetEntries() && iResult>=0; i++) {
	254	argument=((TObjString*)pTokens->At(i))->GetString();
	255	if (argument.IsNull()) continue;
	256
	257	if (argument.CompareTo("-fitTracksToVertex")==0) {
	258	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	259	HLTInfo("fitTracksToVertex is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
608cfbda	260	fFitTracksToVertex=((TObjString*)pTokens->At(i))->GetString().Atoi();
4d5ee3db	261	continue;
	262	}
	263	else if (argument.CompareTo("-plotHistograms")==0) {
	264	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	265	HLTInfo("plotHistograms is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	266	fPlotHistograms=((TObjString*)pTokens->At(i))->GetString().Atoi();
	267	continue;
	268	}
4d5ee3db	269	else if (argument.CompareTo("-constrainedTrackDeviation")==0) {
	270	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	271	HLTInfo("constrainedTrackDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	272	fConstrainedTrackDeviation=((TObjString*)pTokens->At(i))->GetString().Atof();
	273	continue;
	274	}
	275	else if (argument.CompareTo("-v0DaughterPrimDeviation")==0) {
	276	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	277	HLTInfo("v0DaughterPrimDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	278	fV0DaughterPrimDeviation=((TObjString*)pTokens->At(i))->GetString().Atof();
	279	continue;
	280	}
	281	else if (argument.CompareTo("-v0PrimDeviation")==0) {
	282	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	283	HLTInfo("v0PrimDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	284	fV0PrimDeviation=((TObjString*)pTokens->At(i))->GetString().Atof();
	285	continue;
	286	}
	287	else if (argument.CompareTo("-v0Chi")==0) {
	288	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	289	HLTInfo("v0Chi is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	290	fV0Chi=((TObjString*)pTokens->At(i))->GetString().Atof();
	291	continue;
	292	}
	293	else if (argument.CompareTo("-v0DecayLengthInSigmas")==0) {
	294	if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
	295	HLTInfo("v0DecayLengthInSigmas is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
	296	fV0DecayLengthInSigmas=((TObjString*)pTokens->At(i))->GetString().Atof();
	297	continue;
	298	}
	299	else {
	300	HLTError("unknown argument %s", argument.Data());
	301	iResult=-EINVAL;
	302	break;
	303	}
	304	}
	305	delete pTokens;
	306	}
	307	if (bMissingParam) {
	308	HLTError("missing parameter for argument %s", argument.Data());
	309	iResult=-EINVAL;
	310	}
	311
	312	return iResult;
	313	}
	314
	315	int AliHLTGlobalVertexerComponent::Reconfigure(const char* cdbEntry, const char* chainId)
	316	{
	317	// see header file for class documentation
	318	int iResult=0;
	319	const char* path="HLT/ConfigTPC/KryptonHistoComponent";
	320	const char* defaultNotify="";
	321	if (cdbEntry) {
	322	path=cdbEntry;
	323	defaultNotify=" (default)";
	324	}
	325	if (path) {
	326	HLTInfo("reconfigure from entry %s%s, chain id %s", path, defaultNotify,(chainId!=NULL && chainId[0]!=0)?chainId:"<none>");
	327	AliCDBEntry pEntry = AliCDBManager::Instance()->Get(path/,GetRunNo()*/);
	328	if (pEntry) {
	329	TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject());
	330	if (pString) {
	331	HLTInfo("received configuration object string: \'%s\'", pString->GetString().Data());
	332	iResult=Configure(pString->GetString().Data());
333	} else {
334	HLTError("configuration object \"%s\" has wrong type, required TObjString", path);
335	}
336	} else {
337	HLTError("can not fetch object \"%s\" from CDB", path);
338	}
339	}
340
341	return iResult;
342	}
343
344
345
346	void AliHLTGlobalVertexerComponent::SetESD( AliESDEvent *event )
347	{
348	//* Fill fTrackInfo array
349
350	delete[] fTrackInfos;
351	fESD = event;
352
353	AliKFParticle::SetField( fESD->GetMagneticField() );
354
355	Int_t nESDTracks=event->GetNumberOfTracks();
356	fTrackInfos = new AliESDTrackInfo[ nESDTracks ];
357
358	for (Int_t iTr=0; iTr<nESDTracks; iTr++){
359
360	AliESDTrackInfo &info = fTrackInfos[iTr];
361	info.fOK = 0;
362	info.fPrimUsedFlag = 0;
363
364	//* track quality check
365
366	AliESDtrack *pTrack = event->GetTrack(iTr);
367	if( !pTrack ) continue;
368	if (pTrack->GetKinkIndex(0)>0) continue;
369	if ( !( pTrack->GetStatus()&AliESDtrack::kTPCin ) ) continue;
370
371	//* Construct KFParticle for the track
372
373	info.fParticle = AliKFParticle( *pTrack->GetInnerParam(), 211 );
374	info.fOK = 1;
375	}
376	}
377
378
379	void AliHLTGlobalVertexerComponent::FindPrimaryVertex( )
380	{
381	//* Find event primary vertex
382
383	int nTracks = fESD->GetNumberOfTracks();
384
385	const AliKFParticle *vSelected = new const AliKFParticle[nTracks]; //* Selected particles for vertex fit
386	Int_t vIndex = new int [nTracks]; // Indices of selected particles
387	Bool_t vFlag = new bool [nTracks]; // Flags returned by the vertex finder
388
389	fPrimaryVtx.Initialize();
390	fPrimaryVtx.SetBeamConstraint(fESD->GetDiamondX(),fESD->GetDiamondY(),0,
391	TMath::Sqrt(fESD->GetSigma2DiamondX()),TMath::Sqrt(fESD->GetSigma2DiamondY()),5.3);
392
393	Int_t nSelected = 0;
394	for( Int_t i = 0; i<nTracks; i++){
395	if(!fTrackInfos[i].fOK ) continue;
396	//if( fESD->GetTrack(i)->GetTPCNcls()<60 ) continue;
397	const AliKFParticle &p = fTrackInfos[i].fParticle;
398	Double_t chi = p.GetDeviationFromVertex( fPrimaryVtx );
399	if( chi > fConstrainedTrackDeviation ) continue;
400	vSelected[nSelected] = &(fTrackInfos[i].fParticle);
401	vIndex[nSelected] = i;
402	nSelected++;
403	}
404	fPrimaryVtx.ConstructPrimaryVertex( vSelected, nSelected, vFlag, fConstrainedTrackDeviation );
405	for( Int_t i = 0; i<nSelected; i++){
406	if( vFlag[i] ) fTrackInfos[vIndex[i]].fPrimUsedFlag = 1;
407	}
408
409	for( Int_t i = 0; i<nTracks; i++ ){
410	AliESDTrackInfo &info = fTrackInfos[i];
411	info.fPrimDeviation = info.fParticle.GetDeviationFromVertex( fPrimaryVtx );
412	}
413
414	if( fPrimaryVtx.GetNContributors()>=3 ){
415	AliESDVertex vESD( fPrimaryVtx.Parameters(), fPrimaryVtx.CovarianceMatrix(), fPrimaryVtx.GetChi2(), fPrimaryVtx.GetNContributors() );
416	fESD->SetPrimaryVertexTracks( &vESD );
417
418	// relate the tracks to vertex
419
608cfbda	420	if( fFitTracksToVertex ){
4d5ee3db	421	for( Int_t i = 0; i<nTracks; i++ ){
	422	if( !fTrackInfos[i].fPrimUsedFlag ) continue;
	423	if( fTrackInfos[i].fPrimDeviation > fConstrainedTrackDeviation ) continue;
	424	fESD->GetTrack(i)->RelateToVertex( &vESD, fESD->GetMagneticField(),100. );
	425	}
	426	}
	427
	428	} else {
	429	for( Int_t i = 0; i<nTracks; i++)
	430	fTrackInfos[i].fPrimUsedFlag = 0;
	431	}
	432
	433
	434	delete[] vSelected;
	435	delete[] vIndex;
	436	delete[] vFlag;
	437	}
	438
	439
	440	void AliHLTGlobalVertexerComponent::FindV0s( )
	441	{
	442	//* V0 finder
	443
	444	int nTracks = fESD->GetNumberOfTracks();
	445	//AliKFVertex primVtx( *fESD->GetPrimaryVertexTracks() );
	446	AliKFVertex &primVtx = fPrimaryVtx;
	447	if( primVtx.GetNContributors()<3 ) return;
	448
	449	bool *constrainedV0 = new bool[nTracks];
	450	for( Int_t iTr = 0; iTr<nTracks; iTr++ ){
	451	constrainedV0[iTr] = 0;
	452	}
	453
	454
	455	for( Int_t iTr = 0; iTr<nTracks; iTr++ ){ //* first daughter
	456
	457	AliESDTrackInfo &info = fTrackInfos[iTr];
	458	if( !info.fOK ) continue;
	459	if( info.fParticle.GetQ() >0 ) continue;
	460	if( info.fPrimDeviation < fV0DaughterPrimDeviation ) continue;
	461
	462	for( Int_t jTr = 0; jTr<nTracks; jTr++ ){ //* second daughter
	463
	464	AliESDTrackInfo &jnfo = fTrackInfos[jTr];
	465	if( !jnfo.fOK ) continue;
	466	if( jnfo.fParticle.GetQ() < 0 ) continue;
	467	if( jnfo.fPrimDeviation < fV0DaughterPrimDeviation ) continue;
	468
	469	//* construct V0 mother
	470
	471	AliKFParticle v0( info.fParticle, jnfo.fParticle );
	472
	473	//* check V0 Chi^2
	474
	475	if( v0.GetChi2()<0 \|\| v0.GetChi2() > fV0ChifV0Chiv0.GetNDF() ) continue;
	476
	477	//* subtruct daughters from primary vertex
	478
	479	AliKFVertex primVtxCopy = primVtx;
	480
	481	if( info.fPrimUsedFlag ){
	482	if( primVtxCopy.GetNContributors()<=2 ) continue;
	483	primVtxCopy -= info.fParticle;
	484	}
485	if( jnfo.fPrimUsedFlag ){
486	if( primVtxCopy.GetNContributors()<=2 ) continue;
487	primVtxCopy -= jnfo.fParticle;
488	}
489	//* Check v0 Chi^2 deviation from primary vertex
490
491	if( v0.GetDeviationFromVertex( primVtxCopy ) > fV0PrimDeviation ) continue;
492
493	//* Add V0 to primary vertex to improve the primary vertex resolution
494
495	primVtxCopy += v0;
496
497	//* Set production vertex for V0
498
499	v0.SetProductionVertex( primVtxCopy );
500
501	//* Get V0 decay length with estimated error
502
503	Double_t length, sigmaLength;
504	if( v0.GetDecayLength( length, sigmaLength ) ) continue;
505
506	//* Reject V0 if it decays too close[sigma] to the primary vertex
507
508	if( length < fV0DecayLengthInSigmas*sigmaLength ) continue;
509
510	//* add ESD v0
511
512	AliESDv0 v0ESD( fESD->GetTrack( iTr ), iTr, fESD->GetTrack( jTr ), jTr );
513	fESD->AddV0( &v0ESD );
514
515	// relate the tracks to vertex
516
608cfbda	517	if( fFitTracksToVertex ){
4d5ee3db	518	if( constrainedV0[iTr] \|\| constrainedV0[jTr]
	519	\|\| info.fPrimDeviation < fConstrainedTrackDeviation \|\| jnfo.fPrimDeviation < fConstrainedTrackDeviation ) continue;
	520	AliESDVertex vESD(v0.Parameters(), v0.CovarianceMatrix(), v0.GetChi2(), 2);
	521	fESD->GetTrack(iTr)->RelateToVertex( &vESD, fESD->GetMagneticField(),100. );
	522	fESD->GetTrack(jTr)->RelateToVertex( &vESD, fESD->GetMagneticField(),100. );
	523	constrainedV0[iTr] = 1;
	524	constrainedV0[jTr] = 1;
	525	}
	526	}
	527	}
	528	delete[] constrainedV0;
	529	}
	530