test/vmctest/scripts/cuts.cxx

   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
  18 // Program to generate the plots with energy cuts and corresponding
  19 // cuts in range for all materials per detector.
  20 // The program uses a text file output from geant4_vmc
  21 // (in development version only).
  22 //
  23 // By I.Hrivnacova, IPN Orsay
  24
  25 #include <Rtypes.h>
  26 #include <TPaveText.h>
  27 #include <TString.h>
  28 #include <TGraph.h>
  29 #include <TCanvas.h>
  30 #include <TAxis.h>
  31 #include <TLegend.h>
  32
  33 #include <map>
  34 #include <vector>
  35 #include <string>
  36 #include <fstream>
  37 #include <iostream>
  38
  39 using namespace std;
  40
  41 struct MaterialData {
  42   MaterialData(string name,
  43                Double_t x1, Double_t x2, Double_t x3,
  44                Double_t x4, Double_t x5, Double_t x6)
  45     : matName(name),
  46       rangeGam(x1),
  47       rangeEle(x2),
  48       cutGam(x3),
  49       cutEle(x4),
  50       cutGamLimit(x5),
  51       cutEleLimit(x6) {}
  52   string matName;
  53   Double_t rangeGam;
  54   Double_t rangeEle;
  55   Double_t cutGam;
  56   Double_t cutEle;
  57   Double_t cutGamLimit;
  58   Double_t cutEleLimit;
  59 };
  60
  61 map<string, vector<MaterialData> > materialDataMap;
  62
  63 void readCuts()
  64 {
  65   std::ifstream input("cuts.txt");
  66   while ( ! input.eof() ) {
  67     string matName;
  68     Int_t number;
  69     Double_t x1, x2, x3, x4, x5, x6;
  70     input >> number;
  71     if ( input.eof() ) break;
  72     input  >> matName >> x1 >> x2 >> x3 >> x4 >> x5 >> x6;
  73     // cout << ".. reading " << number << " " << matName << endl;
  74
  75     if ( x1 > 1e10 ) x1 = 1e04;
  76     if ( x2 > 1e10 ) x2 = 1e04;
  77     if ( x3 > 1e10 ) x3 = 1e04;
  78     if ( x4 > 1e10 ) x4 = 1e04;
  79     if ( x5 > 1e10 ) x5 = 1e04;
  80     if ( x6 > 1e10 ) x6 = 1e04;
  81
  82     string detName = matName;
  83     detName.erase(detName.find('_'), detName.size()-detName.find('_'));
  84
  85     std::map<string, std::vector<MaterialData> >::iterator it
  86       = materialDataMap.find(detName);
  87     if ( it == materialDataMap.end() ) {
  88       materialDataMap.insert(
  89         std::pair<string, vector<MaterialData> >(detName, vector<MaterialData>()));
  90     }
  91     materialDataMap[detName].
  92       push_back(MaterialData(matName, x1, x2, x3, x4, x5, x6));
  93   }
  94 }
  95
  96 void plotCuts(const string& detName)
  97 {
  98   cout << "Processing " << detName << " ..." << endl;
  99
 100   //if ( detName == "DIPO" ) return;
 101   //if ( detName == "HALL" ) return;
 102   //if ( detName == "SHIL" ) return;
 103   //if ( detName == "ZDC" ) return;
 104
 105   TPaveText* paveText = new TPaveText(0.1, 0.1, 0.98, 0.98);
 106   paveText->SetTextAlign(11);
 107
 108   std::vector<MaterialData> materialData
 109     = materialDataMap[detName];
 110
 111   Double_t matNumber[300];
 112   Double_t rangeGam[300];
 113   Double_t rangeEle[300];
 114   Double_t cutGam[300];
 115   Double_t cutEle[300];
 116   Double_t cutGamLimit[300];
 117   Double_t cutEleLimit[300];
 118   std::vector<string> matNames;
 119
 120   for (UInt_t i=0; i<materialData.size(); i++ ) {
 121
 122     matNumber[i] = i;
 123     rangeGam[i] = materialData[i].rangeGam;
 124     rangeEle[i] = materialData[i].rangeEle;
 125     cutGam[i] = materialData[i].cutGam;
 126     cutEle[i] = materialData[i].cutEle;
 127     cutGamLimit[i] = materialData[i].cutGamLimit;
 128     cutEleLimit[i] = materialData[i].cutEleLimit;
 129
 130     TString legend;
 131     legend += i;
 132     legend += "  ";
 133     legend += materialData[i].matName.data();
 134     paveText->AddText(legend.Data());
 135   }
 136
 137   TGraph* gr1 = new TGraph(materialData.size(), matNumber, rangeGam);
 138   TGraph* gr2 = new TGraph(materialData.size(), matNumber, cutGam);
 139   TGraph* gr3 = new TGraph(materialData.size(), matNumber, cutGamLimit);
 140   TGraph* gr4 = new TGraph(materialData.size(), matNumber, rangeEle);
 141   TGraph* gr5 = new TGraph(materialData.size(), matNumber, cutEle);
 142   TGraph* gr6 = new TGraph(materialData.size(), matNumber, cutEleLimit);
 143
 144   // gamma range cut
 145   gr1->SetMarkerColor(kBlue);
 146   gr1->SetMarkerStyle(22);
 147   gr1->SetTitle("Range cut for gamma");
 148   gr1->GetXaxis()->SetTitle("Material number");
 149   gr1->GetYaxis()->SetTitle("Range [mm]");
 150   gr1->SetMinimum(1e-03);
 151   gr1->SetMaximum(1e+04);
 152
 153   // gamma energy threshold
 154   gr2->SetMarkerColor(kBlue);
 155   gr2->SetMarkerStyle(22);
 156   gr2->SetTitle("Energy threshold for gamma");
 157   gr2->GetXaxis()->SetTitle("Material number");
 158   gr2->GetYaxis()->SetTitle("Energy [MeV]");
 159   gr2->SetMinimum(1e-04);
 160   gr2->SetMaximum(1e+04);
 161
 162   // gamma user limit
 163   gr3->SetMarkerColor(kBlack);
 164   gr3->SetMarkerStyle(23);
 165   gr3->SetTitle("User limit for gamma (GUTGAM)");
 166   gr3->GetXaxis()->SetTitle("Material number");
 167   gr3->GetYaxis()->SetTitle("Energy [MeV]");
 168   gr3->SetMinimum(1e-04);
 169   gr3->SetMaximum(1e+04);
 170
 171   // e- range cut
 172   gr4->SetMarkerColor(kRed);
 173   gr4->SetMarkerStyle(22);
 174   gr4->SetTitle("Range cut for e-");
 175   gr4->GetXaxis()->SetTitle("Material number");
 176   gr4->GetYaxis()->SetTitle("Range [mm]");
 177   gr4->GetYaxis()->SetRange(1e-03, 1e+04);
 178   gr4->SetMinimum(1e-03);
 179   gr4->SetMaximum(1e+04);
 180
 181   // e- energy threshold
 182   gr5->SetMarkerColor(kRed);
 183   gr5->SetMarkerStyle(22);
 184   gr5->SetTitle("Energy threshold for e-");
 185   gr5->GetXaxis()->SetTitle("Material number");
 186   gr5->GetYaxis()->SetTitle("Energy [MeV]");
 187   gr5->SetMinimum(1e-04);
 188   gr5->SetMaximum(1e+04);
 189
 190   // e- user limit
 191   gr6->SetMarkerColor(kBlack);
 192   gr6->SetMarkerStyle(23);
 193   gr6->SetTitle("User limit for e- (CUTELE)");
 194   gr6->GetXaxis()->SetTitle("Material number");
 195   gr6->GetYaxis()->SetTitle("Energy [MeV]");
 196   gr6->SetMinimum(1e-04);
 197   gr6->SetMaximum(1e+04);
 198
 199   //TLegend* leg = new TLegend(0.1,0.7,0.48,0.9);
 200   TLegend* leg = new TLegend(0.1, 0.1, 0.98, 0.98);
 201   leg->SetHeader("Energy threshold/User cuts");
 202   leg->AddEntry(gr2,"Energy threshold for gamma","P");
 203   leg->AddEntry(gr3,"User cut for gamma (CUGAM)","P");
 204   leg->AddEntry(gr5,"Energy threshold for e-","P");
 205   leg->AddEntry(gr6,"User cut for e- (CUTELE)","P");
 206
 207   TCanvas* canvas = new TCanvas("Range cuts", "Range cuts", 1200, 800);
 208   canvas->Divide(3,2);
 209
 210   // Draw graph using the function Draw()
 211   canvas->cd(1);
 212   canvas->cd(1)->SetLogy();
 213   gr1->Draw("AP");
 214
 215   canvas->cd(2);
 216   canvas->cd(2)->SetLogy();
 217   gr2->Draw("AP");
 218
 219   //canvas->cd(3);
 220   //canvas->cd(3)->SetLogy();
 221   gr3->Draw("P");
 222
 223   canvas->cd(3);
 224   leg->Draw();
 225
 226   canvas->cd(4);
 227   canvas->cd(4)->SetLogy();
 228   gr4->Draw("AP");
 229
 230   canvas->cd(5);
 231   canvas->cd(5)->SetLogy();
 232   gr5->Draw("AP");
 233
 234   //canvas->cd(6);
 235   //canvas->cd(6)->SetLogy();
 236   gr6->Draw("P");
 237
 238   canvas->cd(6);
 239
 240   //TCanvas* canvas2 = new TCanvas("Materials", "Materials", 200, 800);
 241   //canvas2->cd();
 242   paveText->Draw();
 243
 244   string outName = detName;
 245   outName += "_cuts.gif";
 246   canvas->SaveAs(outName.data());
 247
 248   //string outName2 = detName;
 249   //outName2 += "_cuts_legend.gif";
 250   //canvas2->SaveAs(outName2.data());
 251 }
 252
 253 int main() {
 254
 255   readCuts();
 256
 257   std::map<string, std::vector<MaterialData> >::iterator it;
 258   for ( it= materialDataMap.begin(); it != materialDataMap.end(); it++ )
 259     plotCuts(it->first);
 260
 261   return 0;
 262 }
 263
 264
 265