Welcome Page Qweak
Region 3 Drift Chambers
Construction Software DAQ Internal Physics
  Software related Topics  Not logged in Qweak
Entry Michael Gericke; Cerenkov Detector Light Transmission Simulation
   Reply Michael Gericke; Screenshot of the handlebar detector
Message ID: 93     Entry time: 03/14/2005 02:06:47 PM     Reply to this: 100
 Author: Michael Gericke 
 Type: Simulation 
 Category: Geant4 
 Subject: Cerenkov Detector Light Transmission Simulation 
 Record date: 03/17/2005 10:31:11 AM 
I have run the GEANT4 QWeak MC with Cerenkov light enabled, for a straight detector, for a
V-Shaped detector, as well as for the new handle shaped detector. For the V-Shape the MC 
was run with several different tilt angles between 0 and 25 degrees w.r.t the vertical.

Attachment 1 below shows the results for the straight bar which, I believe, essentially 
agree with what Neven has seen. 

Attachment 2 shows the results for a V-Shaped detector with 0 tilt (vertical).

Attachment 3 shows the results for the handle-shaped detector at 0 tilt angle.

Attachment 4 shows the light transmission for each PMT and how it changes with tilt
angle (this is for the V-shaped detector). Each plot shows the entire length of the 
detector, so there are no dead spots, but the transmission from the end opposite to 
the PMT gets worse as the tilt angle increases.

The plots included in Attachments 1,2, and 3 show, in order from left to right and top 
to bottom, 

            1) Left PMT number of PEs vs. electron hit location in the detector
            2) Right PMT number of PEs vs. electron hit location in the detector
            3) Left + Right PMT number of PEs vs. electron hit location in the detector
            4) 2D electron hit location
            5) Left PMT number of PEs distribution
            6) Right PMT number of PEs distribution
            7) Distribution of primary electron Q2 weighted by Left + Right PMT number of PEs 
               for the event
            8) Distribution of primary electron Q2 without weight

I suppose the main significance for the last two plots (7 and 8) is that the distribution
of the number of PEs as a function of electron hit location is flat enough to not change
the Q2 distribution very much (for a vertical detector).  

For reference, the ROOT scripts that were used to create these plots are attached below.
Attachment 1: CerenkovLightStraight.gif  27 kB  | Hide | Hide all
CerenkovLightStraight.gif
Attachment 2: CerenkovLightVshape.gif  27 kB  | Hide | Hide all
CerenkovLightVshape.gif
Attachment 3: CerenkovLightHandle.gif  27 kB  | Hide | Hide all
CerenkovLightHandle.gif
Attachment 4: CerenkovLightVshapeTilt.gif  9 kB  | Hide | Hide all
CerenkovLightVshapeTilt.gif
Attachment 5: CerenkovTest5.C  5 kB  | Hide | Hide all
void CerenkovTest5()
{
//   gStyle->SetOptStat(kFALSE);
//   TFile *f = new TFile("QweakSim_03_16_2000_tilt20.8.root","READ");
  TFile *f = new TFile("QweakSim_03_17_2000_tilt0.0.root","READ");
  TTree *tree = (TTree *) gROOT->FindObject("Test_NT");
  Int_t NEvents = tree->GetEntries();
  Double_t PI            = 3.14159265359;  
  Double_t *XPos         = new Double_t[NEvents];
  Double_t *YPos         = new Double_t[NEvents];
  Double_t *LeftPMTHits  = new Double_t[NEvents];
  Double_t *RightPMTHits = new Double_t[NEvents];

  Double_t angle         = 11.31*PI/180;
  Double_t DetWingLength = 55.0/cos(angle);
  Double_t DetWingWidth  =  16.0;
  Double_t yoffs         = 330.0;

  Double_t ax            =  -50.0 - DetWingLength*cos(angle) - DetWingWidth*sin(angle);
  Double_t ay            =  yoffs + DetWingWidth/2.0 - DetWingLength*sin(angle) - DetWingWidth*sin(angle)*tan(angle);
  Double_t bx            =  -50.0;
  Double_t by            =  yoffs + DetWingWidth/2.0;
  Double_t cx            =  50.0;
  Double_t cy            =  yoffs + DetWingWidth/2.0;
  Double_t dx            =  50.0 + DetWingLength*cos(angle) + DetWingWidth*sin(angle);
  Double_t dy            =  yoffs + DetWingWidth/2.0 - DetWingLength*sin(angle) - DetWingWidth*sin(angle)*tan(angle);
  Double_t ex            =  50.0 + DetWingLength*cos(angle);
  Double_t ey            =  yoffs - DetWingWidth/2.0 - DetWingLength*sin(angle);
  Double_t fx            =  50.0;
  Double_t fy            =  yoffs - DetWingWidth/2.0;
  Double_t gx            =  -50.0;
  Double_t gy            =  yoffs - DetWingWidth/2.0;
  Double_t hx            =  -50.0 - DetWingLength*cos(angle);
  Double_t hy            =  yoffs - DetWingWidth/2.0 - DetWingLength*sin(angle);

  TLine *l1 = new TLine(ax,ay,bx,by);
  TLine *l2 = new TLine(bx,by,cx,cy);
  TLine *l3 = new TLine(cx,cy,dx,dy);
  TLine *l4 = new TLine(dx,dy,ex,ey);
  TLine *l5 = new TLine(ex,ey,fx,fy);
  TLine *l6 = new TLine(fx,fy,gx,gy);
  TLine *l7 = new TLine(gx,gy,hx,hy);
  TLine *l8 = new TLine(hx,hy,ax,ay);

  TCanvas *c1 = new TCanvas("c1","",600,400);

  tree->Draw("Cerenkov.Detector.HitGlobalPositionX:Cerenkov.Detector.HitGlobalPositionY","","", 3000, 0);
  for(int i = 0; i < NEvents; i++) { XPos[i] = tree->GetV1()[i];}
  for(int i = 0; i < NEvents; i++) { YPos[i] = tree->GetV2()[i];}
  TGraph *hDetHits = new TGraph(NEvents,XPos,YPos);
  hDetHits->GetYaxis()->SetRangeUser(280,380);
  hDetHits->SetMarkerColor(kBlue);
  hDetHits->SetMarkerStyle(6);
  hDetHits->GetXaxis()->SetTitle("[cm]");
  hDetHits->GetYaxis()->SetTitle("[cm]");

  tree->Draw("Cerenkov.PMT.PMTLeftNbOfHits","","", 3000, 0);
  TH1F *hPMTLHits = new TH1F(*(TH1F*)htemp->Clone());
  for(int i = 0; i < NEvents; i++) {LeftPMTHits[i] = tree->GetV1()[i];}
  hPMTLHits->GetXaxis()->SetTitle("~ Number of pe's");

  tree->Draw("Cerenkov.PMT.PMTRightNbOfHits","","", 3000, 0);
  TH1F *hPMTRHits = new TH1F(*(TH1F*)htemp->Clone());
  for(int i = 0; i < NEvents; i++) {RightPMTHits[i] = tree->GetV1()[i];}
  hPMTRHits->GetXaxis()->SetTitle("~ Number of pe's");

  tree->Draw("Cerenkov.Detector.PrimaryQ2","","", 3000, 0);
  TH1F *hQ2 = new TH1F(*(TH1F*)htemp->Clone());
  hQ2->SetTitle("Primary Q2");
  hQ2->GetXaxis()->SetTitle("Q2");

  tree->Draw("Cerenkov.Detector.PrimaryQ2","(Cerenkov.PMT.PMTLeftNbOfHits + Cerenkov.PMT.PMTRightNbOfHits)","", 3000, 0);
  TH1F *hWQ2 = new TH1F(*(TH1F*)htemp->Clone());
  hWQ2->SetTitle("PE Weighted Primary Q2");
  hWQ2->GetXaxis()->SetTitle("Q2");
 
  delete c1;

  TProfile *prfl = new TProfile("prfl","Left PMT Hits vs. X",30,-108,108);
  TF1 *fit1 = new TF1("fit1","pol2",-108,108);
  fit1->SetLineColor(kRed);
  fit1->SetLineWidth(2.0);
  for(int i = 0; i < NEvents; i++) {prfl->Fill(XPos[i],LeftPMTHits[i]);}
  prfl->GetYaxis()->SetRangeUser(0,prfl->GetMaximum()+50);
  prfl->GetXaxis()->SetTitle("Electron hit location along x [cm]");
  prfl->GetYaxis()->SetTitle("~ Number of pe's");
  prfl->SetStats(kFALSE);
  
  TProfile *prfr = new TProfile("prfr","Right PMT Hits vs. X",30,-108,108);
  TF1 *fit2 = new TF1("fit2","pol2",-108,108);
  fit2->SetLineColor(kRed);
  fit2->SetLineWidth(2.0);
  for(int i = 0; i < NEvents; i++) {prfr->Fill(XPos[i],RightPMTHits[i]);}
  prfr->GetYaxis()->SetRangeUser(0,prfr->GetMaximum()+50);
  prfr->GetXaxis()->SetTitle("Electron hit location along x [cm]");
  prfr->GetYaxis()->SetTitle("~ Number of pe's");
  prfr->SetStats(kFALSE);

  TProfile *prft = new TProfile("prft","Total PMT Hits vs. X",30,-108,108);
  TF1 *fit3 = new TF1("fit3","pol2",-108,108);
  fit3->SetLineColor(kRed);
  fit3->SetLineWidth(2.0);
  for(int i = 0; i < NEvents; i++) {prft->Fill(XPos[i],RightPMTHits[i]+LeftPMTHits[i]);}
  prft->GetYaxis()->SetRangeUser(0,prft->GetMaximum()+50);
  prft->GetXaxis()->SetTitle("Electron hit location along x [cm]");
  prft->GetYaxis()->SetTitle("~ Total Number of pe's");
  prft->SetStats(kFALSE);

  TCanvas *c1 = new TCanvas("c1","",1200,1000);
  c1->Divide(2,4);
  c1->SetFillColor(kWhite);
  c1->cd(1);
  gPad->SetFillColor(kWhite);
  prfl->Fit(fit1->GetName(),"R+");
  prfl->Draw();
  c1->cd(2);
  gPad->SetFillColor(kWhite);
  prfr->Fit(fit2->GetName(),"R+");
  prfr->Draw();
  c1->cd(3);
  gPad->SetFillColor(kWhite);
  prft->Fit(fit3->GetName(),"R+");
  prft->Draw();
  c1->cd(4);
  gPad->SetFillColor(kWhite);
  hDetHits->Draw("ap");
  l1->Draw();
  l2->Draw();
  l3->Draw();
  l4->Draw();
  l5->Draw();
  l6->Draw();
  l7->Draw();
  l8->Draw();
  c1->Update();
  c1->cd(5);
  gPad->SetFillColor(kWhite);
  hPMTLHits->Draw();
  c1->cd(6);
  gPad->SetFillColor(kWhite);
  hPMTRHits->Draw();
  c1->cd(7);
  gPad->SetFillColor(kWhite);
  hWQ2->Draw();
  c1->cd(8);
  gPad->SetFillColor(kWhite);
  hQ2->Draw();
}
Attachment 6: CerenkovTest3.C  2 kB  | Hide | Hide all
void CerenkovTest3()
{

  const char *const filenames[5] = {
    "QweakSim_03_09_2000_kink.root",
    "QweakSim_03_11_3000_kink_tilt5.8.root",
    "QweakSim_03_09_2000_kink_tilt10.8.root",
    "QweakSim_03_09_2000_kink_tilt20.8.root", 
    "QweakSim_03_11_3000_kink_tilt25.8.root", 
//     "QweakSim_03_09_2000_kink_tilt15.8.root",
  }

  TCanvas *c1 = new TCanvas("c1","",600,400);
  TProfile *prfl[5];
  TProfile *prfr[5];
  
  for(int j = 0; j < 5; j++){
    prfl[j] = new TProfile("prfl","Left PMT Hits vs. X",30,-108,108);
    prfl[j]->SetDirectory(0);
    prfr[j] = new TProfile("prfr","Right PMT Hits vs. X",30,-108,108);
    prfr[j]->SetDirectory(0);
    TFile *f = new TFile(filenames[j],"READ");
    TTree *tree = (TTree *) gROOT->FindObject("Test_NT");
    Int_t NEvents     = tree->GetEntries(); 
    tree->Draw("Cerenkov.Detector.HitLocalPositionX:Cerenkov.PMT.PMTLeftNbOfHits:Cerenkov.PMT.PMTRightNbOfHits","","", 10000, 0);
    for(int i = 0; i < NEvents; i++){
      prfl[j]->Fill(tree->GetV1()[i],tree->GetV2()[i]);
      prfr[j]->Fill(tree->GetV1()[i],tree->GetV3()[i]);
    }
    f->Close();
    delete f;
  }  

  Int_t color[5] = {kBlack,kMagenta,kRed,kGreen,kBlue};
  TCanvas *c2 = new TCanvas("c1","",1200,250);
  gStyle->SetOptStat(kFALSE);
  c2->Divide(2,1);
  c2->SetFillColor(kWhite);
  prfl[0]->SetLineColor(color[0]);
  prfr[0]->SetLineColor(color[0]);
  prfl[0]->GetYaxis()->SetRangeUser(0,200);
  prfr[0]->GetYaxis()->SetRangeUser(0,200);
  c2->cd(1);
  gPad->SetFillColor(kWhite);
  prfl[0]->Draw();
  c2->cd(2);
  gPad->SetFillColor(kWhite);
  prfr[0]->Draw();

  for(int j = 1; j < 5; j++){
    prfl[j]->SetLineColor(color[j]);
    prfr[j]->SetLineColor(color[j]);
    c2->cd(1);
    prfl[j]->Draw("same");
    c2->cd(2);
    prfr[j]->Draw("same");
  }

  const char *const legendentry[5] = {
    "Tilt = 0 degrees",
    "Tilt = 5.8",
    "Tilt = 10.8",
    "Tilt = 20.8",
    "Tilt = 25.8",
  }

  c2->cd(1);
  TLegend *leg1 = new TLegend(0.58,0.82,0.98,0.98);
  for(int j = 0; j < 5; j++) leg1->AddEntry(prfl[j],legendentry[j],"l");
  leg1->SetTextFont(62);
  leg1->SetTextSize(0.030);
  leg1->Draw();
}
ELOG V2.6.0-beta