[u/mrichter/AliRoot.git] / MUON / AliMUONDigitizerv1.cxx



#include <Riostream.h>
#include <TDirectory.h>
#include <TFile.h>
#include <TObjArray.h>
#include <TPDGCode.h>
#include <TTree.h> 

#include "AliMUON.h"
#include "AliMUONChamber.h"
#include "AliMUONConstants.h"
#include "AliMUONDigit.h"
#include "AliMUONDigitizerv1.h"
#include "AliMUONHit.h"
#include "AliMUONHitMapA1.h"
#include "AliMUONPadHit.h"
#include "AliMUONTransientDigit.h"
#include "AliRun.h"
#include "AliRunDigitizer.h"

ClassImp(AliMUONDigitizerv1)

//___________________________________________
AliMUONDigitizerv1::AliMUONDigitizerv1() :AliDigitizer()
{
// Default ctor - don't use it
  fHitMap = 0;
  fTDList = 0;
  if (GetDebug()>2) 
    cerr<<"AliMUONDigitizerv1::AliMUONDigitizerv1"
	<<"(AliRunDigitizer* manager) was processed"<<endl;
}

//___________________________________________
AliMUONDigitizerv1::AliMUONDigitizerv1(AliRunDigitizer* manager) 
    :AliDigitizer(manager)
{
// ctor which should be used
  fHitMap  = 0;
  fTDList  = 0;
  fDebug   = 0; 
  fHits = new TClonesArray("AliMUONHit",1000);
  if (GetDebug()>2) 
    cerr<<"AliMUONDigitizerv1::AliMUONDigitizerv1"
	<<"(AliRunDigitizer* manager) was processed"<<endl;
}

//------------------------------------------------------------------------
AliMUONDigitizerv1::~AliMUONDigitizerv1()
{
// Destructor
  delete fHits;
}

//------------------------------------------------------------------------
void AliMUONDigitizerv1::AddTransientDigit(AliMUONTransientDigit * mTD)
{
  // Choosing the maping of the cathode plane of the chamber:
  Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
  fTDList->AddAtAndExpand(mTD, fTDCounter);
  fHitMap[iNchCpl]->SetHit( mTD->PadX(), mTD->PadY(), fTDCounter);
  fTDCounter++;
}

//------------------------------------------------------------------------
Bool_t AliMUONDigitizerv1::ExistTransientDigit(AliMUONTransientDigit * mTD)
{
  // Choosing the maping of the cathode plane of the chamber:
  Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
  return( fHitMap[iNchCpl]->TestHit(mTD->PadX(), mTD->PadY()) );
}

//------------------------------------------------------------------------
Bool_t AliMUONDigitizerv1::Init()
{
// Initialization 
  return kTRUE;
}

//------------------------------------------------------------------------
void AliMUONDigitizerv1::UpdateTransientDigit(Int_t track, AliMUONTransientDigit * mTD)
{
  // Choosing the maping of the cathode plane of the chamber:
  Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
  AliMUONTransientDigit *pdigit = 
    static_cast<AliMUONTransientDigit*>(fHitMap[iNchCpl]->GetHit(mTD->PadX(),mTD->PadY()));
  // update charge
  //
  Int_t iqpad  = mTD->Signal();        // charge per pad
  pdigit->AddSignal(iqpad);
  pdigit->AddPhysicsSignal(iqpad);		
  // update list of tracks
  //
  Int_t charge;    
  track=+ fMask;
  if (fSignal) charge = iqpad;
  else         charge = kBgTag;
  pdigit->UpdateTrackList(track,charge);
}


//--------------------------------------------------------------------------
void AliMUONDigitizerv1::MakeTransientDigit(Int_t track, Int_t iHit, AliMUONHit * mHit)
{
  AliMUON *pMUON  = (AliMUON *) gAlice->GetModule("MUON");
  if (!pMUON) { 
    cerr<<"AliMUONDigitizerv1::Digitize Error:"
	<<" module MUON not found in the input file"<<endl;
  }
  if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit starts"<<endl;
  Int_t   ichamber = mHit->Chamber()-1;
  AliMUONChamber & chamber = pMUON->Chamber(ichamber);
  Float_t xhit = mHit->X();
  Float_t yhit = mHit->Y();
  Float_t zhit = mHit->Z();
  Float_t eloss= mHit->Eloss(); 
  Float_t tof  = mHit->Age();
  // Variables for chamber response from AliMUONChamber::DisIntegration
  Float_t newdigit[6][500];  // Pad information
  Int_t nnew=0;              // Number of touched Pads per hit
  Int_t digits[6];
  
  //
  // Calls the charge disintegration method of the current chamber 
  if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit calling AliMUONChamber::DisIngtegration starts"<<endl;
  chamber.DisIntegration(eloss, tof, xhit, yhit, zhit, nnew, newdigit);
  // Creating a new TransientDigits from hit
  for(Int_t iTD=0; iTD<nnew; iTD++) {
    digits[0] = Int_t(newdigit[1][iTD]);  // Padx of the Digit
    digits[1] = Int_t(newdigit[2][iTD]);  // Pady of the Digit
    digits[2] = Int_t(newdigit[5][iTD]);  // Cathode plane
    digits[3] = Int_t(newdigit[3][iTD]);  // Induced charge in the Pad
    if (fSignal) digits[4] = Int_t(newdigit[3][iTD]);
    else         digits[4] = 0;
    digits[5] = iHit+fMask;    // Hit number in the list
    if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit " <<
			"PadX "<< digits[0] << " " <<
			"PadY "<< digits[1] << " " <<
			"Plane " << digits[2] << " " <<
			"Charge " << digits[3] <<" " <<
			"Hit " << digits[5] << endl;
    // list of tracks
    Int_t charge;    
    track += fMask;
    if (fSignal) charge = digits[3];
    else         charge = kBgTag;
    if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit Creating AliMUONTransientDigit"<<endl;
    AliMUONTransientDigit * mTD = new AliMUONTransientDigit(ichamber, digits);
    mTD->AddToTrackList(track,charge);

    if (!ExistTransientDigit(mTD)) {
      AddTransientDigit(mTD);
      if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit Adding TransientDigit"<<endl;
    }
    else {
      if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit updating TransientDigit"<<endl;
      UpdateTransientDigit(track, mTD);
      delete mTD;
    }
  }
}
//-----------------------------------------------------------------------
void AliMUONDigitizerv1::Exec(Option_t* option)
{
  TString optionString = option;
  if (optionString.Data() == "deb") {
    cout<<"AliMUONDigitizerv1::Exec: called with option deb "<<endl;
    fDebug = 3;
  }

  AliMUONChamber*   chamber;
  AliSegmentation*  c1Segmentation; //Cathode plane c1 of the chamber
  AliSegmentation*  c2Segmentation; //Cathode place c2 of the chamber

  if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::Digitize() starts"<<endl;
  fTDList = new TObjArray;

  // Getting Module MUON
  AliMUON *pMUON  = (AliMUON *) gAlice->GetModule("MUON");
  if (!pMUON) {
    cerr<<"AliMUONDigitizerv1::Digitize Error:"
	<<" module MUON not found in the input file"<<endl;
    return;
  }
  // New branch for MUON digit in the tree of digits
  pMUON->MakeBranchInTreeD(fManager->GetTreeD());

  // Array of pointer of the AliMUONHitMapA1:
  //  two HitMaps per chamber, or one HitMap per cahtode plane
  fHitMap= new AliMUONHitMapA1* [2*AliMUONConstants::NCh()];

  //Loop over chambers for the definition AliMUONHitMap
  for (Int_t i=0; i<AliMUONConstants::NCh(); i++) {
    chamber = &(pMUON->Chamber(i));
    c1Segmentation = chamber->SegmentationModel(1); // Cathode plane 1
    fHitMap[i] = new AliMUONHitMapA1(c1Segmentation, fTDList);
    c2Segmentation = chamber->SegmentationModel(2); // Cathode plane 2
    fHitMap[i+AliMUONConstants::NCh()] = new AliMUONHitMapA1(c2Segmentation, fTDList);
  }

// Loop over files to merge and to digitize
    fSignal = kTRUE;
    for (Int_t inputFile=0; inputFile<fManager->GetNinputs(); inputFile++) {
      // Connect MUON Hit branch
      if (inputFile > 0 ) fSignal = kFALSE;
      TBranch *branchHits = 0;
      TTree *treeH = fManager->GetInputTreeH(inputFile);
      if (GetDebug()>2) {
	cerr<<"AliMUONDigitizerv1::Exec inputFile is "<<inputFile<<" "<<endl;
	cerr<<"AliMUONDigitizerv1::Exec treeH, fHits "<<treeH<<" "<<fHits<<endl;
      }
      if (treeH && fHits) {
	branchHits = treeH->GetBranch("MUON");
	if (branchHits) {
	  fHits->Delete();
	  branchHits->SetAddress(&fHits);
	}
	else
	  Error("Exec","branch MUON was not found");
      }
      if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::Exec branchHits = "<<branchHits<<endl;

      fMask = fManager->GetMask(inputFile);
      //
      // Loop over tracks
      Int_t itrack;
      Int_t ntracks = (Int_t) treeH->GetEntries();
      for (itrack = 0; itrack < ntracks; itrack++) {
	if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::Exec itrack = "<<itrack<<endl;
	fHits->Delete();
	branchHits->GetEntry(itrack);
	//
	//  Loop over hits
	Int_t ihit, ichamber;
	AliMUONHit* mHit;
	for(ihit = 0; ihit < fHits->GetEntriesFast(); ihit++) {
	  mHit = static_cast<AliMUONHit*>(fHits->At(ihit));
	  ichamber = mHit->Chamber()-1;  // chamber number
	  if (ichamber > AliMUONConstants::NCh()-1) {
	    cerr<<"AliMUONDigitizer: ERROR: "
		<<"fNch > AliMUONConstants::NCh()-1, fNch, NCh(): "
		<<ichamber<<", "<< AliMUONConstants::NCh()<<endl;
	    return;
	  }
	  chamber = &(pMUON->Chamber(ichamber));
	  //
	  //Dumping Hit content:
	  if (GetDebug()>2) {
	    cerr<<"AliMuonDigitizerv1::Exec ihit, ichamber, x, y, z, eloss " <<
	      ihit << " " << 
	      mHit->Chamber() << " " <<
	      mHit->X() << " " <<
	      mHit->Y() << " " <<
	      mHit->Z() << " " <<
	      mHit->Eloss() << " " << endl;
	  }
	  // 
	  // Inititializing Correlation
	  chamber->ChargeCorrelationInit();
	  if (ichamber < AliMUONConstants::NTrackingCh()) {
	    // Tracking Chamber
	    // Initialize hit position (cursor) in the segmentation model 
	    chamber->SigGenInit(mHit->X(), mHit->Y(), mHit->Z());
	  } else {
	    // Trigger Chamber 
	  }
	  MakeTransientDigit(itrack, ihit, mHit);
	} // hit loop
      } // track loop
    } // end file loop
    if (GetDebug()>2) cerr<<"AliMUONDigitizer::Exec End of hits, track and file loops"<<endl;

    //
    // Filling Digit List
    Int_t tracks[kMAXTRACKS];
    Int_t charges[kMAXTRACKS];
    Int_t nentries = fTDList->GetEntriesFast();
    Int_t digits[6];
    for (Int_t nent = 0; nent < nentries; nent++) {
      AliMUONTransientDigit *address = (AliMUONTransientDigit*)fTDList->At(nent);
      if (address == 0) continue; 
      Int_t ich = address->Chamber();
      Int_t   q = address->Signal(); 
      chamber = &(pMUON->Chamber(ich));
      //
      //  Digit Response (noise, threshold, saturation, ...)
      AliMUONResponse * response = chamber->ResponseModel();
      q = response->DigitResponse(q,address);
      
      if (!q) continue;
      
      digits[0] = address->PadX();
      digits[1] = address->PadY();
      digits[2] = address->Cathode();
      digits[3] = q;
      digits[4] = address->Physics();
      digits[5] = address->Hit();
      
      Int_t nptracks = address->GetNTracks();
      
      if (nptracks > kMAXTRACKS) {
	if (GetDebug() >0) {
	  cerr<<"AliMUONDigitizer:Exec  nptracks > 10 "<<nptracks;
	  cerr<<"reset to max value "<<kMAXTRACKS<<endl;
	}
	nptracks = kMAXTRACKS;
      }
      if (nptracks > 2 && GetDebug() >2) {
	cerr<<"AliMUONDigitizer::Exec  nptracks > 2 "<<nptracks<<endl;
	//	printf("cat,ich,ix,iy,q %d %d %d %d %d \n",icat,ich,digits[0],digits[1],q);
      }
      for (Int_t tr = 0; tr < nptracks; tr++) {
	tracks[tr]   = address->GetTrack(tr);
	charges[tr]  = address->GetCharge(tr);
      }      //end loop over list of tracks for one pad
      // Sort list of tracks according to charge
      if (nptracks > 1) {
	SortTracks(tracks,charges,nptracks);
      }
      if (nptracks < kMAXTRACKS ) {
	for (Int_t i = nptracks; i < kMAXTRACKS; i++) {
	  tracks[i]  = 0;
	  charges[i] = 0;
	}
      }
      
      // fill digits
      if (GetDebug()>2) cerr<<"AliMUONDigitzerv1::Exex TransientDigit to Digit"<<endl;
      pMUON->AddDigits(ich,tracks,charges,digits);
    }
    fManager->GetTreeD()->Fill();
    pMUON->ResetDigits();  //
    fTDList->Delete();
    
    for(Int_t ii = 0; ii < 2*AliMUONConstants::NCh(); ++ii) {
      if (fHitMap[ii]) {
	delete fHitMap[ii];
	fHitMap[ii] = 0;
      }
    }
    
    if (GetDebug()>2) 
      cerr<<"AliMUONDigitizer::Exec: writing the TreeD: "
	  <<fManager->GetTreeD()->GetName()<<endl;
    fManager->GetTreeD()->Write(0,TObject::kOverwrite);
    delete [] fHitMap;
    delete fTDList;
    
    if (fHits)    fHits->Clear();
}

//------------------------------------------------------------------------
void AliMUONDigitizerv1::SortTracks(Int_t *tracks,Int_t *charges,Int_t ntr)
{
  //
  // Sort the list of tracks contributing to a given digit
  // Only the 3 most significant tracks are acctually sorted
  //
  
  //
  //  Loop over signals, only 3 times
  //
  
  Int_t qmax;
  Int_t jmax;
  Int_t idx[3] = {-2,-2,-2};
  Int_t jch[3] = {-2,-2,-2};
  Int_t jtr[3] = {-2,-2,-2};
  Int_t i,j,imax;
  
  if (ntr<3) imax=ntr;
  else imax=3;
  for(i=0;i<imax;i++){
    qmax=0;
    jmax=0;
    
    for(j=0;j<ntr;j++){
      
      if((i == 1 && j == idx[i-1]) 
	 ||(i == 2 && (j == idx[i-1] || j == idx[i-2]))) continue;
      
      if(charges[j] > qmax) {
	qmax = charges[j];
	jmax=j;
      }       
    } 
    
    if(qmax > 0) {
      idx[i]=jmax;
      jch[i]=charges[jmax]; 
      jtr[i]=tracks[jmax]; 
    }
    
  } 
  
  for(i=0;i<3;i++){
    if (jtr[i] == -2) {
      charges[i]=0;
      tracks[i]=0;
    } else {
      charges[i]=jch[i];
      tracks[i]=jtr[i];
    }
  }
}
Commit	Line	Data
8dbbc4e3	1
	2
	3	#include <Riostream.h>
	4	#include <TDirectory.h>
	5	#include <TFile.h>
	6	#include <TObjArray.h>
	7	#include <TPDGCode.h>
	8	#include <TTree.h>
	9
	10	#include "AliMUON.h"
	11	#include "AliMUONChamber.h"
	12	#include "AliMUONConstants.h"
	13	#include "AliMUONDigit.h"
	14	#include "AliMUONDigitizerv1.h"
	15	#include "AliMUONHit.h"
	16	#include "AliMUONHitMapA1.h"
	17	#include "AliMUONPadHit.h"
	18	#include "AliMUONTransientDigit.h"
	19	#include "AliRun.h"
	20	#include "AliRunDigitizer.h"
	21
	22	ClassImp(AliMUONDigitizerv1)
	23
	24	//___________________________________________
	25	AliMUONDigitizerv1::AliMUONDigitizerv1() :AliDigitizer()
	26	{
	27	// Default ctor - don't use it
	28	fHitMap = 0;
	29	fTDList = 0;
	30	if (GetDebug()>2)
	31	cerr<<"AliMUONDigitizerv1::AliMUONDigitizerv1"
	32	<<"(AliRunDigitizer* manager) was processed"<<endl;
	33	}
	34
	35	//___________________________________________
	36	AliMUONDigitizerv1::AliMUONDigitizerv1(AliRunDigitizer* manager)
	37	:AliDigitizer(manager)
	38	{
	39	// ctor which should be used
	40	fHitMap = 0;
	41	fTDList = 0;
	42	fDebug = 0;
	43	fHits = new TClonesArray("AliMUONHit",1000);
	44	if (GetDebug()>2)
	45	cerr<<"AliMUONDigitizerv1::AliMUONDigitizerv1"
	46	<<"(AliRunDigitizer* manager) was processed"<<endl;
	47	}
	48
	49	//------------------------------------------------------------------------
	50	AliMUONDigitizerv1::~AliMUONDigitizerv1()
	51	{
	52	// Destructor
	53	delete fHits;
	54	}
	55
	56	//------------------------------------------------------------------------
	57	void AliMUONDigitizerv1::AddTransientDigit(AliMUONTransientDigit * mTD)
	58	{
	59	// Choosing the maping of the cathode plane of the chamber:
	60	Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
	61	fTDList->AddAtAndExpand(mTD, fTDCounter);
	62	fHitMap[iNchCpl]->SetHit( mTD->PadX(), mTD->PadY(), fTDCounter);
	63	fTDCounter++;
	64	}
65
66	//------------------------------------------------------------------------
67	Bool_t AliMUONDigitizerv1::ExistTransientDigit(AliMUONTransientDigit * mTD)
68	{
69	// Choosing the maping of the cathode plane of the chamber:
70	Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
71	return( fHitMap[iNchCpl]->TestHit(mTD->PadX(), mTD->PadY()) );
72	}
73
74	//------------------------------------------------------------------------
75	Bool_t AliMUONDigitizerv1::Init()
76	{
77	// Initialization
78	return kTRUE;
79	}
80
81	//------------------------------------------------------------------------
82	void AliMUONDigitizerv1::UpdateTransientDigit(Int_t track, AliMUONTransientDigit * mTD)
83	{
84	// Choosing the maping of the cathode plane of the chamber:
85	Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
86	AliMUONTransientDigit *pdigit =
87	static_cast<AliMUONTransientDigit*>(fHitMap[iNchCpl]->GetHit(mTD->PadX(),mTD->PadY()));
88	// update charge
89	//
90	Int_t iqpad = mTD->Signal(); // charge per pad
91	pdigit->AddSignal(iqpad);
92	pdigit->AddPhysicsSignal(iqpad);
93	// update list of tracks
94	//
95	Int_t charge;
96	track=+ fMask;
97	if (fSignal) charge = iqpad;
98	else charge = kBgTag;
99	pdigit->UpdateTrackList(track,charge);
100	}
101
102
103	//--------------------------------------------------------------------------
104	void AliMUONDigitizerv1::MakeTransientDigit(Int_t track, Int_t iHit, AliMUONHit * mHit)
105	{
106	AliMUON pMUON = (AliMUON ) gAlice->GetModule("MUON");
107	if (!pMUON) {
108	cerr<<"AliMUONDigitizerv1::Digitize Error:"
109	<<" module MUON not found in the input file"<<endl;
110	}
111	if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit starts"<<endl;
112	Int_t ichamber = mHit->Chamber()-1;
113	AliMUONChamber & chamber = pMUON->Chamber(ichamber);
114	Float_t xhit = mHit->X();
115	Float_t yhit = mHit->Y();
116	Float_t zhit = mHit->Z();
117	Float_t eloss= mHit->Eloss();
118	Float_t tof = mHit->Age();
119	// Variables for chamber response from AliMUONChamber::DisIntegration
120	Float_t newdigit[6][500]; // Pad information
121	Int_t nnew=0; // Number of touched Pads per hit
122	Int_t digits[6];
123
124	//
125	// Calls the charge disintegration method of the current chamber
126	if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit calling AliMUONChamber::DisIngtegration starts"<<endl;
127	chamber.DisIntegration(eloss, tof, xhit, yhit, zhit, nnew, newdigit);
128	// Creating a new TransientDigits from hit
129	for(Int_t iTD=0; iTD<nnew; iTD++) {
130	digits[0] = Int_t(newdigit[1][iTD]); // Padx of the Digit
131	digits[1] = Int_t(newdigit[2][iTD]); // Pady of the Digit
132	digits[2] = Int_t(newdigit[5][iTD]); // Cathode plane
133	digits[3] = Int_t(newdigit[3][iTD]); // Induced charge in the Pad
134	if (fSignal) digits[4] = Int_t(newdigit[3][iTD]);
135	else digits[4] = 0;
136	digits[5] = iHit+fMask; // Hit number in the list
137	if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit " <<
138	"PadX "<< digits[0] << " " <<
139	"PadY "<< digits[1] << " " <<
140	"Plane " << digits[2] << " " <<
141	"Charge " << digits[3] <<" " <<
142	"Hit " << digits[5] << endl;
143	// list of tracks
144	Int_t charge;
145	track += fMask;
146	if (fSignal) charge = digits[3];
147	else charge = kBgTag;
148	if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit Creating AliMUONTransientDigit"<<endl;
149	AliMUONTransientDigit * mTD = new AliMUONTransientDigit(ichamber, digits);
150	mTD->AddToTrackList(track,charge);
151
152	if (!ExistTransientDigit(mTD)) {
153	AddTransientDigit(mTD);
154	if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit Adding TransientDigit"<<endl;
155	}
156	else {
157	if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit updating TransientDigit"<<endl;
158	UpdateTransientDigit(track, mTD);
159	delete mTD;
160	}
161	}
162	}
163	//-----------------------------------------------------------------------
164	void AliMUONDigitizerv1::Exec(Option_t* option)
165	{
166	TString optionString = option;
167	if (optionString.Data() == "deb") {
168	cout<<"AliMUONDigitizerv1::Exec: called with option deb "<<endl;
169	fDebug = 3;
170	}
171
172	AliMUONChamber* chamber;
173	AliSegmentation* c1Segmentation; //Cathode plane c1 of the chamber
174	AliSegmentation* c2Segmentation; //Cathode place c2 of the chamber
175
176	if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::Digitize() starts"<<endl;
177	fTDList = new TObjArray;
178
179	// Getting Module MUON
180	AliMUON pMUON = (AliMUON ) gAlice->GetModule("MUON");
181	if (!pMUON) {
182	cerr<<"AliMUONDigitizerv1::Digitize Error:"
183	<<" module MUON not found in the input file"<<endl;
184	return;
185	}
186	// New branch for MUON digit in the tree of digits
187	pMUON->MakeBranchInTreeD(fManager->GetTreeD());
188
189	// Array of pointer of the AliMUONHitMapA1:
190	// two HitMaps per chamber, or one HitMap per cahtode plane
191	fHitMap= new AliMUONHitMapA1* [2*AliMUONConstants::NCh()];
192
193	//Loop over chambers for the definition AliMUONHitMap
194	for (Int_t i=0; i<AliMUONConstants::NCh(); i++) {
195	chamber = &(pMUON->Chamber(i));
196	c1Segmentation = chamber->SegmentationModel(1); // Cathode plane 1
197	fHitMap[i] = new AliMUONHitMapA1(c1Segmentation, fTDList);
198	c2Segmentation = chamber->SegmentationModel(2); // Cathode plane 2
199	fHitMap[i+AliMUONConstants::NCh()] = new AliMUONHitMapA1(c2Segmentation, fTDList);
200	}
201
202	// Loop over files to merge and to digitize
203	fSignal = kTRUE;
204	for (Int_t inputFile=0; inputFile<fManager->GetNinputs(); inputFile++) {
205	// Connect MUON Hit branch
206	if (inputFile > 0 ) fSignal = kFALSE;
207	TBranch *branchHits = 0;
208	TTree *treeH = fManager->GetInputTreeH(inputFile);
209	if (GetDebug()>2) {
210	cerr<<"AliMUONDigitizerv1::Exec inputFile is "<<inputFile<<" "<<endl;
211	cerr<<"AliMUONDigitizerv1::Exec treeH, fHits "<<treeH<<" "<<fHits<<endl;
212	}
213	if (treeH && fHits) {
214	branchHits = treeH->GetBranch("MUON");
215	if (branchHits) {
216	fHits->Delete();
217	branchHits->SetAddress(&fHits);
218	}
219	else
220	Error("Exec","branch MUON was not found");
221	}
222	if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::Exec branchHits = "<<branchHits<<endl;
223
224	fMask = fManager->GetMask(inputFile);
225	//
226	// Loop over tracks
227	Int_t itrack;
228	Int_t ntracks = (Int_t) treeH->GetEntries();
229	for (itrack = 0; itrack < ntracks; itrack++) {
230	if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::Exec itrack = "<<itrack<<endl;
231	fHits->Delete();
232	branchHits->GetEntry(itrack);
233	//
234	// Loop over hits
235	Int_t ihit, ichamber;
236	AliMUONHit* mHit;
237	for(ihit = 0; ihit < fHits->GetEntriesFast(); ihit++) {
238	mHit = static_cast<AliMUONHit*>(fHits->At(ihit));
239	ichamber = mHit->Chamber()-1; // chamber number
240	if (ichamber > AliMUONConstants::NCh()-1) {
241	cerr<<"AliMUONDigitizer: ERROR: "
242	<<"fNch > AliMUONConstants::NCh()-1, fNch, NCh(): "
243	<<ichamber<<", "<< AliMUONConstants::NCh()<<endl;
244	return;
245	}
246	chamber = &(pMUON->Chamber(ichamber));
247	//
248	//Dumping Hit content:
249	if (GetDebug()>2) {
250	cerr<<"AliMuonDigitizerv1::Exec ihit, ichamber, x, y, z, eloss " <<
251	ihit << " " <<
252	mHit->Chamber() << " " <<
253	mHit->X() << " " <<
254	mHit->Y() << " " <<
255	mHit->Z() << " " <<
256	mHit->Eloss() << " " << endl;
257	}
258	//
259	// Inititializing Correlation
260	chamber->ChargeCorrelationInit();
261	if (ichamber < AliMUONConstants::NTrackingCh()) {
262	// Tracking Chamber
263	// Initialize hit position (cursor) in the segmentation model
264	chamber->SigGenInit(mHit->X(), mHit->Y(), mHit->Z());
265	} else {
266	// Trigger Chamber
267	}
268	MakeTransientDigit(itrack, ihit, mHit);
269	} // hit loop
270	} // track loop
271	} // end file loop
272	if (GetDebug()>2) cerr<<"AliMUONDigitizer::Exec End of hits, track and file loops"<<endl;
273
274	//
275	// Filling Digit List
276	Int_t tracks[kMAXTRACKS];
277	Int_t charges[kMAXTRACKS];
278	Int_t nentries = fTDList->GetEntriesFast();
279	Int_t digits[6];
280	for (Int_t nent = 0; nent < nentries; nent++) {
281	AliMUONTransientDigit address = (AliMUONTransientDigit)fTDList->At(nent);
282	if (address == 0) continue;
283	Int_t ich = address->Chamber();
284	Int_t q = address->Signal();
285	chamber = &(pMUON->Chamber(ich));
286	//
287	// Digit Response (noise, threshold, saturation, ...)
288	AliMUONResponse * response = chamber->ResponseModel();
289	q = response->DigitResponse(q,address);
290
291	if (!q) continue;
292
293	digits[0] = address->PadX();
294	digits[1] = address->PadY();
295	digits[2] = address->Cathode();
296	digits[3] = q;
297	digits[4] = address->Physics();
298	digits[5] = address->Hit();
299
300	Int_t nptracks = address->GetNTracks();
301
302	if (nptracks > kMAXTRACKS) {
303	if (GetDebug() >0) {
304	cerr<<"AliMUONDigitizer:Exec nptracks > 10 "<<nptracks;
305	cerr<<"reset to max value "<<kMAXTRACKS<<endl;
306	}
307	nptracks = kMAXTRACKS;
308	}
309	if (nptracks > 2 && GetDebug() >2) {
310	cerr<<"AliMUONDigitizer::Exec nptracks > 2 "<<nptracks<<endl;
311	// printf("cat,ich,ix,iy,q %d %d %d %d %d \n",icat,ich,digits[0],digits[1],q);
312	}
313	for (Int_t tr = 0; tr < nptracks; tr++) {
314	tracks[tr] = address->GetTrack(tr);
315	charges[tr] = address->GetCharge(tr);
316	} //end loop over list of tracks for one pad
317	// Sort list of tracks according to charge
318	if (nptracks > 1) {
319	SortTracks(tracks,charges,nptracks);
320	}
321	if (nptracks < kMAXTRACKS ) {
322	for (Int_t i = nptracks; i < kMAXTRACKS; i++) {
323	tracks[i] = 0;
324	charges[i] = 0;
325	}
326	}
327
328	// fill digits
329	if (GetDebug()>2) cerr<<"AliMUONDigitzerv1::Exex TransientDigit to Digit"<<endl;
330	pMUON->AddDigits(ich,tracks,charges,digits);
331	}
332	fManager->GetTreeD()->Fill();
333	pMUON->ResetDigits(); //
334	fTDList->Delete();
335
336	for(Int_t ii = 0; ii < 2*AliMUONConstants::NCh(); ++ii) {
337	if (fHitMap[ii]) {
338	delete fHitMap[ii];
339	fHitMap[ii] = 0;
340	}
341	}
342
343	if (GetDebug()>2)
344	cerr<<"AliMUONDigitizer::Exec: writing the TreeD: "
345	<<fManager->GetTreeD()->GetName()<<endl;
346	fManager->GetTreeD()->Write(0,TObject::kOverwrite);
347	delete [] fHitMap;
348	delete fTDList;
349
350	if (fHits) fHits->Clear();
351	}
352
353	//------------------------------------------------------------------------
354	void AliMUONDigitizerv1::SortTracks(Int_t tracks,Int_t charges,Int_t ntr)
355	{
356	//
357	// Sort the list of tracks contributing to a given digit
358	// Only the 3 most significant tracks are acctually sorted
359	//
360
361	//
362	// Loop over signals, only 3 times
363	//
364
365	Int_t qmax;
366	Int_t jmax;
367	Int_t idx[3] = {-2,-2,-2};
368	Int_t jch[3] = {-2,-2,-2};
369	Int_t jtr[3] = {-2,-2,-2};
370	Int_t i,j,imax;
371
372	if (ntr<3) imax=ntr;
373	else imax=3;
374	for(i=0;i<imax;i++){
375	qmax=0;
376	jmax=0;
377
378	for(j=0;j<ntr;j++){
379
380	if((i == 1 && j == idx[i-1])
381	\|\|(i == 2 && (j == idx[i-1] \|\| j == idx[i-2]))) continue;
382
383	if(charges[j] > qmax) {
384	qmax = charges[j];
385	jmax=j;
386	}
387	}
388
389	if(qmax > 0) {
390	idx[i]=jmax;
391	jch[i]=charges[jmax];
392	jtr[i]=tracks[jmax];
393	}
394
395	}
396
397	for(i=0;i<3;i++){
398	if (jtr[i] == -2) {
399	charges[i]=0;
400	tracks[i]=0;
401	} else {
402	charges[i]=jch[i];
403	tracks[i]=jtr[i];
404	}
405	}
406	}