#include "TEveCaloData.h"
#include "TEveCalo.h"
#include "TGLSelectRecord.h"
#include "TAxis.h"
#include "THStack.h"
#include "TH2.h"
#include "TMath.h"
#include "TList.h"
#include <cassert>
#include <algorithm>
#include <set>
void TEveCaloData::CellGeom_t::Dump() const
{
printf("%f, %f %f, %f \n", fEtaMin, fEtaMax, fPhiMin, fPhiMax);
}
void TEveCaloData::CellGeom_t::Configure(Float_t etaMin, Float_t etaMax, Float_t phiMin, Float_t phiMax)
{
fEtaMin = etaMin;
fEtaMax = etaMax;
fPhiMin = phiMin;
fPhiMax = phiMax;
if (fPhiMin < - TMath::TwoPi() || fPhiMin > TMath::TwoPi() ||
fPhiMax < - TMath::TwoPi() || fPhiMax > TMath::TwoPi())
{
::Error("TEveCaloData::CellGeom_t::Configure", "phiMin and phiMax should be between -2*pi and 2*pi (min=%f, max=%f). RhoZ projection will be wrong.",
fPhiMin, fPhiMax);
}
fThetaMin = EtaToTheta(fEtaMax);
fThetaMax = EtaToTheta(fEtaMin);
}
Float_t TEveCaloData::CellData_t::Value(Bool_t isEt) const
{
if (isEt)
return fValue;
else
return TMath::Abs(fValue/TMath::Sin(Theta()));
}
void TEveCaloData::CellData_t::Dump() const
{
printf("%f, %f %f, %f \n", fEtaMin, fEtaMax, fPhiMin, fPhiMax);
}
Float_t* TEveCaloData::RebinData_t::GetSliceVals(Int_t bin)
{
if (fBinData[bin] == -1)
{
fBinData[bin] = fSliceData.size();
for (Int_t i=0; i<fNSlices; i++)
fSliceData.push_back(0.f);
}
return &fSliceData[fBinData[bin]];
}
ClassImp(TEveCaloData);
TEveCaloData::TEveCaloData(const char* n, const char* t):
TEveElement(),
TNamed(n, t),
fEtaAxis(0),
fPhiAxis(0),
fWrapTwoPi(kTRUE),
fMaxValEt(0),
fMaxValE(0),
fEps(0)
{
}
void TEveCaloData::UnSelected()
{
fCellsSelected.clear();
}
void TEveCaloData::UnHighlighted()
{
fCellsHighlighted.clear();
}
TString TEveCaloData::GetHighlightTooltip()
{
if (fCellsHighlighted.empty()) return "";
CellData_t cellData;
Bool_t single = fCellsHighlighted.size() == 1;
Float_t sum = 0;
TString s;
for (vCellId_i i = fCellsHighlighted.begin(); i!=fCellsHighlighted.end(); ++i)
{
GetCellData(*i, cellData);
s += TString::Format("%s %.2f (%.3f, %.3f)",
fSliceInfos[i->fSlice].fName.Data(), cellData.fValue,
cellData.Eta(), cellData.Phi());
if (single) return s;
s += "\n";
sum += cellData.fValue;
}
s += TString::Format("Sum = %.2f", sum);
return s;
}
void TEveCaloData::FillImpliedSelectedSet(Set_t& impSelSet)
{
for (List_ci i=fChildren.begin(); i!=fChildren.end(); ++i)
{
impSelSet.insert(*i);
}
}
void TEveCaloData::PrintCellsSelected()
{
printf("%d Selected selected cells:\n", (Int_t)fCellsSelected.size());
CellData_t cellData;
for (vCellId_i i = fCellsSelected.begin(); i != fCellsSelected.end(); ++i)
{
GetCellData(*i, cellData);
printf("Tower [%d] Slice [%d] Value [%.2f] ", i->fTower, i->fSlice, cellData.fValue);
printf("Eta:(%f, %f) Phi(%f, %f)\n", cellData.fEtaMin, cellData.fEtaMax, cellData.fPhiMin, cellData.fPhiMax);
}
}
void TEveCaloData::ProcessSelection(vCellId_t& sel_cells, TGLSelectRecord& rec)
{
typedef std::set<CellId_t> sCellId_t;
typedef std::set<CellId_t>::iterator sCellId_i;
struct helper
{
static void fill_cell_set(sCellId_t& cset, vCellId_t& cvec)
{
for (vCellId_i i = cvec.begin(); i != cvec.end(); ++i)
cset.insert(*i);
}
static void fill_cell_vec(vCellId_t& cvec, sCellId_t& cset)
{
for (sCellId_i i = cset.begin(); i != cset.end(); ++i)
cvec.push_back(*i);
}
};
vCellId_t& cells = rec.GetHighlight() ? fCellsHighlighted : fCellsSelected;
if (cells.empty())
{
if (!sel_cells.empty())
{
cells.swap(sel_cells);
rec.SetSecSelResult(TGLSelectRecord::kEnteringSelection);
}
}
else
{
if (!sel_cells.empty())
{
if (rec.GetMultiple())
{
sCellId_t cs;
helper::fill_cell_set(cs, cells);
for (vCellId_i i = sel_cells.begin(); i != sel_cells.end(); ++i)
{
std::set<CellId_t>::iterator csi = cs.find(*i);
if (csi == cs.end())
cs.insert(*i);
else
cs.erase(csi);
}
cells.clear();
if (cs.empty())
{
rec.SetSecSelResult(TGLSelectRecord::kLeavingSelection);
}
else
{
helper::fill_cell_vec(cells, cs);
rec.SetSecSelResult(TGLSelectRecord::kModifyingInternalSelection);
}
}
else
{
Bool_t differ = kFALSE;
if (cells.size() == sel_cells.size())
{
sCellId_t cs;
helper::fill_cell_set(cs, cells);
for (vCellId_i i = sel_cells.begin(); i != sel_cells.end(); ++i)
{
if (cs.find(*i) == cs.end())
{
differ = kTRUE;
break;
}
}
}
else
{
differ = kTRUE;
}
if (differ)
{
cells.swap(sel_cells);
rec.SetSecSelResult(TGLSelectRecord::kModifyingInternalSelection);
}
}
}
else
{
if (!rec.GetMultiple())
{
cells.clear();
rec.SetSecSelResult(TGLSelectRecord::kLeavingSelection);
}
}
}
if (rec.GetSecSelResult() != TGLSelectRecord::kNone)
{
CellSelectionChanged();
}
}
void TEveCaloData::SetSliceThreshold(Int_t slice, Float_t val)
{
fSliceInfos[slice].fThreshold = val;
InvalidateUsersCellIdCache();
}
Float_t TEveCaloData::GetSliceThreshold(Int_t slice) const
{
return fSliceInfos[slice].fThreshold;
}
void TEveCaloData::SetSliceColor(Int_t slice, Color_t col)
{
fSliceInfos[slice].fColor = col;
for (List_ci i=fChildren.begin(); i!=fChildren.end(); ++i)
{
(*i)->AddStamp(TEveElement::kCBObjProps);
}
}
Color_t TEveCaloData::GetSliceColor(Int_t slice) const
{
return fSliceInfos[slice].fColor;
}
void TEveCaloData::SetSliceTransparency(Int_t slice, Char_t t)
{
fSliceInfos[slice].fTransparency = t;
for (List_ci i=fChildren.begin(); i!=fChildren.end(); ++i)
{
(*i)->AddStamp(TEveElement::kCBObjProps);
}
}
Char_t TEveCaloData::GetSliceTransparency(Int_t slice) const
{
return fSliceInfos[slice].fTransparency;
}
void TEveCaloData::InvalidateUsersCellIdCache()
{
TEveCaloViz* calo;
for (List_ci i=fChildren.begin(); i!=fChildren.end(); ++i)
{
calo = dynamic_cast<TEveCaloViz*>(*i);
calo->InvalidateCellIdCache();
calo->StampObjProps();
}
}
void TEveCaloData::DataChanged()
{
TEveCaloViz* calo;
for (List_ci i=fChildren.begin(); i!=fChildren.end(); ++i)
{
calo = dynamic_cast<TEveCaloViz*>(*i);
calo->DataChanged();
calo->StampObjProps();
}
}
void TEveCaloData::CellSelectionChanged()
{
TEveCaloViz* calo;
for (List_ci i=fChildren.begin(); i!=fChildren.end(); ++i)
{
calo = dynamic_cast<TEveCaloViz*>(*i);
calo->CellSelectionChanged();
calo->StampColorSelection();
}
}
Float_t TEveCaloData::EtaToTheta(Float_t eta)
{
using namespace TMath;
if (eta < 0)
return Pi() - 2*ATan(Exp(- Abs(eta)));
else
return 2*ATan(Exp(- Abs(eta)));
}
ClassImp(TEveCaloDataVec);
TEveCaloDataVec::TEveCaloDataVec(Int_t nslices):
TEveCaloData(),
fTower(0),
fEtaMin( 1e3),
fEtaMax(-1e3),
fPhiMin( 1e3),
fPhiMax(-1e3)
{
fSliceInfos.assign(nslices, SliceInfo_t());
fSliceVec.assign(nslices, std::vector<Float_t> ());
}
TEveCaloDataVec::~TEveCaloDataVec()
{
if (fEtaAxis) delete fEtaAxis;
if (fPhiAxis) delete fPhiAxis;
}
Int_t TEveCaloDataVec::AddSlice()
{
fSliceInfos.push_back(SliceInfo_t());
fSliceVec.push_back(std::vector<Float_t> ());
fSliceVec.back().resize(fGeomVec.size(), 0.f);
return fSliceInfos.size() - 1;
}
Int_t TEveCaloDataVec::AddTower(Float_t etaMin, Float_t etaMax, Float_t phiMin, Float_t phiMax)
{
assert (etaMin < etaMax);
assert (phiMin < phiMax);
fGeomVec.push_back(CellGeom_t(etaMin, etaMax, phiMin, phiMax));
for (vvFloat_i it=fSliceVec.begin(); it!=fSliceVec.end(); ++it)
(*it).push_back(0);
if (etaMin < fEtaMin) fEtaMin = etaMin;
if (etaMax > fEtaMax) fEtaMax = etaMax;
if (phiMin < fPhiMin) fPhiMin = phiMin;
if (phiMax > fPhiMax) fPhiMax = phiMax;
fTower = fGeomVec.size() - 1;
return fTower;
}
void TEveCaloDataVec::FillSlice(Int_t slice, Float_t val)
{
fSliceVec[slice][fTower] = val;
}
void TEveCaloDataVec::FillSlice(Int_t slice, Int_t tower, Float_t val)
{
fSliceVec[slice][tower] = val;
}
void TEveCaloDataVec::GetCellList(Float_t eta, Float_t etaD,
Float_t phi, Float_t phiD,
TEveCaloData::vCellId_t &out) const
{
using namespace TMath;
Float_t etaMin = eta - etaD*0.5;
Float_t etaMax = eta + etaD*0.5;
Float_t phiMin = phi - phiD*0.5;
Float_t phiMax = phi + phiD*0.5;
Int_t nS = fSliceVec.size();
Int_t tower = 0;
Float_t fracx=0, fracy=0, frac;
Float_t minQ, maxQ;
for(vCellGeom_ci i=fGeomVec.begin(); i!=fGeomVec.end(); i++)
{
const CellGeom_t &cg = *i;
fracx = TEveUtil::GetFraction(etaMin, etaMax, cg.fEtaMin, cg.fEtaMax);
if (fracx > 1e-3)
{
minQ = cg.fPhiMin;
maxQ = cg.fPhiMax;
if (fWrapTwoPi)
{
if (maxQ < phiMin)
{
minQ += TwoPi(); maxQ += TwoPi();
}
else if (minQ > phiMax)
{
minQ -= TwoPi(); maxQ -= TwoPi();
}
}
if (maxQ >= phiMin && minQ <= phiMax)
{
fracy = TEveUtil::GetFraction(phiMin, phiMax, minQ, maxQ);
if (fracy > 1e-3)
{
frac = fracx*fracy;
for (Int_t s=0; s<nS; s++)
{
if (fSliceVec[s][tower] > fSliceInfos[s].fThreshold)
out.push_back(CellId_t(tower, s, frac));
}
}
}
}
tower++;
}
}
void TEveCaloDataVec::Rebin(TAxis* ax, TAxis* ay, vCellId_t &ids, Bool_t et, RebinData_t& rdata) const
{
rdata.fNSlices = GetNSlices();
rdata.fBinData.assign((ax->GetNbins()+2)*(ay->GetNbins()+2), -1);
CellData_t cd;
for (vCellId_i it = ids.begin(); it != ids.end(); ++it)
{
GetCellData(*it, cd);
Int_t iMin = ax->FindBin(cd.EtaMin());
Int_t iMax = ax->FindBin(cd.EtaMax());
Int_t jMin = ay->FindBin(cd.PhiMin());
Int_t jMax = ay->FindBin(cd.PhiMax());
for (Int_t i = iMin; i <= iMax; ++i)
{
if (i < 0 || i > ax->GetNbins()) continue;
for (Int_t j = jMin; j <= jMax; ++j)
{
if (j < 0 || j > ay->GetNbins()) continue;
Double_t ratio = TEveUtil::GetFraction(ax->GetBinLowEdge(i), ax->GetBinUpEdge(i), cd.EtaMin(), cd.EtaMax())
* TEveUtil::GetFraction(ay->GetBinLowEdge(j), ay->GetBinUpEdge(j), cd.PhiMin(), cd.PhiMax());
if (ratio > 1e-6f)
{
Float_t* slices = rdata.GetSliceVals(i + j*(ax->GetNbins()+2));
slices[(*it).fSlice] += ratio * cd.Value(et);
}
}
}
}
}
void TEveCaloDataVec::GetCellData(const TEveCaloData::CellId_t &id,
TEveCaloData::CellData_t& cellData) const
{
cellData.CellGeom_t::operator=( fGeomVec[id.fTower] );
cellData.fValue = fSliceVec[id.fSlice][id.fTower];
}
void TEveCaloDataVec::DataChanged()
{
using namespace TMath;
fMaxValE = 0;
fMaxValEt = 0;
Float_t sum=0;
for (UInt_t tw=0; tw<fGeomVec.size(); tw++)
{
sum=0;
for (vvFloat_i it=fSliceVec.begin(); it!=fSliceVec.end(); ++it)
sum += (*it)[tw];
if (sum > fMaxValEt ) fMaxValEt=sum;
sum /= Abs(Sin(EtaToTheta(fGeomVec[tw].Eta())));
if (sum > fMaxValE) fMaxValE=sum;
}
TEveCaloData::DataChanged();
}
void TEveCaloDataVec::SetAxisFromBins(Double_t epsX, Double_t epsY)
{
std::vector<Double_t> binX;
std::vector<Double_t> binY;
for(vCellGeom_ci i=fGeomVec.begin(); i!=fGeomVec.end(); i++)
{
const CellGeom_t &ch = *i;
binX.push_back(ch.EtaMin());
binX.push_back(ch.EtaMax());
binY.push_back(ch.PhiMin());
binY.push_back(ch.PhiMax());
}
std::sort(binX.begin(), binX.end());
std::sort(binY.begin(), binY.end());
Int_t cnt = 0;
Double_t sum = 0;
Double_t val;
Double_t dx = binX.back() - binX.front();
epsX *= dx;
std::vector<Double_t> newX;
newX.push_back(binX.front());
Int_t nX = binX.size()-1;
for(Int_t i=0; i<nX; i++)
{
val = (sum +binX[i])/(cnt+1);
if (binX[i+1] -val > epsX)
{
newX.push_back(val);
cnt = 0;
sum = 0;
}
else
{
sum += binX[i];
cnt++;
}
}
newX.push_back(binX.back());
cnt = 0;
sum = 0;
std::vector<Double_t> newY;
Double_t dy = binY.back() - binY.front();
epsY *= dy;
newY.push_back(binY.front());
Int_t nY = binY.size()-1;
for(Int_t i=0 ; i<nY; i++)
{
val = (sum +binY[i])/(cnt+1);
if (binY[i+1] -val > epsY )
{
newY.push_back(val);
cnt = 0;
sum = 0;
}
else
{
sum += binY[i];
cnt++;
}
}
newY.push_back(binY.back());
if (fEtaAxis) delete fEtaAxis;
if (fPhiAxis) delete fPhiAxis;
fEtaAxis = new TAxis(newX.size()-1, &newX[0]);
fPhiAxis = new TAxis(newY.size()-1, &newY[0]);
fEtaAxis->SetNdivisions(510);
fPhiAxis->SetNdivisions(510);
}
ClassImp(TEveCaloDataHist);
TEveCaloDataHist::TEveCaloDataHist():
TEveCaloData(),
fHStack(0)
{
fHStack = new THStack();
fEps = 1e-5;
}
TEveCaloDataHist::~TEveCaloDataHist()
{
delete fHStack;
}
void TEveCaloDataHist::DataChanged()
{
using namespace TMath;
fMaxValE = 0;
fMaxValEt = 0;
if (GetNSlices() < 1) return;
TH2* hist = GetHist(0);
fEtaAxis = hist->GetXaxis();
fPhiAxis = hist->GetYaxis();
for (Int_t ieta = 1; ieta <= fEtaAxis->GetNbins(); ++ieta)
{
Double_t eta = fEtaAxis->GetBinCenter(ieta);
for (Int_t iphi = 1; iphi <= fPhiAxis->GetNbins(); ++iphi)
{
Double_t value = 0;
for (Int_t i = 0; i < GetNSlices(); ++i)
{
hist = GetHist(i);
Int_t bin = hist->GetBin(ieta, iphi);
value += hist->GetBinContent(bin);
}
if (value > fMaxValEt ) fMaxValEt = value;
value /= Abs(Sin(EtaToTheta(eta)));
if (value > fMaxValE) fMaxValE = value;
}
}
TEveCaloData::DataChanged();
}
void TEveCaloDataHist::GetCellList(Float_t eta, Float_t etaD,
Float_t phi, Float_t phiD,
TEveCaloData::vCellId_t &out) const
{
using namespace TMath;
Float_t etaMin = eta - etaD*0.5 -fEps;
Float_t etaMax = eta + etaD*0.5 +fEps;
Float_t phiMin = phi - phiD*0.5 -fEps;
Float_t phiMax = phi + phiD*0.5 +fEps;
Int_t nEta = fEtaAxis->GetNbins();
Int_t nPhi = fPhiAxis->GetNbins();
Int_t nSlices = GetNSlices();
Int_t bin = 0;
Bool_t accept;
for (Int_t ieta = 1; ieta <= nEta; ++ieta)
{
if (fEtaAxis->GetBinLowEdge(ieta) >= etaMin && fEtaAxis->GetBinUpEdge(ieta) <= etaMax)
{
for (Int_t iphi = 1; iphi <= nPhi; ++iphi)
{
if (fWrapTwoPi )
{
accept = TEveUtil::IsU1IntervalContainedByMinMax
(phiMin, phiMax, fPhiAxis->GetBinLowEdge(iphi), fPhiAxis->GetBinUpEdge(iphi));
}
else
{
accept = fPhiAxis->GetBinLowEdge(iphi) >= phiMin && fPhiAxis->GetBinUpEdge(iphi) <= phiMax &&
fPhiAxis->GetBinLowEdge(iphi) >= phiMin && fPhiAxis->GetBinUpEdge(iphi) <= phiMax;
}
if (accept)
{
for (Int_t s = 0; s < nSlices; ++s)
{
TH2F *hist = GetHist(s);
bin = hist->GetBin(ieta, iphi);
if (hist->GetBinContent(bin) > fSliceInfos[s].fThreshold)
out.push_back(TEveCaloData::CellId_t(bin, s));
}
}
}
}
}
}
void TEveCaloDataHist::Rebin(TAxis* ax, TAxis* ay, TEveCaloData::vCellId_t &ids, Bool_t et, RebinData_t &rdata) const
{
rdata.fNSlices = GetNSlices();
rdata.fBinData.assign((ax->GetNbins()+2)*(ay->GetNbins()+2), -1);
TEveCaloData::CellData_t cd;
Float_t *val;
Int_t i, j, w;
Int_t binx, biny;
Int_t bin;
for (vCellId_i it=ids.begin(); it!=ids.end(); ++it)
{
GetCellData(*it, cd);
GetHist(it->fSlice)->GetBinXYZ((*it).fTower, i, j, w);
binx = ax->FindBin(fEtaAxis->GetBinCenter(i));
biny = ay->FindBin(fPhiAxis->GetBinCenter(j));
bin = biny*(ax->GetNbins()+2)+binx;
val = rdata.GetSliceVals(bin);
Double_t ratio = TEveUtil::GetFraction(ax->GetBinLowEdge(binx), ax->GetBinUpEdge(binx), cd.EtaMin(), cd.EtaMax())
* TEveUtil::GetFraction(ay->GetBinLowEdge(biny), ay->GetBinUpEdge(biny), cd.PhiMin(), cd.PhiMax());
val[(*it).fSlice] += cd.Value(et)*ratio;
}
}
void TEveCaloDataHist::GetCellData(const TEveCaloData::CellId_t &id,
TEveCaloData::CellData_t& cellData) const
{
TH2F* hist = GetHist(id.fSlice);
Int_t x, y, z;
hist->GetBinXYZ(id.fTower, x, y, z);
cellData.fValue = hist->GetBinContent(id.fTower);
cellData.Configure(hist->GetXaxis()->GetBinLowEdge(x),
hist->GetXaxis()->GetBinUpEdge(x),
hist->GetYaxis()->GetBinLowEdge(y),
hist->GetYaxis()->GetBinUpEdge(y));
}
Int_t TEveCaloDataHist::AddHistogram(TH2F* hist)
{
fHStack->Add(hist);
fSliceInfos.push_back(SliceInfo_t());
fSliceInfos.back().fName = hist->GetName();
fSliceInfos.back().fColor = hist->GetLineColor();
DataChanged();
return fSliceInfos.size() - 1;
}
TH2F* TEveCaloDataHist::GetHist(Int_t slice) const
{
assert(slice >= 0 && slice < fHStack->GetHists()->GetSize());
return (TH2F*) fHStack->GetHists()->At(slice);
}
void TEveCaloDataHist::GetEtaLimits(Double_t &min, Double_t &max) const
{
min = fEtaAxis->GetXmin();
max = fEtaAxis->GetXmax();
}
void TEveCaloDataHist::GetPhiLimits(Double_t &min, Double_t &max) const
{
min = fPhiAxis->GetXmin();
max = fPhiAxis->GetXmax();
}