#include "TROOT.h"
#include "TStyle.h"
#include "TVector.h"
#include "TLinearFitter.h"
#include "TCanvas.h"
#include "TGeoMatrix.h"
#include "TGeoManager.h"
#include "TEveTrans.h"
#include "TEveManager.h"
#include "EveBase/AliEveEventManager.h"
#include "AliEveTRDData.h"
#include "AliEveTRDModuleImp.h"
#include "AliEveTRDLoader.h"
#include "AliEveTRDLoaderImp.h"
#include "AliGeomManager.h"
#include "AliESDtrack.h"
#include "AliLog.h"
#include "AliPID.h"
#include "AliTrackPointArray.h"
#include "AliRieman.h"
#include "AliTRDhit.h"
#include "AliTRDcluster.h"
#include "AliTRDseedV1.h"
#include "AliTRDtrackletMCM.h"
#include "AliTRDtrackletWord.h"
#include "AliTRDmcmSim.h"
#include "AliTRDtrackV1.h"
#include "AliTRDtrackerV1.h"
#include "AliTRDpadPlane.h"
#include "AliTRDdigitsManager.h"
#include "AliTRDmcmSim.h"
#include "AliTRDarrayADC.h"
#include "AliTRDSignalIndex.h"
#include "AliTRDgeometry.h"
#include "AliTRDtransform.h"
#include "AliTRDReconstructor.h"
#include "AliTRDrecoParam.h"
ClassImp(AliEveTRDHits)
ClassImp(AliEveTRDDigits)
ClassImp(AliEveTRDClusters)
ClassImp(AliEveTRDTracklet)
ClassImp(AliEveTRDTrack)
ClassImp(AliEveTRDTrackletOnline)
ClassImp(AliEveTRDmcm)
AliEveTRDDigits::AliEveTRDDigits(AliEveTRDChamber *p)
:TEveQuadSet("digits", "")
,fParent(p)
{
SetOwnIds(kTRUE);
gStyle->SetPalette(1, 0);
SetPalette(new TEveRGBAPalette(0, 512));
Reset(TEveQuadSet::kQT_RectangleYZ, kFALSE, 32);
}
AliEveTRDDigits::~AliEveTRDDigits()
{
}
void AliEveTRDDigits::SetData(AliTRDdigitsManager *digits)
{
Int_t det(fParent->GetID());
AliTRDarrayADC *data = digits->GetDigits(det);
if(!data->GetDim()) return;
data->Expand();
AliTRDSignalIndex *indexes = digits->GetIndexes(det);
if(!indexes->IsAllocated()) digits->BuildIndexes(det);
Double_t scale, dy, dz;
Int_t ly = AliTRDgeometry::GetLayer(det),
stk = AliTRDgeometry::GetStack(det),
sec = AliTRDgeometry::GetSector(det),
vid = AliGeomManager::LayerToVolUID(AliGeomManager::kTRD1 + ly, stk + AliTRDgeometry::Nstack() * sec);
SetNameTitle(Form("digits%03d", det), Form("D-%03d [%02d_%d_%d]", det, sec, stk, ly));
Short_t sig[7]={0,0,0,10,0,0,0};
AliTRDtransform transform(det);
AliTRDpadPlane *pp(fParent->fGeo->GetPadPlane(ly, stk));
Int_t row, col;
AliTRDcluster c;
indexes->ResetCounters();
while (indexes->NextRCIndex(row, col)){
dz = pp->GetRowSize(row);
dy = pp->GetColSize(col);
Short_t *const adc = data->GetDataAddress(row,col);
for (Int_t time(0); time<30; time++){
if(data->IsPadCorrupted(row, col, time)) break;
if(adc[time] <= 1) continue;
new (&c) AliTRDcluster(det, col, row, time, sig, vid);
transform.Transform(&c);
scale = adc[time] < 512 ? adc[time]/512. : 1.;
AddQuad(c.GetY()-0.5*dy, c.GetZ()-0.5*dz*scale, c.GetX(), dy*0.95, dz*scale);
QuadValue(Int_t(adc[time]));
QuadId(new TNamed(Form("ADC %d", adc[time]), Form("det[%3d(%02d_%d_%d)] col[%3d] row[%2d] tb[%2d]", det, sec, stk, ly, col, row, time)));
}
}
RefitPlex();
TEveTrans& t = RefMainTrans();
t.SetRotByAngles((sec+.5)*AliTRDgeometry::GetAlpha(), 0.,0.);
}
AliEveTRDHits::AliEveTRDHits() : TEvePointSet("hits", 20)
{
SetMarkerSize(.1);
SetMarkerColor(kGreen);
SetOwnIds(kTRUE);
}
AliEveTRDHits::~AliEveTRDHits()
{
}
void AliEveTRDHits::PointSelected(Int_t n)
{
AliTRDhit *h = NULL;
if(!(h = dynamic_cast<AliTRDhit*>(GetPointId(n)))) return;
printf("Id[%3d] Det[%3d] Reg[%c] TR[%c] Q[%3d] MC[%d] t[%f]\n",
n, h->GetDetector(),
h->FromAmplification() ? 'A' : 'D',
h->FromTRphoton() ? 'y' : 'n',
h->GetCharge(), h->GetTrack(), h->GetTime());
}
AliEveTRDClusters::AliEveTRDClusters():AliEveTRDHits()
{
SetName("clusters");
SetMarkerSize(.4);
SetMarkerStyle(24);
SetMarkerColor(kGray);
SetOwnIds(kTRUE);
}
void AliEveTRDClusters::PointSelected(Int_t n)
{
AliTRDcluster *c = dynamic_cast<AliTRDcluster*>(GetPointId(n));
if(!c) return;
c->Print();
Emit("PointSelected(Int_t)", n);
}
void AliEveTRDClusters::Print(Option_t *o) const
{
AliTRDcluster *c = NULL;
for(Int_t n = GetN(); n--;){
if(!(c = dynamic_cast<AliTRDcluster*>(GetPointId(n)))) continue;
c->Print(o);
}
}
void AliEveTRDClusters::Load(const Char_t *w) const
{
Int_t typ = -1;
if(strcmp(w, "hit")==0) typ = 0;
else if(strcmp(w, "dig")==0) typ = 1;
else if(strcmp(w, "cls")==0) typ = 2;
else if(strcmp(w, "tlt")==0) typ = 3;
else if(strcmp(w, "all")==0) typ = 4;
else{
AliInfo("The following arguments are accepted:");
AliInfo(" \"hit\" : loading of MC hits");
AliInfo(" \"dig\" : loading of digits");
AliInfo(" \"cls\" : loading of reconstructed clusters");
AliInfo(" \"tlt\" : loading of online tracklets");
AliInfo(" \"all\" : loading everything available");
return;
}
AliTRDcluster *c = NULL;
Int_t n = 0;
while((n = GetN() && !(c = dynamic_cast<AliTRDcluster*>(GetPointId(n))))) n++;
if(!c) return;
Int_t det = c->GetDetector();
AliEveTRDLoader *loader(NULL);
switch(typ){
case 0:
loader = new AliEveTRDLoader("Hits");
if(!loader->Open("TRD.Hits.root")){
delete loader;
return;
}
break;
case 1:
loader = new AliEveTRDLoader("Digits");
if(!loader->Open("TRD.Digits.root")){
delete loader;
return;
}
break;
case 2:
loader = new AliEveTRDLoader("Clusters");
if(!loader->Open("TRD.RecPoints.root")){
delete loader;
return;
}
break;
case 3:
loader = new AliEveTRDLoader("Tracklets");
if(!loader->Open("TRD.Tracklets.root")){
delete loader;
return;
}
break;
case 4:
loader = new AliEveTRDLoaderSim("MC");
if(!loader->Open("galice.root")){
delete loader; loader=NULL;
Int_t failed(0);
loader = new AliEveTRDLoader("Data");
if(!loader->Open("TRD.Digits.root")) failed++;
if(!loader->Open("TRD.RecPoints.root")) failed++;
if(!loader->Open("TRD.Tracklets.root")) failed++;
if(failed==3){delete loader;loader=NULL;}
}else loader->SetDataType(AliEveTRDLoader::kTRDHits | AliEveTRDLoader::kTRDDigits | AliEveTRDLoader::kTRDClusters);
break;
default: return;
}
if(!loader) return;
loader->AddChambers(AliTRDgeometry::GetSector(det),AliTRDgeometry::GetStack(det), AliTRDgeometry::GetLayer(det));
if(!loader->GoToEvent(AliEveEventManager::GetCurrent()->GetEventId())) AliWarning(Form("No data loaded for event %d", AliEveEventManager::GetCurrent()->GetEventId()));
gEve->AddElement(loader->GetChamber(det), *(BeginParents()));
gEve->Redraw3D();
}
AliEveTRDTracklet::AliEveTRDTracklet(AliTRDseedV1 *trklt):TEveLine()
,fClusters(NULL)
{
SetName("tracklet");
if(!gGeoManager) AliEveEventManager::AssertGeometry();
SetUserData(trklt);
Int_t sec = AliTRDgeometry::GetSector(trklt->GetDetector());
Float_t alpha((0.5+sec)*AliTRDgeometry::GetAlpha()),
cphi(TMath::Cos(alpha)),
sphi(TMath::Sin(alpha));
Float_t dx;
Float_t x0 = trklt->GetX0();
Float_t y0 = trklt->GetYref(0);
Float_t z0 = trklt->GetZref(0);
Float_t dydx = trklt->GetYref(1);
Float_t dzdx = trklt->GetZref(1);
Float_t tilt = trklt->GetTilt();
Float_t g[3];
AliTRDcluster *c = NULL;
for(Int_t ic=0; ic<AliTRDseedV1::kNclusters; ic++){
if(!(c = trklt->GetClusters(ic))) continue;
if(!fClusters) AddElement(fClusters = new AliEveTRDClusters());
dx = x0 - c->GetX();
Float_t zt = z0 - dx*dzdx;
Float_t yc = c->GetY();
c->SetY(yc-tilt*(c->GetZ()-zt));
Float_t loc[3] = { c->GetX(), c->GetY(), c->GetZ() };
g[0] = loc[0] * cphi - loc[1] * sphi;
g[1] = loc[0] * sphi + loc[1] * cphi;
g[2] = loc[2];
Int_t id = fClusters->SetNextPoint(g[0], g[1], g[2]);
c->SetY(yc);
fClusters->SetPointId(id, new AliTRDcluster(*c));
}
if(fClusters){
fClusters->SetName("TC clusters");
fClusters->SetTitle(Form("N[%d]", trklt->GetN()));
fClusters->SetMarkerColor(kMagenta);
}
SetTitle(Form("Det[%d] RC[%c] Layer[%d] P[%7.3f]", trklt->GetDetector(), trklt->IsRowCross()?'y':'n', trklt->GetPlane(), trklt->GetMomentum()));
SetLineColor(kRed);
y0 = trklt->GetYfit(0);
dydx = trklt->GetYfit(1);
Double_t xg(x0 * cphi - y0 * sphi),
yg(x0 * sphi + y0 * cphi);
SetPoint(0, xg, yg, z0);
dx = .5*AliTRDgeometry::CamHght()+AliTRDgeometry::CdrHght();
x0 -= dx;
y0 -= dydx*dx,
z0 -= dzdx*dx;
xg = x0 * cphi - y0 * sphi;
yg = x0 * sphi + y0 * cphi;
SetPoint(1, xg, yg, z0);
}
AliEveTRDTracklet::~AliEveTRDTracklet()
{
}
void AliEveTRDTracklet::Print(Option_t *o) const
{
AliTRDseedV1 *tracklet = (AliTRDseedV1*)GetUserData();
if(!tracklet) return;
tracklet->Print(o);
}
AliEveTRDTrack::AliEveTRDTrack(AliTRDtrackV1 *trk)
:TEveLine()
,fTrackState(0)
,fESDStatus(0)
,fAlpha(0.)
,fPoints(NULL)
,fRim(NULL)
{
SetUserData(trk);
SetName("");
fRim = new AliRieman(trk->GetNumberOfClusters());
AliTRDseedV1 *tracklet = NULL;
for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++){
if(!(tracklet = trk->GetTracklet(il))) continue;
if(!tracklet->IsOK()) continue;
AddElement(new AliEveTRDTracklet(tracklet));
}
if(trk->GetNumberOfTracklets()>1) fRim->Update();
SetStatus(fTrackState);
}
AliEveTRDTrack::~AliEveTRDTrack()
{
if(fPoints) delete [] fPoints; fPoints = NULL;
}
void AliEveTRDTrack::Print(Option_t *o) const
{
AliTRDtrackV1 *track = (AliTRDtrackV1*)GetUserData();
if(!track) return;
track->Print(o);
}
void AliEveTRDTrack::SetStatus(UChar_t s)
{
if(fPoints && fTrackState == s) return;
const Int_t nc = AliTRDtrackV1::kMAXCLUSTERSPERTRACK;
AliTRDtrackV1 *trk(NULL);
if(!(trk=static_cast<AliTRDtrackV1*>(GetUserData()))) {
AliError("Failed casting data to TRD track.");
return;
}
Bool_t BUILD = kFALSE;
if(!fPoints){
fPoints = new AliTrackPoint[nc];
Double_t xmin = -1., xmax = -1.;
Int_t det = 0;
AliTRDseedV1 *trklt = NULL;
for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
if(!(trklt = trk->GetTracklet(ily))) continue;
if(xmin<0.) xmin = trklt->GetX0() - AliTRDgeometry::CamHght() - AliTRDgeometry::CdrHght();
if(trklt->GetX0()>xmax) xmax = trklt->GetX0();
det = trklt->GetDetector();
}
Int_t sec = det/AliTRDgeometry::kNdets;
fAlpha = AliTRDgeometry::GetAlpha() * (sec<9 ? sec + .5 : sec - 17.5);
Double_t dx =(xmax - xmin)/nc;
for(Int_t ip=0; ip<nc; ip++){
fPoints[ip].SetXYZ(xmin, 0., 0.);
xmin+=dx;
}
BUILD = kTRUE;
}
if(BUILD || ((s&12) != (fTrackState&12))){
if(TESTBIT(s, kTrackCosmics)){
AliTRDtrackerV1::FitLine(trk, NULL, kFALSE, nc, fPoints);
} else {
if(TESTBIT(s, kTrackModel)){
if(trk->GetNumberOfTracklets() >=4) AliTRDtrackerV1::FitKalman(trk, NULL, kFALSE, nc, fPoints);
} else {
Float_t x = 0.;
for(Int_t ip = nc; ip--;){
x = fPoints[ip].GetX();
fPoints[ip].SetXYZ(x, fRim->GetYat(x), fRim->GetZat(x));
}
}
}
Float_t global[3];
for(Int_t ip=0; ip<nc; ip++){
fPoints[ip].Rotate(-fAlpha).GetXYZ(global);
SetPoint(ip, global[0], global[1], global[2]);
}
SetSmooth(kTRUE);
}
if(BUILD || ((s&3) != (fTrackState&3))){
if(TESTBIT(s, kSource)){
if(fESDStatus&AliESDtrack::kTRDin){
SetMarkerColor(kGreen);
SetLineColor(kGreen);
} else {
SetMarkerColor(kMagenta);
SetLineColor(kMagenta);
}
} else {
if(TESTBIT(s, kPID) == AliTRDpidUtil::kLQ){
} else {
}
Int_t species = 0; Float_t pid = 0.;
for(Int_t is=0; is<AliPID::kSPECIES; ++is)
if(trk->GetPID(is) > pid){
pid = trk->GetPID(is);
species = (AliPID::EParticleType) is;
}
switch(species){
case AliPID::kElectron:
SetMarkerColor(kRed);
SetLineColor(kRed);
break;
default:
SetMarkerColor(kBlue);
SetLineColor(kBlue);
break;
}
}
SetLineWidth(2);
}
const Char_t *model = "line";
if(!TESTBIT(s, kTrackCosmics)){
if(TESTBIT(s, kTrackModel)) model = "kalman";
else model = "rieman";
}
Int_t species = 0; Float_t pid = 0.;
for(Int_t is=0; is<AliPID::kSPECIES; is++)
if(trk->GetPID(is) > pid){
pid = trk->GetPID(is);
species = is;
}
SetTitle(Form(
"Tracklets[%d] Clusters[%d]\n"
"Reconstruction Source[%s]\n"
"PID[%4.1f %4.1f %4.1f %4.1f %4.1f]\n"
"MC[%d]", trk->GetNumberOfTracklets(), trk->GetNumberOfClusters(), fESDStatus&AliESDtrack::kTRDin ? "barrel" : "sa",
1.E2*trk->GetPID(0), 1.E2*trk->GetPID(1),
1.E2*trk->GetPID(2), 1.E2*trk->GetPID(3), 1.E2*trk->GetPID(4), trk->GetLabel()));
if(GetName()){
char id[6]; snprintf(id, 6, "%s", GetName());
SetName(Form("%s %s", id, AliPID::ParticleName(species)));
}
fTrackState = s;
}
void AliEveTRDTrack::Load(const Char_t *what) const
{
const AliEveTRDTracklet* trklt(NULL);
TEveElement::List_ci itrklt=BeginChildren();
while(itrklt!=EndChildren()){
if((trklt = dynamic_cast<const AliEveTRDTracklet*>(*itrklt))) trklt->Load(what);
itrklt++;
}
}
AliEveTRDTrackletOnline::AliEveTRDTrackletOnline(AliTRDtrackletMCM *tracklet) :
TEveLine(),
fDetector(-1),
fROB(-1),
fMCM(-1)
{
AliTRDtrackletMCM *trkl = new AliTRDtrackletMCM(*tracklet);
SetUserData(trkl);
fDetector = trkl->GetDetector();
fROB = trkl->GetROB();
fMCM = trkl->GetMCM();
SetName("sim. tracklet");
SetTitle(Form("Det: %i, ROB: %i, MCM: %i, Label: %i\n0x%08x",
trkl->GetDetector(), trkl->GetROB(), trkl->GetMCM(), trkl->GetLabel(),
trkl->GetTrackletWord()));
SetLineColor(kGreen);
SetLineWidth(3);
AliTRDgeometry geo;
TGeoHMatrix *matrix = geo.GetClusterMatrix(trkl->GetDetector());
fDetector = trkl->GetDetector();
fROB = trkl->GetROB();
fMCM = trkl->GetMCM();
Float_t length = 3.;
Double_t x[3];
Double_t p[3];
Double_t p2[3];
x[0] = AliTRDgeometry::AnodePos();
x[1] = trkl->GetY();
x[2] = trkl->GetLocalZ();
matrix->LocalToMaster(x, p);
geo.RotateBack(trkl->GetDetector(), p, p2);
SetPoint(0, p2[0], p2[1], p2[2]);
x[0] -= length;
x[1] -= length * trkl->GetdYdX();
matrix->LocalToMaster(x, p);
p[2] *= p[0] / (p[0] + length);
geo.RotateBack(trkl->GetDetector(), p, p2);
SetPoint(1, p2[0], p2[1], p2[2]);
}
AliEveTRDTrackletOnline::AliEveTRDTrackletOnline(AliTRDtrackletWord *tracklet) :
TEveLine(),
fDetector(-1),
fROB(-1),
fMCM(-1)
{
AliTRDtrackletWord *trkl = new AliTRDtrackletWord(*tracklet);
SetUserData(trkl);
fDetector = trkl->GetDetector();
fROB = trkl->GetROB();
fMCM = trkl->GetMCM();
SetName("raw tracklet");
SetTitle(Form("Det: %i, ROB: %i, MCM: %i, Label: %i\n0x%08x",
trkl->GetDetector(), fROB, fMCM, -1,
trkl->GetTrackletWord()));
SetLineColor(kRed);
SetLineWidth(3);
AliTRDgeometry geo;
TGeoHMatrix *matrix = geo.GetClusterMatrix(trkl->GetDetector());
Float_t length = 3.;
Double_t x[3];
Double_t p[3];
Double_t p2[3];
x[0] = AliTRDgeometry::AnodePos();
x[1] = trkl->GetY();
x[2] = trkl->GetLocalZ();
matrix->LocalToMaster(x, p);
geo.RotateBack(trkl->GetDetector(), p, p2);
SetPoint(0, p2[0], p2[1], p2[2]);
x[0] -= length;
x[1] -= length * trkl->GetdYdX();
matrix->LocalToMaster(x, p);
p[2] *= p[0] / (p[0] + length);
geo.RotateBack(trkl->GetDetector(), p, p2);
SetPoint(1, p2[0], p2[1], p2[2]);
}
AliEveTRDTrackletOnline::~AliEveTRDTrackletOnline()
{
delete ((AliTRDtrackletBase*) GetUserData());
SetUserData(0x0);
}
void AliEveTRDTrackletOnline::ShowMCM(Option_t *opt) const
{
if (fDetector < 0 || fROB < 0 || fMCM < 0)
return;
AliEveTRDmcm *evemcm = new AliEveTRDmcm();
evemcm->Init(fDetector, fROB, fMCM);
evemcm->LoadDigits();
evemcm->Draw(opt);
TEveElementList *mcmlist = NULL;
if (gEve->GetCurrentEvent())
mcmlist = (TEveElementList*) gEve->GetCurrentEvent()->FindChild("TRD MCMs");
if (!mcmlist) {
mcmlist = new TEveElementList("TRD MCMs");
gEve->AddElement(mcmlist);
}
gEve->AddElement(evemcm, mcmlist);
}
AliEveTRDmcm::AliEveTRDmcm() :
TEveElement(),
TNamed(),
fMCM(new AliTRDmcmSim())
{
SetName("MCM");
SetTitle("Unknown MCM");
}
AliEveTRDmcm::~AliEveTRDmcm()
{
delete fMCM;
}
Bool_t AliEveTRDmcm::Init(Int_t det, Int_t rob, Int_t mcm)
{
SetName(Form("MCM: Det. %i, ROB %i, MCM %i", det, rob, mcm));
SetTitle(Form("MCM: Det. %i, ROB %i, MCM %i", det, rob, mcm));
fMCM->Init(det, rob, mcm);
fMCM->Reset();
return kTRUE;
}
Bool_t AliEveTRDmcm::LoadDigits()
{
AliRunLoader *rl = AliEveEventManager::AssertRunLoader();
return fMCM->LoadMCM(rl, fMCM->GetDetector(),
fMCM->GetRobPos(), fMCM->GetMcmPos());
}
Bool_t AliEveTRDmcm::Filter()
{
fMCM->Filter();
return kTRUE;
}
Bool_t AliEveTRDmcm::Tracklet()
{
fMCM->Tracklet();
return kTRUE;
}
void AliEveTRDmcm::Draw(Option_t* option)
{
const char *mcmname = Form("mcm_%i_%i_%i", fMCM->GetDetector(),
fMCM->GetRobPos(), fMCM->GetMcmPos());
TCanvas *c = dynamic_cast<TCanvas*> (gROOT->FindObject(mcmname));
if (!c)
c = gEve->AddCanvasTab("TRD MCM");
c->SetTitle(Form("MCM %i on ROB %i of det. %i",
fMCM->GetMcmPos(), fMCM->GetRobPos(), fMCM->GetDetector()));
c->SetName(mcmname);
c->cd();
fMCM->Draw(option);
}
Bool_t AliEveTRDmcm::AssignPointer(const char* ptrname)
{
gROOT->ProcessLine(Form("AliTRDmcmSim* %s = (AliTRDmcmSim *)%p", ptrname, (void*)fMCM));
return kTRUE;
}