Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-08-05 08:12:44 

0001 
0002 
0003 #include "SaveCanvas.C"
0004 #include "sPhenixStyle.C"
0005 
0006 #include <TCanvas.h>
0007 #include <TDatabasePDG.h>
0008 #include <TFile.h>
0009 #include <TH1D.h>
0010 #include <TH1F.h>
0011 #include <TH2F.h>
0012 #include <TH3F.h>
0013 #include <TLegend.h>
0014 #include <TLine.h>
0015 #include <TLorentzVector.h>
0016 #include <TRandom3.h>
0017 #include <TString.h>
0018 #include <TTree.h>
0019 #include <TVirtualFitter.h>
0020 
0021 #include <assert.h>
0022 
0023 #include <cmath>
0024 #include <cstddef>
0025 #include <iostream>
0026 
0027 using namespace std;
0028 
0029 TFile *_file0 = NULL;
0030 TFile *_file_trigger = NULL;
0031 TString description;
0032 
0033 typedef vector<vector<TH1 *> > chipMultiplicitySet_vec;
0034 
0035 chipMultiplicitySet_vec MakeChipMultiplicitySet(TFile *file)
0036 {
0037   assert(file);
0038 
0039   TH3 *LayerChipMultiplicity = (TH3 *) file->GetObjectChecked("LayerChipMultiplicity", "TH3");
0040 
0041   TCanvas *c1 = new TCanvas("MakeChipMultiplicitySet", "MakeChipMultiplicitySet", 1800, 960);
0042   c1->Divide(9, 3, 0, 0);
0043   int idx = 1;
0044   TPad *p;
0045 
0046   chipMultiplicitySet_vec chipMultiplicitySet;
0047 
0048   for (int ilayer = 1; ilayer <= LayerChipMultiplicity->GetNbinsX(); ++ilayer)
0049   {
0050     vector<TH1 *> chipMultiplicityLayer;
0051 
0052     for (int ichip = 1; ichip <= LayerChipMultiplicity->GetNbinsY(); ++ichip)
0053     {
0054       p = (TPad *) c1->cd(idx++);
0055       c1->Update();
0056       p->SetLogy();
0057 
0058       LayerChipMultiplicity->GetXaxis()->SetRange(ilayer, ilayer);
0059       LayerChipMultiplicity->GetYaxis()->SetRange(ichip, ichip);
0060 
0061       TString histname(Form("ChipMultiplicity_Layer%dChip%d", ilayer - 1, ichip - 1));
0062       TH1 *h = LayerChipMultiplicity->Project3D("z");
0063       h->SetTitle(";Chip Multiplicity [pixel];Count");
0064 
0065       h->SetName(histname);
0066 
0067       h->Draw();
0068 
0069       h->SetDirectory(NULL);
0070       h->ComputeIntegral(true);
0071 
0072       chipMultiplicityLayer.push_back(h);
0073     }
0074 
0075     chipMultiplicitySet.push_back(chipMultiplicityLayer);
0076   }
0077 
0078   SaveCanvas(c1, TString(file->GetName()) + TString(c1->GetName()), false);
0079 
0080   c1 = new TCanvas("MakeChipMultiplicitySetChip4", "MakeChipMultiplicitySetChip4", 1200, 960);
0081   gPad->SetLogy();
0082   TH1 *LayerMultiplicityLayer0 = chipMultiplicitySet[0][4];
0083   TH1 *LayerMultiplicityLayer1 = chipMultiplicitySet[1][4];
0084   TH1 *LayerMultiplicityLayer2 = chipMultiplicitySet[2][4];
0085 
0086   LayerMultiplicityLayer0->SetLineColor(kRed + 2);
0087   LayerMultiplicityLayer1->SetLineColor(kBlue + 2);
0088   LayerMultiplicityLayer2->SetLineColor(kGreen + 2);
0089 
0090   LayerMultiplicityLayer0->Draw();
0091   LayerMultiplicityLayer1->Draw("same");
0092   LayerMultiplicityLayer2->Draw("same");
0093 
0094   LayerMultiplicityLayer0->SetTitle(";Chip multiplicity [Pixel];Count");
0095 
0096   TLegend *leg = new TLegend(.5, .5, .93, .93);
0097   leg->AddEntry("", "#it{#bf{sPHENIX}} Simulation", "");
0098   leg->AddEntry("", description, "");
0099   leg->AddEntry(LayerMultiplicityLayer0, Form("MVTX Layer0 Chip4, <hit> = %.1f", LayerMultiplicityLayer0->GetMean()), "l");
0100   leg->AddEntry(LayerMultiplicityLayer1, Form("MVTX Layer1 Chip4, <hit> = %.1f", LayerMultiplicityLayer1->GetMean()), "l");
0101   leg->AddEntry(LayerMultiplicityLayer2, Form("MVTX Layer2 Chip4, <hit> = %.1f", LayerMultiplicityLayer2->GetMean()), "l");
0102   leg->Draw();
0103 
0104   SaveCanvas(c1, TString(file->GetName()) + TString(c1->GetName()), false);
0105 
0106   return chipMultiplicitySet;
0107 }
0108 
0109 void Check(TFile *file)
0110 {
0111   assert(file);
0112 
0113   TCanvas *c1 = new TCanvas("Check", "Check", 1800, 960);
0114   c1->Divide(2, 2);
0115   int idx = 1;
0116   TPad *p;
0117 
0118   p = (TPad *) c1->cd(idx++);
0119   c1->Update();
0120   //  p->SetLogz();
0121 
0122   TH1 *hNormalization = (TH1 *) file->GetObjectChecked("hNormalization", "TH1");
0123   hNormalization->Draw();
0124 
0125   p = (TPad *) c1->cd(idx++);
0126   c1->Update();
0127   //  p->SetLogy();
0128 
0129   TH1 *hNChEta = (TH1 *) file->GetObjectChecked("hNChEta", "TH1");
0130   assert(hNChEta);
0131   double norm = hNChEta->GetBinWidth(1) * hNormalization->GetBinContent(2);
0132   hNChEta->Scale(1. / norm);
0133   hNChEta->Draw();
0134 
0135   hNChEta->GetYaxis()->SetTitle("dN_{Ch}/d#eta");
0136   //  hNormalization->Draw()
0137 
0138   p = (TPad *) c1->cd(idx++);
0139   c1->Update();
0140   TH1 *hVertexZ = (TH1 *) file->GetObjectChecked("hVertexZ", "TH1");
0141   hVertexZ->Draw();
0142 
0143   p = (TPad *) c1->cd(idx++);
0144   c1->Update();
0145   p->SetLogy();
0146 
0147   TH2 *LayerMultiplicity = (TH2 *) file->GetObjectChecked("LayerMultiplicity", "TH2");
0148   TH1 *LayerMultiplicityLayer0 =
0149       LayerMultiplicity->ProjectionY("LayerMultiplicityLayer0", 1, 1);
0150   TH1 *LayerMultiplicityLayer1 =
0151       LayerMultiplicity->ProjectionY("LayerMultiplicityLayer1", 2, 2);
0152   TH1 *LayerMultiplicityLayer2 =
0153       LayerMultiplicity->ProjectionY("LayerMultiplicityLayer2", 3, 3);
0154 
0155   LayerMultiplicityLayer0->SetLineColor(kRed + 2);
0156   LayerMultiplicityLayer1->SetLineColor(kBlue + 2);
0157   LayerMultiplicityLayer2->SetLineColor(kGreen + 2);
0158 
0159   LayerMultiplicityLayer0->Draw();
0160   LayerMultiplicityLayer1->Draw("same");
0161   LayerMultiplicityLayer2->Draw("same");
0162 
0163   LayerMultiplicityLayer0->SetTitle(";Chip multiplicity [Pixel];Count");
0164 
0165   TLegend *leg = new TLegend(.5, .5, .93, .93);
0166   leg->AddEntry("", "#it{#bf{sPHENIX}} Simulation", "");
0167   leg->AddEntry("", description, "");
0168   leg->AddEntry(LayerMultiplicityLayer0, Form("MVTX Layer0, <hit> = %.1f", LayerMultiplicityLayer0->GetMean()), "l");
0169   leg->AddEntry(LayerMultiplicityLayer1, Form("MVTX Layer1, <hit> = %.1f", LayerMultiplicityLayer1->GetMean()), "l");
0170   leg->AddEntry(LayerMultiplicityLayer2, Form("MVTX Layer2, <hit> = %.1f", LayerMultiplicityLayer2->GetMean()), "l");
0171   leg->Draw();
0172 
0173   SaveCanvas(c1, TString(file->GetName()) + TString(c1->GetName()), false);
0174 }
0175 
0176 TH1 *MakeCDF(TH1 *h)
0177 {
0178   assert(h);
0179 
0180   TH1 *hCDF = (TH1 *) h->Clone(TString("CDF") + h->GetName());
0181   hCDF->SetDirectory(NULL);
0182 
0183   double integral = 0;
0184 
0185   for (int bin = h->GetNbinsX() + 1; bin >= 0; --bin)
0186   {
0187     integral += h->GetBinContent(bin);
0188     hCDF->SetBinContent(bin, integral);
0189   }
0190 
0191   for (int bin = h->GetNbinsX(); bin >= 1; --bin)
0192   {
0193     hCDF->SetBinContent(bin, hCDF->GetBinContent(bin) / integral);
0194   }
0195 
0196   return hCDF;
0197 }
0198 
0199 chipMultiplicitySet_vec TriggerMultiplicity(chipMultiplicitySet_vec cm_MB, chipMultiplicitySet_vec cm_Trigger, const double mu_MB, const int n_Trigger, const double mu_Noise)
0200 {
0201   //init
0202   chipMultiplicitySet_vec cm;
0203 
0204   const static int MaxMultiplicity(2000);
0205 
0206   for (auto &layer : cm_MB)
0207   {
0208     vector<TH1 *> chipMultiplicityLayer;
0209 
0210     for (auto &chip : layer)
0211     {
0212       TH1 *h = new TH1D(
0213           TString("Trigger") + chip->GetName(), TString("Trigger") + chip->GetName(),
0214           MaxMultiplicity, -.5, MaxMultiplicity - .5);
0215       //          (TH1 *) chip->Clone(TString("Trigger") + chip->GetName());
0216 
0217       h->SetDirectory(NULL);
0218       h->GetXaxis()->SetTitle(chip->GetXaxis()->GetTitle());
0219       h->GetYaxis()->SetTitle(chip->GetYaxis()->GetTitle());
0220 
0221       assert(h);
0222       chipMultiplicityLayer.push_back(h);
0223     }
0224     cm.push_back(chipMultiplicityLayer);
0225   }
0226 
0227   //composition for trigger
0228   TRandom3 rnd(1234);
0229   assert(mu_MB >= 0);
0230   assert(n_Trigger >= 0 && n_Trigger <= 1);
0231   assert(mu_Noise >= 0);
0232 
0233   //  const int NSample = 1000000;
0234   const int NSample = 10000000;
0235   for (int i = 0; i < NSample; ++i)
0236   {
0237     if (i % (NSample / 10) == 0)
0238     {
0239       cout << "TriggerMultiplicity - " << i << endl;
0240     }
0241 
0242     int n_MB = rnd.Poisson(mu_MB);
0243     int n_Noise = rnd.Poisson(mu_Noise);
0244 
0245     int ilayer = 0;
0246     for (auto &layer : cm)
0247     {
0248       int ichip = 0;
0249       for (auto &chip : layer)
0250       {
0251         assert(chip);
0252         int n_hit = n_Noise;
0253 
0254         assert(ilayer < cm_Trigger.size());
0255         assert(ichip < cm_Trigger[ilayer].size());
0256         assert(cm_Trigger[ilayer][ichip]);
0257         if (n_Trigger)
0258         {
0259           double rndHit = cm_Trigger[ilayer][ichip]->GetRandom();
0260           if (rndHit < .5) rndHit = 0;
0261 
0262           //          cout << "n_hit += " << cm_Trigger[ilayer][ichip]->GetRandom() << endl;
0263           n_hit += rndHit;
0264         }
0265 
0266         assert(ilayer < cm_MB.size());
0267         assert(ichip < cm_MB[ilayer].size());
0268         assert(cm_MB[ilayer][ichip]);
0269         for (int iMB = 0; iMB < n_MB; ++iMB)
0270         {
0271           double rndHit = cm_MB[ilayer][ichip]->GetRandom();
0272           if (rndHit < .5) rndHit = 0;
0273           //          cout << "n_hit +=" << cm_MB[ilayer][ichip]->GetRandom() << endl;
0274           n_hit += rndHit;
0275         }
0276 
0277         //        cout << "chip->Fill(n_hit);" << n_hit << ", n_Noise = " << n_Noise << endl;
0278         chip->Fill(n_hit);
0279 
0280         ++ichip;
0281       }
0282       ++ilayer;
0283     }
0284   }
0285 
0286   TString UniqueName(Form("_muMB%.0f_nTrig%d_muNoise%.0f", mu_MB * 100, n_Trigger, mu_Noise));
0287 
0288   //save - plot
0289   TCanvas *c1 = new TCanvas("TriggerMultiplicity" + UniqueName, "TriggerMultiplicity" + UniqueName, 1800, 960);
0290   c1->Divide(9, 3, 0, 0);
0291   int idx = 1;
0292   TPad *p;
0293 
0294   for (auto &layer : cm)
0295   {
0296     vector<TH1 *> chipMultiplicityLayer;
0297 
0298     for (auto &chip : layer)
0299     {
0300       p = (TPad *) c1->cd(idx++);
0301       c1->Update();
0302       p->SetLogy();
0303 
0304       chip->SetTitle(";Chip Multiplicity [pixel];Count");
0305       chip->Draw();
0306     }
0307   }
0308 
0309   SaveCanvas(c1, TString(_file0->GetName()) + TString(c1->GetName()), false);
0310 
0311   c1 = new TCanvas("TriggerMultiplicityChip4" + UniqueName, "TriggerMultiplicityChip4" + UniqueName, 1900, 750);
0312   c1->Divide(2, 1);
0313   idx = 1;
0314 
0315   p = (TPad *) c1->cd(idx++);
0316   c1->Update();
0317   p->SetLogy();
0318 
0319   TH1 *LayerMultiplicityLayer0 = cm[0][4];
0320   TH1 *LayerMultiplicityLayer1 = cm[1][4];
0321   TH1 *LayerMultiplicityLayer2 = cm[2][4];
0322 
0323   LayerMultiplicityLayer0->SetLineColor(kRed + 2);
0324   LayerMultiplicityLayer1->SetLineColor(kBlue + 2);
0325   LayerMultiplicityLayer2->SetLineColor(kGreen + 2);
0326 
0327   LayerMultiplicityLayer0->Draw();
0328   LayerMultiplicityLayer1->Draw("same");
0329   LayerMultiplicityLayer2->Draw("same");
0330 
0331   LayerMultiplicityLayer0->SetTitle(";Chip multiplicity [Pixel];Count");
0332 
0333   TLegend *leg = new TLegend(.45, .55, .93, .93);
0334   leg->AddEntry("", "#it{#bf{sPHENIX}} Simulation", "");
0335   leg->AddEntry("", description, "");
0336   leg->AddEntry("", Form("#mu_{MB} = %.1f, N_{Trig} = %d, #mu_{Noise} = %.1f", mu_MB, n_Trigger, mu_Noise), "");
0337   leg->AddEntry(LayerMultiplicityLayer0, Form("MVTX Layer0 Chip4, <hit> = %.1f", LayerMultiplicityLayer0->GetMean()), "l");
0338   leg->AddEntry(LayerMultiplicityLayer1, Form("MVTX Layer1 Chip4, <hit> = %.1f", LayerMultiplicityLayer1->GetMean()), "l");
0339   leg->AddEntry(LayerMultiplicityLayer2, Form("MVTX Layer2 Chip4, <hit> = %.1f", LayerMultiplicityLayer2->GetMean()), "l");
0340   leg->Draw();
0341 
0342   p = (TPad *) c1->cd(idx++);
0343   c1->Update();
0344   p->SetLogy();
0345 
0346   p->DrawFrame(0, 1. / NSample, MaxMultiplicity, 1)->SetTitle(";Chip multiplicity [Pixel];CCDF");
0347   TLine *line = new TLine(0, 1e-3, MaxMultiplicity, 1e-3);
0348   line->SetLineColor(kGray);
0349   line->SetLineWidth(3);
0350   line->Draw("same");
0351 
0352   MakeCDF(LayerMultiplicityLayer0)->Draw("same");
0353   MakeCDF(LayerMultiplicityLayer1)->Draw("same");
0354   MakeCDF(LayerMultiplicityLayer2)->Draw("same");
0355 
0356   SaveCanvas(c1, TString(_file0->GetName()) + TString(c1->GetName()), false);
0357 
0358   return cm;
0359 }
0360 
0361 void OccupancySim(                                                                                                                                                     //
0362                                                                                                                                                                        //    const TString infile = "/sphenix/user/jinhuang/HF-jet/MVTX_Multiplicity/pp200MB_30cmVZ_Iter5/pp200MB_30cmVZ_Iter5_SUM.cfg_HFMLTriggerOccupancy.root",          //
0363                                                                                                                                                                        //    const TString infile_trigger = "/sphenix/user/jinhuang/HF-jet/MVTX_Multiplicity/pp200MB_30cmVZ_Iter5/pp200MB_30cmVZ_Iter5_SUM.cfg_HFMLTriggerOccupancy.root",  //
0364                                                                                                                                                                        //    const TString disc = "p+p MB, #sqrt{s} = 200 GeV"                                                                                                              //
0365     const TString infile = "/sphenix/user/jinhuang/HF-jet/MVTX_Multiplicity/AuAu200MB_30cmVZ_Iter5/AuAu200MB_30cmVZ_Iter5_SUM.xml_HFMLTriggerOccupancy.root",          //
0366     const TString infile_trigger = "/sphenix/user/jinhuang/HF-jet/MVTX_Multiplicity/AuAu200MB_10cmVZ_Iter5/AuAu200MB_10cmVZ_Iter5_SUM.xml_HFMLTriggerOccupancy.root",  //
0367     const TString disc = "Au+Au MB, #sqrt{s_{NN}} = 200 GeV"                                                                                                           //
0368 )
0369 {
0370   SetsPhenixStyle();
0371   TVirtualFitter::SetDefaultFitter("Minuit2");
0372 
0373   description = disc;
0374   _file0 = new TFile(infile);
0375   assert(_file0->IsOpen());
0376   _file_trigger = new TFile(infile_trigger);
0377   assert(_file_trigger->IsOpen());
0378 
0379   _file0->cd();
0380   Check(_file0);
0381   chipMultiplicitySet_vec chipMultiplicity = MakeChipMultiplicitySet(_file0);
0382 
0383   _file_trigger->cd();
0384   Check(_file_trigger);
0385   chipMultiplicitySet_vec chipMultiplicityTrigger = MakeChipMultiplicitySet(_file_trigger);
0386 
0387   // AuAu
0388   TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0389                       200e3 * 10e-6, 1, 1024 * 512 * 1e-6 + 2);
0390 
0391   TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0392                       200e3 * 15e-6, 0, 1024 * 512 * 1e-6 + 2);
0393 
0394   TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0395                       200e3 * 10e-6, 0, 1024 * 512 * 1e-6 + 2);
0396 
0397   TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0398                       200e3 * 5e-6, 0, 1024 * 512 * 1e-6 + 2);
0399 
0400   TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0401                       100e3 * 5e-6, 0, 1024 * 512 * 1e-6 + 2);
0402 
0403   TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0404                       50e3 * 5e-6, 0, 1024 * 512 * 1e-6 + 2);
0405 
0406   TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0407                       200e3 * 15e-6, 0, 1024 * 512 * 1e-5 + 2);
0408 
0409   TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0410                       200e3 * 10e-6, 0, 1024 * 512 * 1e-5 + 2);
0411 
0412   TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0413                       200e3 * 5e-6, 0, 1024 * 512 * 1e-5 + 2);
0414 
0415   //  //   pp
0416   //  TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0417   //                      13e6 * 10e-6, 1, 1024 * 512 * 1e-6);
0418   //
0419   //  TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0420   //                      13e6 * 15e-6, 0, 1024 * 512 * 1e-6);
0421   //  TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0422   //                      13e6 * 10e-6, 0, 1024 * 512 * 1e-6);
0423   //  TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0424   //                      13e6 * 5e-6, 0, 1024 * 512 * 1e-6);
0425   //
0426   //  TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0427   //                      13e6 * 15e-6, 0, 1024 * 512 * 1e-5);
0428   //  TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0429   //                      13e6 * 10e-6, 0, 1024 * 512 * 1e-5);
0430   //  TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0431   //                      13e6 * 5e-6, 0, 1024 * 512 * 1e-5);
0432   //  TriggerMultiplicity(chipMultiplicity, chipMultiplicityTrigger,
0433   //                      4e6 * 5e-6, 0, 1024 * 512 * 1e-5);
0434 
0435   //
0436 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!