diff --git a/examples/compclustepr.C b/examples/compclustepr.C new file mode 100644 index 0000000000000000000000000000000000000000..97bc64a8c0a3eb80930398aaa5d8904b206b031e --- /dev/null +++ b/examples/compclustepr.C @@ -0,0 +1,57 @@ +void compclustepr(Int_t run) +{ + TFile* f = new TFile(Form("hodtest_%d.root",run)); + cout << "hcana root file " << Form("hodtest_%d.root",run) << endl; + TH1F* h = eprmax; + + TFile* f1 = new TFile(Form("%d_hbk.root",run)); + cout << "Engine root file " << Form("%d_hbk.root",run) << endl; + TH1F* h1; + switch (run) { + case 50017 : + // h1 = h212; //A+ + break; + default : + h1 = h415; //edep + } + + TCanvas *c1 = new TCanvas("c1", "Shower Largest cluster Eprsh", 1000, 667); + + gPad->SetLogy(); + + h1->SetFillColor(kGreen); + h1->SetLineColor(kGreen); + h1->Draw(); + + h->SetLineColor(kBlue); + h->SetFillStyle(0); + h->SetLineWidth(2); + h->Draw("same"); + + TLatex l; + l.SetTextSize(0.04); + Float_t maxy = h1->GetBinContent(h1->GetMaximumBin()); + Float_t xmin = h1->GetXaxis()->GetXmin(); + Float_t xmax = h1->GetXaxis()->GetXmax(); + Float_t xt = xmin + 0.67*(xmax-xmin); + + l.SetTextColor(kGreen); + l.DrawLatex(xt,0.095*maxy,"Engine"); + l.SetTextColor(kBlue); + l.DrawLatex(xt,0.045*maxy,"hcana"); + + // Difference between the histograms. + + TCanvas *c2 = new TCanvas("c2", "Eprsh differences", 1000, 667); + + TH1F* dif = h->Clone(); + + dif->Add(h,h1,1.,-1.); + + dif->SetTitle("Eprsh Difference"); + dif->SetFillColor(kRed); + dif->SetLineColor(kRed); + dif->SetLineWidth(1); + dif->SetFillStyle(1111); + dif->Draw(); +} diff --git a/examples/compclustx.C b/examples/compclustx.C new file mode 100644 index 0000000000000000000000000000000000000000..1fe40d6aaf11421876cbb3b812ac10afe26638ad --- /dev/null +++ b/examples/compclustx.C @@ -0,0 +1,57 @@ +void compclustx(Int_t run) +{ + TFile* f = new TFile(Form("hodtest_%d.root",run)); + cout << "hcana root file " << Form("hodtest_%d.root",run) << endl; + TH1F* h = xmax; + + TFile* f1 = new TFile(Form("%d_hbk.root",run)); + cout << "Engine root file " << Form("%d_hbk.root",run) << endl; + TH1F* h1; + switch (run) { + case 50017 : + // h1 = h212; //A+ + break; + default : + h1 = h416; //X + } + + TCanvas *c1 = new TCanvas("c1", "Shower Largest cluster X", 1000, 667); + + gPad->SetLogy(); + + h1->SetFillColor(kGreen); + h1->SetLineColor(kGreen); + h1->Draw(); + + h->SetLineColor(kBlue); + h->SetFillStyle(0); + h->SetLineWidth(2); + h->Draw("same"); + + TLatex l; + l.SetTextSize(0.04); + Float_t maxy = h1->GetBinContent(h1->GetMaximumBin()); + Float_t xmin = h1->GetXaxis()->GetXmin(); + Float_t xmax = h1->GetXaxis()->GetXmax(); + Float_t xt = xmin + 0.67*(xmax-xmin); + + l.SetTextColor(kGreen); + l.DrawLatex(xt,0.095*maxy,"Engine"); + l.SetTextColor(kBlue); + l.DrawLatex(xt,0.045*maxy,"hcana"); + + // Difference between the histograms. + + TCanvas *c2 = new TCanvas("c2", "X differences", 1000, 667); + + TH1F* dif = h->Clone(); + + dif->Add(h,h1,1.,-1.); + + dif->SetTitle("X Difference"); + dif->SetFillColor(kRed); + dif->SetLineColor(kRed); + dif->SetLineWidth(1); + dif->SetFillStyle(1111); + dif->Draw(); +} diff --git a/src/THcShower.cxx b/src/THcShower.cxx index b49a109136009ef94aa54603477bf3e7cf2588a3..07dd8576a52e5ef07e1ba60864b4292af4e2ac6a 100644 --- a/src/THcShower.cxx +++ b/src/THcShower.cxx @@ -173,15 +173,15 @@ Int_t THcShower::ReadDatabase( const TDatime& date ) BlockThick = new Double_t [fNLayers]; fNBlocks = new Int_t [fNLayers]; fNLayerZPos = new Double_t [fNLayers]; - XPos = new Double_t [2*fNLayers]; + YPos = new Double_t [2*fNLayers]; for(Int_t i=0;i<fNLayers;i++) { DBRequest list[]={ {Form("cal_%s_thick",fLayerNames[i]), &BlockThick[i], kDouble}, {Form("cal_%s_nr",fLayerNames[i]), &fNBlocks[i], kInt}, {Form("cal_%s_zpos",fLayerNames[i]), &fNLayerZPos[i], kDouble}, - {Form("cal_%s_left",fLayerNames[i]), &XPos[2*i], kDouble}, - {Form("cal_%s_right",fLayerNames[i]), &XPos[2*i+1], kDouble}, + {Form("cal_%s_right",fLayerNames[i]), &YPos[2*i], kDouble}, + {Form("cal_%s_left",fLayerNames[i]), &YPos[2*i+1], kDouble}, {0} }; gHcParms->LoadParmValues((DBRequest*)&list, prefix); @@ -189,11 +189,11 @@ Int_t THcShower::ReadDatabase( const TDatime& date ) //Caution! Z positions (fronts) are off in hcal.param! Correct later on. - YPos = new Double_t* [fNLayers]; + XPos = new Double_t* [fNLayers]; for(Int_t i=0;i<fNLayers;i++) { - YPos[i] = new Double_t [fNBlocks[i]]; + XPos[i] = new Double_t [fNBlocks[i]]; DBRequest list[]={ - {Form("cal_%s_top",fLayerNames[i]),YPos[i], kDouble, fNBlocks[i]}, + {Form("cal_%s_top",fLayerNames[i]),XPos[i], kDouble, fNBlocks[i]}, {0} }; gHcParms->LoadParmValues((DBRequest*)&list, prefix); @@ -204,10 +204,10 @@ Int_t THcShower::ReadDatabase( const TDatime& date ) cout << " Block thickness: " << BlockThick[i] << endl; cout << " NBlocks : " << fNBlocks[i] << endl; cout << " Z Position : " << fNLayerZPos[i] << endl; - cout << " X Positions : " << XPos[2*i] << ", " << XPos[2*i+1] << endl; - cout << " Y Positions :"; + cout << " Y Positions : " << YPos[2*i] << ", " << YPos[2*i+1] <<endl; + cout << " X Positions :"; for(Int_t j=0; j<fNBlocks[i]; j++) { - cout << " " << YPos[i][j]; + cout << " " << XPos[i][j]; } cout << endl; @@ -392,8 +392,9 @@ Int_t THcShower::DefineVariables( EMode mode ) // { "asum_c", "Sum of calibrated ADCs", "fAsum_c" }, { "nclust", "Number of clusters", "fNclust" }, { "emax", "Energy (MeV) of largest cluster", "fE" }, - // { "x", "x-position (cm) of largest cluster", "fX" }, - // { "y", "y-position (cm) of largest cluster", "fY" }, + { "eprmax", "Preshower Energy (MeV) of largest cluster", "fEpr" }, + { "xmax", "x-position (cm) of largest cluster", "fX" }, + // { "z", "z-position (cm) of largest cluster", "fZ" }, { "mult", "Multiplicity of largest cluster", "fMult" }, // { "nblk", "Numbers of blocks in main cluster", "fNblk" }, // { "eblk", "Energies of blocks in main cluster", "fEblk" }, @@ -430,7 +431,7 @@ void THcShower::DeleteArrays() delete [] fNBlocks; fNBlocks = NULL; delete [] fNLayerZPos; fNLayerZPos = NULL; delete [] XPos; XPos = NULL; - delete [] YPos; YPos = NULL; + delete [] ZPos; ZPos = NULL; //delete [] fSpacing; fSpacing = NULL; //delete [] fCenter; fCenter = NULL; // This 2D. What is correct way to delete? } @@ -452,6 +453,8 @@ void THcShower::Clear(Option_t* opt) fNclust = 0; fMult = 0; fE = 0.; + fEpr = 0.; + fX = -75.; //out of acceptance } //_____________________________________________________________________________ @@ -541,9 +544,9 @@ Int_t THcShower::CoarseProcess( TClonesArray& ) //tracks // // Float_t Edep = fPlanes[j]->GetEmean(i); // if (Edep > 0.) { //hit - // Float_t y = YPos[j][i] + BlockThick[j]/2.; //top + thick/2 + // Float_t x = YPos[j][i] + BlockThick[j]/2.; //top + thick/2 // Float_t z = fNLayerZPos[j] + BlockThick[j]/2.; //front + thick/2 - // THcShowerHit* hit = new THcShowerHit(i,j,y,z,Edep); + // THcShowerHit* hit = new THcShowerHit(i,j,x,z,Edep); //ENGINE way. // @@ -554,13 +557,13 @@ Int_t THcShower::CoarseProcess( TClonesArray& ) //tracks fPlanes[j]->GetAneg(i) - fPlanes[j]->GetNegPed(i) > fPlanes[j]->GetNegThr(i) - fPlanes[j]->GetNegPed(i)) { //hit Float_t Edep = fPlanes[j]->GetEmean(i); - Float_t y = YPos[j][i] + BlockThick[j]/2.; //top + thick/2 + Float_t x = XPos[j][i] + BlockThick[j]/2.; //top + thick/2 Float_t z = fNLayerZPos[j] + BlockThick[j]/2.; //front + thick/2 - THcShowerHit* hit = new THcShowerHit(i,j,y,z,Edep); + THcShowerHit* hit = new THcShowerHit(i,j,x,z,Edep); HitList.push_back(hit); - //cout << "Hit: Edep = " << Edep << " Y = " << y << " Z = " << z << + //cout << "Hit: Edep = " << Edep << " X = " << x << " Z = " << z << // " Block " << i << " Layer " << j << endl; }; @@ -591,7 +594,7 @@ Int_t THcShower::CoarseProcess( TClonesArray& ) //tracks cout << "Cluster #" << i <<": E=" << (*cluster).clE() << " Epr=" << (*cluster).clEpr() - << " Y=" << (*cluster).clY() + << " X=" << (*cluster).clX() << " Z=" << (*cluster).clZ() << " size=" << (*cluster).clSize() << endl; @@ -624,13 +627,17 @@ Int_t THcShower::CoarseProcess( TClonesArray& ) //tracks } fE = (*MaxCluster).clE(); + fEpr = (*MaxCluster).clEpr(); + fX = (*MaxCluster).clX(); } + cout << fEpr << " " << fE << " " << fX << " PrSh" << endl; + /* cout << "Cluster #" << i <<": E=" << (*cluster).clE() << " Epr=" << (*cluster).clEpr() - << " Y=" << (*cluster).clY() + << " X=" << (*cluster).clX() << " Z=" << (*cluster).clZ() << " size=" << (*cluster).clSize() << endl; diff --git a/src/THcShower.h b/src/THcShower.h index 24f4c5601147cb1a26a37d3a08c61b9162093de0..f6ebd34cd8400a59d0ca02a774aeb164bd150e00 100644 --- a/src/THcShower.h +++ b/src/THcShower.h @@ -93,8 +93,9 @@ protected: Float_t fAsum_c; // Sum of calibrated ADCs Int_t fNclust; // Number of clusters Float_t fE; // Energy (MeV) of largest cluster + Float_t fEpr; // Preshower Energy (MeV) of largest cluster Float_t fX; // x-position (cm) of largest cluster - Float_t fY; // y-position (cm) of largest cluster + Float_t fZ; // z-position (cm) of largest cluster Int_t fMult; // Multiplicity of largest cluster Int_t fNblk; // Number of blocks in main cluster Float_t* fEblk; // Energies of blocks in main cluster @@ -110,8 +111,9 @@ protected: Double_t* BlockThick; // Thickness of shower counter blocks, blocks Int_t* fNBlocks; // Number of shower counter blocks per layer Int_t fNtotBlocks; // Total number of shower counter blocks - Double_t** YPos; //X,Y positions of shower counter blocks - Double_t* XPos; + Double_t** XPos; //X,Y,Z positions of shower counter blocks + Double_t* YPos; + Double_t* ZPos; Int_t fNegCols; //number of columns with PMTTs on the negative side only. Double_t fSlop; //Track to cluster vertical slop distance. Int_t fvTest; //fiducial volume test flag diff --git a/src/THcShowerCluster.h b/src/THcShowerCluster.h index 513e71c9860d0d47b8e3b5abb02099c05db0d40e..61213e5f9793cf851582046a526f5b5f4804c046 100644 --- a/src/THcShowerCluster.h +++ b/src/THcShowerCluster.h @@ -44,18 +44,18 @@ class THcShowerCluster : THcShowerHitList { (*(THcShowerHitList::begin()+num))->show(); } - //Y coordinate of cluster's center of gravity. + //X coordinate of cluster's center of gravity. // - float clY() { - float y_sum=0.; + float clX() { + float x_sum=0.; float Etot=0.; for (THcShowerHitIt it=THcShowerHitList::begin(); it!=THcShowerHitList::end(); it++) { - y_sum += (*it)->hitY() * (*it)->hitE(); + x_sum += (*it)->hitX() * (*it)->hitE(); Etot += (*it)->hitE(); } - // cout << "y_sum=" << y_sum << " Etot=" << Etot << endl; - return y_sum/Etot; + // cout << "x_sum=" << x_sum << " Etot=" << Etot << endl; + return (Etot != 0. ? x_sum/Etot : -75.); } //Z coordinate for a cluster, calculated as a weighted by energy average. diff --git a/src/THcShowerHit.h b/src/THcShowerHit.h index d2504d23868767147fa1494df9cea82529ed4cbf..1437699903ca603a2f633ce211a79e541ba223e5 100644 --- a/src/THcShowerHit.h +++ b/src/THcShowerHit.h @@ -13,22 +13,22 @@ class THcShowerHit { //HMS calorimeter hit class private: unsigned int fCol, fRow; //hit colomn and row - float fY, fZ; //hit Y (vert.) and Z (along spect.axis) coordinates + float fX, fZ; //hit X (vert.) and Z (along spect.axis) coordinates float fE; //hit energy deposition public: THcShowerHit() { //default constructor fCol=fRow=0; - fY=fZ=0.; + fX=fZ=0.; fE=0.; } - THcShowerHit(unsigned int hRow, unsigned int hCol, float hY, float hZ, + THcShowerHit(unsigned int hRow, unsigned int hCol, float hX, float hZ, float hE) { fRow=hRow; fCol=hCol; - fY=hY; + fX=hX; fZ=hZ; fE=hE; } @@ -45,8 +45,8 @@ class THcShowerHit { //HMS calorimeter hit class return fRow; } - float hitY() { - return fY; + float hitX() { + return fX; } float hitZ() { @@ -67,7 +67,7 @@ class THcShowerHit { //HMS calorimeter hit class // void show() { cout << "row=" << fRow << " column=" << fCol << - " y=" << fY << " z=" << fZ << " E=" << fE << endl; + " x=" << fX << " z=" << fZ << " E=" << fE << endl; } };