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

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * 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.                  *
 **************************************************************************/

/* $Id$ */

// --------------------------
// Class AliMUONDetElement
// --------------------------
// Detection element object containing information for combined 
// cluster / track finder in the MUON arm 
// Author: Alexander Zinchenko, JINR Dubna
 
#include "AliMUONDetElement.h"
#include "AliMUON.h"
#include "AliMUONSegmentation.h"
#include "AliMUONDigit.h"
#include "AliMUONHitMapA1.h"
#include "AliMUONData.h"
#include "AliMUONRawCluster.h"
#include "AliMUONHitForRec.h"
#include "AliMUONClusterInput.h"
#include "AliMUONClusterFinderAZ.h"
#include "AliMUONGeometryModuleTransformer.h" //AZ
#include "AliMUONVGeometryDESegmentation.h" //AZ
#include "AliMpVSegmentation.h" //AZ

#include "AliRun.h"
#include "AliLog.h"

#include <TObjArray.h>
#include <TClonesArray.h>
#include <TVector2.h> //AZ

/// \cond CLASSIMP
ClassImp(AliMUONDetElement) // Class implementation in ROOT context
/// \endcond
  FILE *lun = 0x0; //fopen("hitmap.dat","w");

//_________________________________________________________________________
AliMUONDetElement::AliMUONDetElement()
  : TObject(),
    fidDE(0),
    fIndex(0),
    fChamber(0),
    fZ(0.),
    fNHitsForRec(0),
    fRawClus(0x0),
    fHitsForRec(0x0),
    fRecModel(0x0)
{
  /// Default constructor
  for (Int_t i = 0; i < 2; i++) {
    fHitMap[i] = NULL;
    fDigits[i] = NULL;
    fSeg[i] = NULL;
  }
} 

//_________________________________________________________________________
AliMUONDetElement::AliMUONDetElement(Int_t idDE, AliMUONDigit *dig, AliMUONClusterFinderAZ *recModel) 
  : TObject(),
    fidDE(idDE),
    fIndex(0),
    fChamber(idDE / 100 - 1),
    fZ(0.),
    fNHitsForRec(0),
    fRawClus(0x0),
    fHitsForRec(new TObjArray(10)),
    fRecModel(recModel)
{
  /// Constructor
  fDigits[0] = new TObjArray(10);
  fDigits[1] = new TObjArray(10);
  fRawClus = new TObjArray(10);
  AliMUON *pMUON = (AliMUON*) gAlice->GetModule("MUON");
  AliMUONSegmentation *pSegmentation = pMUON->GetSegmentation();
  fSeg[0] = pSegmentation->GetModuleSegmentationByDEId(idDE, 0);
  fSeg[1] = pSegmentation->GetModuleSegmentationByDEId(idDE, 1);
  /*
  Float_t x, y, z;
  fSeg[dig->Cathode()]->GetPadC(fidDE, dig->PadX(), dig->PadY(), x, y, z);
  fZ = z;
  */
  AddDigit(dig);
}

//_________________________________________________________________________
AliMUONDetElement::~AliMUONDetElement()
{
  /// Destructor

  for (Int_t i = 0; i < 2; i++) {
    delete fHitMap[i]; fHitMap[i] = NULL; 
    delete fDigits[i]; fDigits[i] = NULL; 
    //fDigits[i]->Delete(); delete fDigits[i]; fDigits[i] = NULL; 
  }
  //fHitsForRec->Clear(); delete fHitsForRec; fHitsForRec = 0;
  if (fRawClus) { fRawClus->Delete(); delete fRawClus; fRawClus = 0; }
  //if (fRawClus) { delete fRawClus; fRawClus = 0; }
  for (Int_t i = 0; i < fNHitsForRec; i++) {
    Int_t ntracks = ((AliMUONHitForRec*)fHitsForRec->UncheckedAt(i))->GetNTrackHits();
    //cout << ntracks << endl;
    if (ntracks) fHitsForRec->RemoveAt(i);
  }
  //fHitsForRec->Delete(); 
  delete fHitsForRec; fHitsForRec = 0;
  //fHitsForRec->SetOwner(kFALSE);
  //fHitsForRec->Clear(); cout << " here " << fHitsForRec->GetEntriesFast() << endl; delete fHitsForRec; fHitsForRec = 0;
}

//_________________________________________________________________________
Int_t AliMUONDetElement::Compare(const TObject* detElem) const
{
  /// "Compare" function to sort in Z (towards interaction point)
  /// Returns -1 (0, +1) if charge of current pixel
  /// is greater than (equal to, less than) charge of pixel
  if (fZ > ((AliMUONDetElement*)detElem)->Z()) return(+1);
  else if (fZ == ((AliMUONDetElement*)detElem)->Z()) return( 0);
  else return(-1);
}

//_________________________________________________________________________
void AliMUONDetElement::Fill(AliMUONData */*data*/)
{
  /// Fill hit maps
  fLeft[0] = fDigits[0]->GetEntriesFast();
  fLeft[1] = fDigits[1]->GetEntriesFast();

  Int_t  npx0  = fSeg[0]->Npx(fidDE)+1;
  Int_t  npy0  = fSeg[0]->Npy(fidDE)+1;
  fHitMap[0] = new AliMUONHitMapA1(npx0, npy0, fDigits[0]);
  fHitMap[0]->FillHits();

  Int_t  npx1  = fSeg[1]->Npx(fidDE)+1;
  Int_t  npy1  = fSeg[1]->Npy(fidDE)+1;
  fHitMap[1] = new AliMUONHitMapA1(npx1, npy1, fDigits[1]);
  fHitMap[1]->FillHits();

  // The part below is just for debugging (fill rec. points already found)
  /*
  fLeft[0] = fLeft[1] = 0;
  TClonesArray *rawClus = data->RawClusters(fChamber);
  cout << rawClus << " " << rawClus->GetEntriesFast() << endl;
  for (Int_t i = 0; i < rawClus->GetEntriesFast(); i++) {
    AliMUONRawCluster *recP = (AliMUONRawCluster*) rawClus->UncheckedAt(i);
    cout << fChamber << " " << recP->GetZ(0) << " " << recP->GetZ(1) << " " << fZ << endl;
    if (TMath::Abs(recP->GetZ(0)-fZ) > 0.5) continue;
    if (!Inside(recP->GetX(0), recP->GetY(0), recP->GetZ(0))) continue;
    AddHitForRec(recP); // add hit for rec.
    rawClus->RemoveAt(i); // remove
  }
  cout << fHitsForRec->GetEntriesFast() << endl;
  rawClus->Compress();
  */
}

//_________________________________________________________________________
void AliMUONDetElement::AddDigit(AliMUONDigit *dig)
{
  /// Add digit

  //fDigits[dig->Cathode()]->Add(dig);
  fDigits[dig->Cathode()]->Add(new AliMUONDigit(*dig));
  Float_t x, y, z;
  fSeg[dig->Cathode()]->GetPadC(fidDE, dig->PadX(), dig->PadY(), x, y, z);
  fZ += z;
  /*
  Int_t ndig = fDigits[dig->Cathode()]->GetEntriesFast();
  new((*fDigits[dig->Cathode()])[ndig]) AliMUONDigit(*dig);
  */
}

//_________________________________________________________________________
Bool_t AliMUONDetElement::Inside(Double_t x, Double_t y, Double_t z, Double_t dx, Double_t dy) const
{
  /// Check if point is inside detection element

  dx = TMath::Max (dx, 1.);
  dy = TMath::Max (dy, 1.);
  Int_t ix, iy;
  ix = iy = 0;
  Double_t xl, yl, zl;
  for (Int_t i = 0; i < 2; i++) {
    if (!fSeg[i]) continue;
    if (i == 0 || fSeg[0] == 0x0) {
      fSeg[i]->GetTransformer()->Global2Local(fidDE, x, y, z, xl, yl, zl);
      //cout << xl << " " << yl << " " << zl << endl;
      if (TMath::Abs (zl) > 0.5) return kFALSE;
      //if (fChamber < 4 && (xl < 0 || yl < 0)) return kFALSE;
    }
    fSeg[i]->GetPadI(fidDE, x, y, z, ix, iy);
    //cout << x << " " << y << " " << z << " " << fChamber << " " << ix << " " << iy << " " << fSeg[i]->Npx(fidDE) << " " << fSeg[i]->Npy(fidDE) /*<< " " << fSeg[i]->GetPadI(fidDE, x, y, z, ix, iy)*/ << " " << fSeg[i]->HasPad(fidDE, (Float_t)x, (Float_t)y, (Float_t)z) << endl; 
    //if (ix > 0 && iy > 0 && ix <= fSeg[i]->Npx(fidDE) && iy <= fSeg[i]->Npy(fidDE)) return kTRUE;
    if (fSeg[i]->HasPad(fidDE, (Float_t)x, (Float_t)y, (Float_t)z)) return kTRUE;
  }
  //Double_t xl, yl, zl;
  //fSeg[0]->GetTransformer()->Global2Local(fidDE, x, y, z, xl, yl, zl);
  //cout << xl << " " << yl << " " << zl << endl;

  Int_t cath = fSeg[0] ? 0 : 1;
  const AliMpVSegmentation *mpSeg = fSeg[cath]->GetDESegmentation(fidDE)->GetMpSegmentation();
  TVector2 dim = mpSeg->Dimensions(); 
  if (fChamber < 4) dim *= 2.; // quadrants
  //cout << fChamber << " " << dim.X() << " " << dim.Y() << endl;
  if (TMath::Abs (xl) - dim.X() > dx) return kFALSE;
  if (TMath::Abs (yl) - dim.Y() > dy) return kFALSE;

  Int_t flagX = 0, flagY = 0;
  if (TMath::Abs (xl) - dim.X() < -2) flagX = 1;
  if (TMath::Abs (yl) - dim.Y() < -2) flagY = 1;
  if (flagX && flagY) flagX = flagY = 0; // both coordinates "inside" det. elem. limits

  // Check for edge effect (extrapolated track "right outside" det. elem. boundaries 
  // (+- dx * cath in X and Y)
  dx /= 2;
  dy /= 2;
  for (Int_t i = 0; i < 2; i++) {
    if (!fSeg[i]) continue;
    for (Int_t idx = -1; idx < 2; idx++) {
      if (flagX && idx) continue;
      Double_t x1 = x + dx * idx * (i + 1);
      for (Int_t idy = -1; idy < 2; idy++) {
	if (idx == 0 && idy == 0) continue;
	if (flagY && idy) continue;
	Double_t y1 = y + dy * idy * (i + 1);
	fSeg[i]->GetPadI(fidDE, x1, y1, z, ix, iy);
	//cout << x1 << " " << y1 << " " << z << " " << fChamber << " " << ix << " " << iy << " " << fSeg[i]->Npx(fidDE) << " " << fSeg[i]->Npy(fidDE) /*<< " " << fSeg[i]->GetPadI(fidDE, x, y, z, ix, iy)*/ << " " << fSeg[i]->HasPad(fidDE, (Float_t)x1, (Float_t)y1, (Float_t)z) << endl; 
	//if (ix > 0 && iy > 0 && ix <= fSeg[i]->Npx(fidDE) && iy <= fSeg[i]->Npy(fidDE)) return kTRUE;
	if (fSeg[i]->HasPad(fidDE, (Float_t)x1, (Float_t)y1, (Float_t)z)) return kTRUE;
      }
    }
  }
  return kFALSE;
}

//_________________________________________________________________________
void AliMUONDetElement::ClusterReco(Double_t xTrack, Double_t yTrack, Double_t dx, Double_t dy)
{
  /// Run cluster reconstruction around point (x,y) inside window (dx,dy)

  if (fLeft[0] == 0 && fLeft[1] == 0) return; // all digits have been used 
  AliMUONClusterInput::Instance()->SetDigits(fChamber, fidDE, 
             (TClonesArray*)fDigits[0], (TClonesArray*)fDigits[1]);

  // Mark used pads
  for (Int_t cath = 0; cath < 2; cath++) {
    Int_t ndig = fDigits[cath]->GetEntriesFast();
    if (ndig == 0) continue; // empty cathode
    //cout << fidDE << " " << cath << " " << ndig << endl; 
    for (Int_t i = 0; i < ndig; i++) {
      if (fLeft[cath] == 0) { fRecModel->SetUsed(cath,i); continue; }
      AliMUONDigit *dig = (AliMUONDigit*) fDigits[cath]->UncheckedAt(i);
      //cout << i << " " << dig << " " /*<< dig->PadX() << " " << dig->PadY() << " " << fHitMap[cath]->TestHit(dig->PadX(), dig->PadY()) */<< endl;
      if (fHitMap[cath]->TestHit(dig->PadX(), dig->PadY()) == kUsed) fRecModel->SetUsed(cath,i);
      else fRecModel->SetUnused(cath,i);
    }
  }

  fRecModel->ResetRawClusters();

  for (Int_t cath = 0; cath < 2; cath++) {
    Int_t ndig = fDigits[cath]->GetEntriesFast();
    if (ndig == 0) continue; // empty cathode
    // Loop over pads
    for (fSeg[cath]->FirstPad(fidDE, xTrack, yTrack, fZ, dx, dy);
	 fSeg[cath]->MorePads(fidDE);
	 fSeg[cath]->NextPad(fidDE)) {
      if (fLeft[cath] == 0) break;
      //cout << cath << " " << fSeg[cath]->Ix() << " " << fSeg[cath]->Iy() << " " << fSeg[cath]->DetElemId() << " " << fHitMap[cath]->TestHit(fSeg[cath]->Ix(), fSeg[cath]->Iy()) << endl;
      if (fHitMap[cath]->TestHit(fSeg[cath]->Ix(), fSeg[cath]->Iy()) == kEmpty ||
	  fHitMap[cath]->TestHit(fSeg[cath]->Ix(), fSeg[cath]->Iy()) == kUsed) continue;

      Int_t iOK = 0;
      // Set starting pad
      for (Int_t j = 0; j < ndig; j++) {
	AliMUONDigit *dig = (AliMUONDigit*) fDigits[cath]->UncheckedAt(j);
	if (dig->PadX() != fSeg[cath]->Ix() || dig->PadY() != fSeg[cath]->Iy()) continue;
	//cout << fidDE << " " << j << " " << fDigits[cath]->GetEntriesFast() << " " << fSeg[cath]->Ix() << " " << fSeg[cath]->Iy() << endl;
	fRecModel->SetStart(cath, j);
	iOK = 1;
	break;
      }
      if (!iOK) continue;

      fRecModel->FindRawClusters();
      Int_t nClusEnd = fRecModel->GetRawClusters()->GetEntriesFast();
      //cout << " ***nclus: " << nClusEnd << endl;
      for (Int_t i = 0; i < nClusEnd; i++) {
	AliMUONRawCluster *clus = (AliMUONRawCluster*) fRecModel->GetRawClusters()->UncheckedAt(i);
	AddHitForRec(clus); // add hit for rec.
	//cout << clus->GetX(0) << " " << clus->GetY(0) << endl;
      }
      // Mark used pads
      for (Int_t cath1 = 0; cath1 < 2; cath1++) {
	Int_t ndig1 = fDigits[cath1]->GetEntriesFast();
	for (Int_t j = 0; j < ndig1; j++) {
	  if (fLeft[cath1] == 0) break;
	  AliMUONDigit *dig = (AliMUONDigit*) fDigits[cath1]->UncheckedAt(j);
	  Float_t x, y, z;
	  fSeg[cath1]->GetPadC(fidDE,dig->PadX(),dig->PadY(),x,y,z);
	  //cout << "clus " << cath1 << " " << fLeft[cath1] << " " << dig->PadX() << " " << dig->PadY() << " " << x << " " << y << " " << z << " " << fRecModel->GetUsed(cath1,j) << endl;
	  if (!fRecModel->GetUsed(cath1,j)) continue;
	  if (fHitMap[cath1]->TestHit(dig->PadX(), dig->PadY()) == kUsed) continue;
	  fHitMap[cath1]->FlagHit(dig->PadX(), dig->PadY());
	  if (lun) fprintf(lun," %d %d %d %d \n", cath1, fidDE, dig->PadX(), dig->PadY());
	  fLeft[cath1]--;
	}
      }
    } // for (fSeg[cath]->FirstPad(...
  } // for (Int_t cath = 0;
}

//_________________________________________________________________________
void AliMUONDetElement::AddHitForRec(AliMUONRawCluster *clus)
{
  /// Make HitForRec from raw cluster (rec. point)

  fRawClus->Add(new AliMUONRawCluster(*clus));
  //AliMUONHitForRec *hitForRec = 
  //new ((*fHitsForRec)[fNHitsForRec++]) AliMUONHitForRec(clus);
  AliMUONHitForRec *hitForRec = new AliMUONHitForRec(clus);
  fHitsForRec->Add(hitForRec);
  fNHitsForRec++;

  // more information into HitForRec
  //  resolution: info should be already in raw cluster and taken from it ????
  //hitForRec->SetBendingReso2(-1); //fBendingResolution * fBendingResolution);
  //hitForRec->SetNonBendingReso2(-1); //fNonBendingResolution * fNonBendingResolution);
  hitForRec->SetBendingReso2(clus->GetErrY() * clus->GetErrY());
  hitForRec->SetNonBendingReso2(clus->GetErrX() * clus->GetErrX());
  //  original raw cluster
  hitForRec->SetChamberNumber(fChamber);
  hitForRec->SetZ(clus->GetZ(0));
  //hitForRec->SetHitNumber(-(fIndex+1)*100000-fNHitsForRec+1);
  hitForRec->SetHitNumber(-(fIndex+1)*100000-fRawClus->GetEntriesFast()+1);
  //delete clus; // for now
}

/*
//_________________________________________________________________________
Int_t AliMUONDetElement::GetMapElem(AliMUONDigit *digit)
{
  Int_t cath = digit->Cathode();
  return 0;

}

//_________________________________________________________________________
void AliMUONDetElement::SetMapElem(const AliMUONDigit *digit, Int_t flag)
{
  Int_t cath = digit->Cathode();
}
*/
Commit	Line	Data
cc87ebcd	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	7	* Permission to use, copy, modify and distribute this software and its *
	8	* documentation strictly for non-commercial purposes is hereby granted *
	9	* without fee, provided that the above copyright notice appears in all *
	10	* copies and that both the copyright notice and this permission notice *
	11	* appear in the supporting documentation. The authors make no claims *
	12	* about the suitability of this software for any purpose. It is *
	13	* provided "as is" without express or implied warranty. *
	14	**************************************************************************/
	15
	16	/* $Id$ */
	17
d19b6003	18	// --------------------------
	19	// Class AliMUONDetElement
	20	// --------------------------
	21	// Detection element object containing information for combined
03ee5d34	22	// cluster / track finder in the MUON arm
d19b6003	23	// Author: Alexander Zinchenko, JINR Dubna
d19b6003	24
cc87ebcd	25	#include "AliMUONDetElement.h"
	26	#include "AliMUON.h"
	27	#include "AliMUONSegmentation.h"
	28	#include "AliMUONDigit.h"
	29	#include "AliMUONHitMapA1.h"
03ee5d34	30	#include "AliMUONData.h"
cc87ebcd	31	#include "AliMUONRawCluster.h"
	32	#include "AliMUONHitForRec.h"
	33	#include "AliMUONClusterInput.h"
	34	#include "AliMUONClusterFinderAZ.h"
03ee5d34	35	#include "AliMUONGeometryModuleTransformer.h" //AZ
	36	#include "AliMUONVGeometryDESegmentation.h" //AZ
	37	#include "AliMpVSegmentation.h" //AZ
	38
cc87ebcd	39	#include "AliRun.h"
	40	#include "AliLog.h"
	41
03ee5d34	42	#include <TObjArray.h>
	43	#include <TClonesArray.h>
	44	#include <TVector2.h> //AZ
	45
7945aae7	46	/// \cond CLASSIMP
cc87ebcd	47	ClassImp(AliMUONDetElement) // Class implementation in ROOT context
7945aae7	48	/// \endcond
03ee5d34	49	FILE *lun = 0x0; //fopen("hitmap.dat","w");
cc87ebcd	50
	51	//_________________________________________________________________________
	52	AliMUONDetElement::AliMUONDetElement()
54d7ba50	53	: TObject(),
	54	fidDE(0),
	55	fIndex(0),
	56	fChamber(0),
	57	fZ(0.),
	58	fNHitsForRec(0),
	59	fRawClus(0x0),
	60	fHitsForRec(0x0),
	61	fRecModel(0x0)
cc87ebcd	62	{
03ee5d34	63	/// Default constructor
cc87ebcd	64	for (Int_t i = 0; i < 2; i++) {
	65	fHitMap[i] = NULL;
	66	fDigits[i] = NULL;
	67	fSeg[i] = NULL;
	68	}
cc87ebcd	69	}
	70
	71	//_________________________________________________________________________
	72	AliMUONDetElement::AliMUONDetElement(Int_t idDE, AliMUONDigit dig, AliMUONClusterFinderAZ recModel)
54d7ba50	73	: TObject(),
	74	fidDE(idDE),
	75	fIndex(0),
	76	fChamber(idDE / 100 - 1),
	77	fZ(0.),
	78	fNHitsForRec(0),
	79	fRawClus(0x0),
03ee5d34	80	fHitsForRec(new TObjArray(10)),
54d7ba50	81	fRecModel(recModel)
cc87ebcd	82	{
03ee5d34	83	/// Constructor
cc87ebcd	84	fDigits[0] = new TObjArray(10);
	85	fDigits[1] = new TObjArray(10);
	86	fRawClus = new TObjArray(10);
cc87ebcd	87	AliMUON pMUON = (AliMUON) gAlice->GetModule("MUON");
cc87ebcd	88	AliMUONSegmentation *pSegmentation = pMUON->GetSegmentation();
03ee5d34	89	fSeg[0] = pSegmentation->GetModuleSegmentationByDEId(idDE, 0);
	90	fSeg[1] = pSegmentation->GetModuleSegmentationByDEId(idDE, 1);
	91	/*
cc87ebcd	92	Float_t x, y, z;
	93	fSeg[dig->Cathode()]->GetPadC(fidDE, dig->PadX(), dig->PadY(), x, y, z);
	94	fZ = z;
03ee5d34	95	*/
cc87ebcd	96	AddDigit(dig);
	97	}
	98
	99	//_________________________________________________________________________
	100	AliMUONDetElement::~AliMUONDetElement()
	101	{
03ee5d34	102	/// Destructor
03ee5d34	103
cc87ebcd	104	for (Int_t i = 0; i < 2; i++) {
	105	delete fHitMap[i]; fHitMap[i] = NULL;
	106	delete fDigits[i]; fDigits[i] = NULL;
03ee5d34	107	//fDigits[i]->Delete(); delete fDigits[i]; fDigits[i] = NULL;
cc87ebcd	108	}
03ee5d34	109	//fHitsForRec->Clear(); delete fHitsForRec; fHitsForRec = 0;
cc87ebcd	110	if (fRawClus) { fRawClus->Delete(); delete fRawClus; fRawClus = 0; }
cc87ebcd	111	//if (fRawClus) { delete fRawClus; fRawClus = 0; }
03ee5d34	112	for (Int_t i = 0; i < fNHitsForRec; i++) {
	113	Int_t ntracks = ((AliMUONHitForRec*)fHitsForRec->UncheckedAt(i))->GetNTrackHits();
	114	//cout << ntracks << endl;
	115	if (ntracks) fHitsForRec->RemoveAt(i);
	116	}
	117	//fHitsForRec->Delete();
cc87ebcd	118	delete fHitsForRec; fHitsForRec = 0;
03ee5d34	119	//fHitsForRec->SetOwner(kFALSE);
03ee5d34	120	//fHitsForRec->Clear(); cout << " here " << fHitsForRec->GetEntriesFast() << endl; delete fHitsForRec; fHitsForRec = 0;
cc87ebcd	121	}
cc87ebcd	122
cc87ebcd	123	//_________________________________________________________________________
	124	Int_t AliMUONDetElement::Compare(const TObject* detElem) const
	125	{
03ee5d34	126	/// "Compare" function to sort in Z (towards interaction point)
	127	/// Returns -1 (0, +1) if charge of current pixel
	128	/// is greater than (equal to, less than) charge of pixel
cc87ebcd	129	if (fZ > ((AliMUONDetElement*)detElem)->Z()) return(+1);
	130	else if (fZ == ((AliMUONDetElement*)detElem)->Z()) return( 0);
	131	else return(-1);
	132	}
	133
	134	//_________________________________________________________________________
	135	void AliMUONDetElement::Fill(AliMUONData /data*/)
	136	{
03ee5d34	137	/// Fill hit maps
cc87ebcd	138	fLeft[0] = fDigits[0]->GetEntriesFast();
	139	fLeft[1] = fDigits[1]->GetEntriesFast();
	140
4e92ede3	141	Int_t npx0 = fSeg[0]->Npx(fidDE)+1;
	142	Int_t npy0 = fSeg[0]->Npy(fidDE)+1;
	143	fHitMap[0] = new AliMUONHitMapA1(npx0, npy0, fDigits[0]);
cc87ebcd	144	fHitMap[0]->FillHits();
4e92ede3	145
	146	Int_t npx1 = fSeg[1]->Npx(fidDE)+1;
	147	Int_t npy1 = fSeg[1]->Npy(fidDE)+1;
	148	fHitMap[1] = new AliMUONHitMapA1(npx1, npy1, fDigits[1]);
cc87ebcd	149	fHitMap[1]->FillHits();
	150
	151	// The part below is just for debugging (fill rec. points already found)
	152	/*
	153	fLeft[0] = fLeft[1] = 0;
	154	TClonesArray *rawClus = data->RawClusters(fChamber);
	155	cout << rawClus << " " << rawClus->GetEntriesFast() << endl;
	156	for (Int_t i = 0; i < rawClus->GetEntriesFast(); i++) {
	157	AliMUONRawCluster recP = (AliMUONRawCluster) rawClus->UncheckedAt(i);
	158	cout << fChamber << " " << recP->GetZ(0) << " " << recP->GetZ(1) << " " << fZ << endl;
	159	if (TMath::Abs(recP->GetZ(0)-fZ) > 0.5) continue;
	160	if (!Inside(recP->GetX(0), recP->GetY(0), recP->GetZ(0))) continue;
	161	AddHitForRec(recP); // add hit for rec.
	162	rawClus->RemoveAt(i); // remove
	163	}
	164	cout << fHitsForRec->GetEntriesFast() << endl;
	165	rawClus->Compress();
	166	*/
	167	}
	168
	169	//_________________________________________________________________________
	170	void AliMUONDetElement::AddDigit(AliMUONDigit *dig)
	171	{
03ee5d34	172	/// Add digit
cc87ebcd	173
03ee5d34	174	//fDigits[dig->Cathode()]->Add(dig);
	175	fDigits[dig->Cathode()]->Add(new AliMUONDigit(*dig));
	176	Float_t x, y, z;
	177	fSeg[dig->Cathode()]->GetPadC(fidDE, dig->PadX(), dig->PadY(), x, y, z);
	178	fZ += z;
	179	/*
	180	Int_t ndig = fDigits[dig->Cathode()]->GetEntriesFast();
	181	new((fDigits[dig->Cathode()])[ndig]) AliMUONDigit(dig);
	182	*/
cc87ebcd	183	}
	184
	185	//_________________________________________________________________________
03ee5d34	186	Bool_t AliMUONDetElement::Inside(Double_t x, Double_t y, Double_t z, Double_t dx, Double_t dy) const
cc87ebcd	187	{
03ee5d34	188	/// Check if point is inside detection element
cc87ebcd	189
03ee5d34	190	dx = TMath::Max (dx, 1.);
03ee5d34	191	dy = TMath::Max (dy, 1.);
cc87ebcd	192	Int_t ix, iy;
03ee5d34	193	ix = iy = 0;
03ee5d34	194	Double_t xl, yl, zl;
cc87ebcd	195	for (Int_t i = 0; i < 2; i++) {
cc87ebcd	196	if (!fSeg[i]) continue;
03ee5d34	197	if (i == 0 \|\| fSeg[0] == 0x0) {
	198	fSeg[i]->GetTransformer()->Global2Local(fidDE, x, y, z, xl, yl, zl);
	199	//cout << xl << " " << yl << " " << zl << endl;
	200	if (TMath::Abs (zl) > 0.5) return kFALSE;
	201	//if (fChamber < 4 && (xl < 0 \|\| yl < 0)) return kFALSE;
	202	}
cc87ebcd	203	fSeg[i]->GetPadI(fidDE, x, y, z, ix, iy);
03ee5d34	204	//cout << x << " " << y << " " << z << " " << fChamber << " " << ix << " " << iy << " " << fSeg[i]->Npx(fidDE) << " " << fSeg[i]->Npy(fidDE) /<< " " << fSeg[i]->GetPadI(fidDE, x, y, z, ix, iy)/ << " " << fSeg[i]->HasPad(fidDE, (Float_t)x, (Float_t)y, (Float_t)z) << endl;
	205	//if (ix > 0 && iy > 0 && ix <= fSeg[i]->Npx(fidDE) && iy <= fSeg[i]->Npy(fidDE)) return kTRUE;
	206	if (fSeg[i]->HasPad(fidDE, (Float_t)x, (Float_t)y, (Float_t)z)) return kTRUE;
cc87ebcd	207	}
03ee5d34	208	//Double_t xl, yl, zl;
	209	//fSeg[0]->GetTransformer()->Global2Local(fidDE, x, y, z, xl, yl, zl);
	210	//cout << xl << " " << yl << " " << zl << endl;
	211
	212	Int_t cath = fSeg[0] ? 0 : 1;
	213	const AliMpVSegmentation *mpSeg = fSeg[cath]->GetDESegmentation(fidDE)->GetMpSegmentation();
	214	TVector2 dim = mpSeg->Dimensions();
	215	if (fChamber < 4) dim *= 2.; // quadrants
	216	//cout << fChamber << " " << dim.X() << " " << dim.Y() << endl;
	217	if (TMath::Abs (xl) - dim.X() > dx) return kFALSE;
	218	if (TMath::Abs (yl) - dim.Y() > dy) return kFALSE;
	219
	220	Int_t flagX = 0, flagY = 0;
	221	if (TMath::Abs (xl) - dim.X() < -2) flagX = 1;
	222	if (TMath::Abs (yl) - dim.Y() < -2) flagY = 1;
	223	if (flagX && flagY) flagX = flagY = 0; // both coordinates "inside" det. elem. limits
	224
	225	// Check for edge effect (extrapolated track "right outside" det. elem. boundaries
	226	// (+- dx * cath in X and Y)
	227	dx /= 2;
	228	dy /= 2;
cc87ebcd	229	for (Int_t i = 0; i < 2; i++) {
	230	if (!fSeg[i]) continue;
	231	for (Int_t idx = -1; idx < 2; idx++) {
03ee5d34	232	if (flagX && idx) continue;
03ee5d34	233	Double_t x1 = x + dx * idx * (i + 1);
cc87ebcd	234	for (Int_t idy = -1; idy < 2; idy++) {
cc87ebcd	235	if (idx == 0 && idy == 0) continue;
03ee5d34	236	if (flagY && idy) continue;
03ee5d34	237	Double_t y1 = y + dy * idy * (i + 1);
cc87ebcd	238	fSeg[i]->GetPadI(fidDE, x1, y1, z, ix, iy);
03ee5d34	239	//cout << x1 << " " << y1 << " " << z << " " << fChamber << " " << ix << " " << iy << " " << fSeg[i]->Npx(fidDE) << " " << fSeg[i]->Npy(fidDE) /<< " " << fSeg[i]->GetPadI(fidDE, x, y, z, ix, iy)/ << " " << fSeg[i]->HasPad(fidDE, (Float_t)x1, (Float_t)y1, (Float_t)z) << endl;
	240	//if (ix > 0 && iy > 0 && ix <= fSeg[i]->Npx(fidDE) && iy <= fSeg[i]->Npy(fidDE)) return kTRUE;
	241	if (fSeg[i]->HasPad(fidDE, (Float_t)x1, (Float_t)y1, (Float_t)z)) return kTRUE;
cc87ebcd	242	}
	243	}
	244	}
	245	return kFALSE;
	246	}
	247
	248	//_________________________________________________________________________
03ee5d34	249	void AliMUONDetElement::ClusterReco(Double_t xTrack, Double_t yTrack, Double_t dx, Double_t dy)
cc87ebcd	250	{
03ee5d34	251	/// Run cluster reconstruction around point (x,y) inside window (dx,dy)
cc87ebcd	252
cc87ebcd	253	if (fLeft[0] == 0 && fLeft[1] == 0) return; // all digits have been used
cc87ebcd	254	AliMUONClusterInput::Instance()->SetDigits(fChamber, fidDE,
	255	(TClonesArray)fDigits[0], (TClonesArray)fDigits[1]);
	256
	257	// Mark used pads
	258	for (Int_t cath = 0; cath < 2; cath++) {
03ee5d34	259	Int_t ndig = fDigits[cath]->GetEntriesFast();
	260	if (ndig == 0) continue; // empty cathode
	261	//cout << fidDE << " " << cath << " " << ndig << endl;
	262	for (Int_t i = 0; i < ndig; i++) {
cc87ebcd	263	if (fLeft[cath] == 0) { fRecModel->SetUsed(cath,i); continue; }
cc87ebcd	264	AliMUONDigit dig = (AliMUONDigit) fDigits[cath]->UncheckedAt(i);
03ee5d34	265	//cout << i << " " << dig << " " /<< dig->PadX() << " " << dig->PadY() << " " << fHitMap[cath]->TestHit(dig->PadX(), dig->PadY()) /<< endl;
cc87ebcd	266	if (fHitMap[cath]->TestHit(dig->PadX(), dig->PadY()) == kUsed) fRecModel->SetUsed(cath,i);
	267	else fRecModel->SetUnused(cath,i);
	268	}
	269	}
	270
	271	fRecModel->ResetRawClusters();
	272
	273	for (Int_t cath = 0; cath < 2; cath++) {
03ee5d34	274	Int_t ndig = fDigits[cath]->GetEntriesFast();
03ee5d34	275	if (ndig == 0) continue; // empty cathode
cc87ebcd	276	// Loop over pads
	277	for (fSeg[cath]->FirstPad(fidDE, xTrack, yTrack, fZ, dx, dy);
	278	fSeg[cath]->MorePads(fidDE);
	279	fSeg[cath]->NextPad(fidDE)) {
	280	if (fLeft[cath] == 0) break;
	281	//cout << cath << " " << fSeg[cath]->Ix() << " " << fSeg[cath]->Iy() << " " << fSeg[cath]->DetElemId() << " " << fHitMap[cath]->TestHit(fSeg[cath]->Ix(), fSeg[cath]->Iy()) << endl;
	282	if (fHitMap[cath]->TestHit(fSeg[cath]->Ix(), fSeg[cath]->Iy()) == kEmpty \|\|
	283	fHitMap[cath]->TestHit(fSeg[cath]->Ix(), fSeg[cath]->Iy()) == kUsed) continue;
	284
03ee5d34	285	Int_t iOK = 0;
cc87ebcd	286	// Set starting pad
03ee5d34	287	for (Int_t j = 0; j < ndig; j++) {
cc87ebcd	288	AliMUONDigit dig = (AliMUONDigit) fDigits[cath]->UncheckedAt(j);
cc87ebcd	289	if (dig->PadX() != fSeg[cath]->Ix() \|\| dig->PadY() != fSeg[cath]->Iy()) continue;
03ee5d34	290	//cout << fidDE << " " << j << " " << fDigits[cath]->GetEntriesFast() << " " << fSeg[cath]->Ix() << " " << fSeg[cath]->Iy() << endl;
cc87ebcd	291	fRecModel->SetStart(cath, j);
03ee5d34	292	iOK = 1;
cc87ebcd	293	break;
cc87ebcd	294	}
03ee5d34	295	if (!iOK) continue;
cc87ebcd	296
	297	fRecModel->FindRawClusters();
	298	Int_t nClusEnd = fRecModel->GetRawClusters()->GetEntriesFast();
	299	//cout << " ***nclus: " << nClusEnd << endl;
	300	for (Int_t i = 0; i < nClusEnd; i++) {
	301	AliMUONRawCluster clus = (AliMUONRawCluster) fRecModel->GetRawClusters()->UncheckedAt(i);
	302	AddHitForRec(clus); // add hit for rec.
	303	//cout << clus->GetX(0) << " " << clus->GetY(0) << endl;
	304	}
	305	// Mark used pads
	306	for (Int_t cath1 = 0; cath1 < 2; cath1++) {
03ee5d34	307	Int_t ndig1 = fDigits[cath1]->GetEntriesFast();
03ee5d34	308	for (Int_t j = 0; j < ndig1; j++) {
cc87ebcd	309	if (fLeft[cath1] == 0) break;
	310	AliMUONDigit dig = (AliMUONDigit) fDigits[cath1]->UncheckedAt(j);
	311	Float_t x, y, z;
	312	fSeg[cath1]->GetPadC(fidDE,dig->PadX(),dig->PadY(),x,y,z);
	313	//cout << "clus " << cath1 << " " << fLeft[cath1] << " " << dig->PadX() << " " << dig->PadY() << " " << x << " " << y << " " << z << " " << fRecModel->GetUsed(cath1,j) << endl;
	314	if (!fRecModel->GetUsed(cath1,j)) continue;
	315	if (fHitMap[cath1]->TestHit(dig->PadX(), dig->PadY()) == kUsed) continue;
	316	fHitMap[cath1]->FlagHit(dig->PadX(), dig->PadY());
03ee5d34	317	if (lun) fprintf(lun," %d %d %d %d \n", cath1, fidDE, dig->PadX(), dig->PadY());
cc87ebcd	318	fLeft[cath1]--;
	319	}
	320	}
	321	} // for (fSeg[cath]->FirstPad(...
	322	} // for (Int_t cath = 0;
	323	}
	324
	325	//_________________________________________________________________________
	326	void AliMUONDetElement::AddHitForRec(AliMUONRawCluster *clus)
	327	{
03ee5d34	328	/// Make HitForRec from raw cluster (rec. point)
cc87ebcd	329
cc87ebcd	330	fRawClus->Add(new AliMUONRawCluster(*clus));
03ee5d34	331	//AliMUONHitForRec *hitForRec =
	332	//new ((*fHitsForRec)[fNHitsForRec++]) AliMUONHitForRec(clus);
	333	AliMUONHitForRec *hitForRec = new AliMUONHitForRec(clus);
	334	fHitsForRec->Add(hitForRec);
	335	fNHitsForRec++;
cc87ebcd	336
	337	// more information into HitForRec
	338	// resolution: info should be already in raw cluster and taken from it ????
03ee5d34	339	//hitForRec->SetBendingReso2(-1); //fBendingResolution * fBendingResolution);
	340	//hitForRec->SetNonBendingReso2(-1); //fNonBendingResolution * fNonBendingResolution);
	341	hitForRec->SetBendingReso2(clus->GetErrY() * clus->GetErrY());
	342	hitForRec->SetNonBendingReso2(clus->GetErrX() * clus->GetErrX());
cc87ebcd	343	// original raw cluster
	344	hitForRec->SetChamberNumber(fChamber);
	345	hitForRec->SetZ(clus->GetZ(0));
	346	//hitForRec->SetHitNumber(-(fIndex+1)*100000-fNHitsForRec+1);
	347	hitForRec->SetHitNumber(-(fIndex+1)*100000-fRawClus->GetEntriesFast()+1);
	348	//delete clus; // for now
	349	}
	350
	351	/*
	352	//_________________________________________________________________________
	353	Int_t AliMUONDetElement::GetMapElem(AliMUONDigit *digit)
	354	{
	355	Int_t cath = digit->Cathode();
	356	return 0;
	357
	358	}
	359
	360	//_________________________________________________________________________
	361	void AliMUONDetElement::SetMapElem(const AliMUONDigit *digit, Int_t flag)
	362	{
	363	Int_t cath = digit->Cathode();
	364	}
	365	*/