RICH/Opticals.h

   1 #ifdef __CINT__
   2 void Opticals()
   3 {
   4   gROOT->Reset();
   5 #endif
   6   int i;
   7   const Int_t kNbins=30;//number of photon energy points
   8
   9   Float_t aPckov[kNbins];  for(i=0;i<kNbins;i++) aPckov[i]=1e-9*AliRICHParam::PhotonEnergy(i);
  10
  11
  12   Float_t aIdxC6F14[kNbins];  //Freon ref index
  13   Float_t aIdxCH4[kNbins];    //Methane ref index
  14   Float_t aIdxSiO2[kNbins];   //Quartz ref index
  15   Float_t aIdxMetal[kNbins];  //Metal ref index always 0
  16   Float_t aIdxGel[kNbins];    //Aerogel ref index
  17
  18   for (i=0;i<kNbins;i++){
  19     aIdxC6F14[i] = AliRICHParam::RefIdxC6F14(aPckov[i]*1e9);
  20     aIdxSiO2[i]  = AliRICHParam::RefIdxSiO2(aPckov[i]*1e9);
  21     aIdxCH4[i]   = AliRICHParam::RefIdxCH4(5.5);
  22     aIdxGel[i]   = AliRICHParam::RefIdxGel(5.5);
  23     aIdxMetal[i]  =0;//metal
  24   }
  25
  26   Float_t aAbsC6F14[kNbins]={//New values from A.DiMauro 28.10.03 total 30
  27     32701.4219, 17996.1141, 10039.7281, 1799.1230, 1799.1231, 1799.1231, 1241.4091, 179.0987, 179.0986, 179.0987,
  28       179.0987,   118.9800,    39.5058,   23.7244,   11.1283,    7.1573,    3.6249,   2.1236,   0.7362,   0.5348,
  29         0.3387,     0.3074,     0.3050,    0.0001,    0.0001,    0.0001,    0.0001,   0.0001,   0.0001,   0.0001};
  30   Float_t aAbsSiO2[kNbins]={//New values from A.DiMauro 28.10.03 total 30
  31        34.4338, 30.5424, 30.2584, 31.4928, 31.7868, 17.8397, 9.3410, 6.4492, 6.1128, 5.8128,
  32         5.5589,  5.2877,  5.0162,  4.7999,  4.5734,  4.2135, 3.7471, 2.6033, 1.5223, 0.9658,
  33         0.4242,  0.2500,  0.1426,  0.0863,  0.0793,  0.0724, 0.0655, 0.0587, 0.0001, 0.0001};
  34
  35   Float_t aAbsOpSiO2[kNbins];
  36   Float_t aAbsCH4[kNbins];
  37   Float_t aAbsGrid[kNbins];
  38   Float_t aAbsCsI[kNbins];
  39   Float_t aAbsGel[kNbins];
  40   Float_t aAbsRef[kNbins];//abs prob on reflector
  41
  42     for (i=0;i<kNbins;i++){
  43       aAbsCH4[i]         =AliRICHParam::AbsCH4(AliRICHParam::PhotonEnergy(i));
  44       aAbsOpSiO2[i]      =1e-5;
  45       aAbsCsI[i]         =1e-4;
  46       aAbsGrid[i]        =1e-4;
  47       aAbsGel[i]         =AliRICHParam::AbsGel(AliRICHParam::PhotonEnergy(i));
  48       aAbsRef[i]         =0.01;//1% reflector abs prob
  49     }
  50
  51   Float_t aQeCsI[kNbins] = {//New values from A.DiMauro 28.10.03 total 31
  52                             0.0002, 0.0006, 0.0007, 0.0010, 0.0049, 0.0073, 0.0104, 0.0519, 0.0936, 0.1299,
  53                             0.1560, 0.1768, 0.1872, 0.1976, 0.2142, 0.2288, 0.2434, 0.2599, 0.2673, 0.2808,
  54                             0.2859, 0.2954, 0.3016, 0.3120, 0.3172, 0.3224, 0.3266, 0.3328, 0.3359, 0.3390};
  55 //                             0.3431};
  56
  57
  58   Float_t aQeAll[kNbins];
  59   for(i=0;i<kNbins;i++){
  60     aQeCsI[i]/= (1.0-Fresnel(AliRICHParam::PhotonEnergy(i),1.0,0)); //FRESNEL LOSS CORRECTION
  61     aQeAll[i]=1; //QE for all other materials except for PC must be 1.
  62   }
  63
  64
  65 #ifdef __CINT__
  66
  67   Float_t aAbsC6F14old[kNbins]={//previous values 26 in total added 0.0001 to the 30
  68       179.0987, 179.0987, 179.0987, 179.0987, 179.0987, 179.0987, 179.0987, 179.0987, 179.0987, 142.9206,
  69        56.6496,  25.5862,  13.9529,  12.0391,  10.4295,   8.8042,   7.0690,   4.4611,   2.0284,   1.2930,
  70         0.5772,   0.4074,   0.3350,   0.0001,   0.0001,   0.0001,   0.0001,   0.0001,   0.0001,   0.0001};
  71   Float_t aAbsSiO2old[kNbins]={//previous values 26 in total added 0.0001 to the 30
  72       105.8000,  65.5200,  48.5800,  42.8500,  35.7900,  31.2620,  28.5980,  27.5270,  25.0070,  22.8150,
  73        21.0040,  19.2660,  17.5250,  15.8780,  14.1770,  11.7190,   9.2820,   6.6200,   4.0930,   2.6010,
  74         1.1490,   0.6670,   0.3630,   0.1920,   0.1500,   0.1090,   0.0001,   0.0001,   0.0001,   0.0001};
  75   Float_t aQeCsIold[kNbins]={
  76     0.0002, 0.0006, 0.0007, 0.0050, 0.0075, 0.0101, 0.0243, 0.0405, 0.0689, 0.1053,
  77     0.1215, 0.1417, 0.1579, 0.1620, 0.1661, 0.1677, 0.1743, 0.1768, 0.1793, 0.1826,
  78     0.1859, 0.1876, 0.1892, 0.1909, 0.2075, 0.2158, 0.0001, 0.0001, 0.0001, 0.0001 };
  79
  80 //Now plot all the thigs
  81 //Markers and colors for different curves:
  82
  83   const Double_t kWidth=0.25,kHeight=0.2;
  84   const Int_t kC6F14Marker=24 , kC6F14Color=kRed;
  85   const Int_t kCH4Marker  =25 , kCH4Color  =kGreen;
  86   const Int_t kSiO2M      =26 , kSiO2Color =kBlue;
  87   const Int_t kCsIMarker  = 2 , kCsIColor  =kMagenta;
  88   const Int_t kGelMarker  =27 , kGelColor  =46;
  89   const Int_t kRefMarker  =28 , kRefColor  =47;
  90
  91   TCanvas *pC=new TCanvas("c1","RICH optics to check",1100,900);
  92
  93   pC->Divide(3,2);
  94 //Ref index
  95   pC->cd(1);
  96   TGraph *pIdxC6F14=new TGraph(kNbins,aPckov,aIdxC6F14);pIdxC6F14->SetMarkerStyle(kC6F14Marker);pIdxC6F14->SetMarkerColor(kC6F14Color);
  97   TGraph *pIdxSiO2 =new TGraph(kNbins,aPckov,aIdxSiO2); pIdxSiO2 ->SetMarkerStyle(kSiO2M)      ;pIdxSiO2 ->SetMarkerColor(kSiO2Color);
  98   TGraph *pIdxCH4  =new TGraph(kNbins,aPckov,aIdxCH4);  pIdxCH4  ->SetMarkerStyle(kCH4Marker)       ;pIdxCH4  ->SetMarkerColor(kCH4Color);
  99   TGraph *pIdxGel  =new TGraph(kNbins,aPckov,aIdxGel);  pIdxGel  ->SetMarkerStyle(kGelMarker)  ;pIdxGel  ->SetMarkerColor(kGelColor);
 100
 101   TMultiGraph *pIdxMG=new TMultiGraph();  TLegend *pIdxLe=new TLegend(0.5,0.21,0.5+kWidth,0.21+kHeight);
 102   pIdxMG->Add(pIdxC6F14); pIdxLe->AddEntry(pIdxC6F14,"C6F14"  ,"p");
 103   pIdxMG->Add(pIdxSiO2) ; pIdxLe->AddEntry(pIdxSiO2 ,"SiO2"   ,"p");
 104   pIdxMG->Add(pIdxCH4)  ; pIdxLe->AddEntry(pIdxCH4  ,"CH4"    ,"p");
 105   pIdxMG->Add(pIdxGel)  ; pIdxLe->AddEntry(pIdxGel  ,"Aerogel","p");
 106   pIdxMG->Draw("APL");
 107   pIdxMG->SetTitle("Refractive index");  pIdxMG->GetXaxis()->SetTitle("energy, GeV");
 108   pIdxMG->Draw("APL");
 109   pIdxLe->Draw();
 110 //Absorbtion
 111   pC->cd(2);
 112   gPad->SetLogy();
 113   TGraph *pAbsC6F14=new TGraph(kNbins,aPckov,aAbsC6F14);pAbsC6F14->SetMarkerStyle(kC6F14Marker); pAbsC6F14->SetMarkerColor(kC6F14Color);
 114   TGraph *pAbsSiO2 =new TGraph(kNbins,aPckov,aAbsSiO2) ;pAbsSiO2 ->SetMarkerStyle(kSiO2M)      ; pAbsSiO2 ->SetMarkerColor(kSiO2Color);
 115   TGraph *pAbsCH4  =new TGraph(kNbins,aPckov,aAbsCH4)  ;pAbsCH4  ->SetMarkerStyle(kCH4Marker)  ; pAbsCH4  ->SetMarkerColor(kCH4Color);
 116   TGraph *pAbsGel  =new TGraph(kNbins,aPckov,aAbsGel)  ;pAbsGel  ->SetMarkerStyle(kGelMarker)  ; pAbsGel  ->SetMarkerColor(kGelColor);
 117   TGraph *pAbsRef  =new TGraph(kNbins,aPckov,aAbsRef)  ;pAbsRef  ->SetMarkerStyle(kRefMarker)  ; pAbsRef  ->SetMarkerColor(kRefColor);
 118
 119   TMultiGraph *pAbsMG=new TMultiGraph();  TLegend *pAbsLe=new TLegend(0.2,0.15,0.2+kWidth,0.15+kHeight);
 120   pAbsMG->Add(pAbsC6F14);      pAbsLe->AddEntry(pAbsC6F14,  "C6F14"    ,"p");
 121   pAbsMG->Add(pAbsSiO2);       pAbsLe->AddEntry(pAbsSiO2 ,  "SiO2"     ,"p");
 122   pAbsMG->Add(pAbsCH4);        pAbsLe->AddEntry(pAbsCH4  ,  "CH4"      ,"p");
 123   pAbsMG->Add(pAbsGel);        pAbsLe->AddEntry(pAbsGel  ,  "Aerogel"  ,"p");
 124   pAbsMG->Add(pAbsRef);        pAbsLe->AddEntry(pAbsRef  ,  "Reflector","p");
 125   pAbsMG->Draw("APL");
 126   pAbsMG->SetTitle("Absorbtion length,cm"); pAbsMG->GetXaxis()->SetTitle("energy, GeV");
 127   pAbsMG->Draw("APL");
 128   pAbsLe->Draw();
 129 //QE
 130   pC->cd(4);
 131   TGraph *pQeCsI=new TGraph(kNbins,aPckov,aQeCsI); pQeCsI->SetMarkerStyle(kCsIMarker); pQeCsI->SetMarkerColor(kCsIColor);
 132   pQeCsI->Draw("APL");
 133   pQeCsI->SetTitle("QE");  pQeCsI->GetXaxis()->SetTitle("energy, GeV");
 134   pQeCsI->Draw("APL");
 135 //transmission
 136   pC->cd(5);
 137   Float_t mm =0.1;Float_t cm=1.;
 138   Float_t aTrC6F14[kNbins],aTrSiO2[kNbins],aTrCH4[kNbins];
 139   Float_t aTotTr[kNbins];
 140   for(Int_t i=0;i<kNbins;i++){//calculate probability for photon to survive during transversing a volume of material with absorption length
 141     aTrC6F14[i]=TMath::Exp(-15*mm/(aAbsC6F14[i]+0.0001));
 142     aTrSiO2[i] =TMath::Exp(-5*mm/(aAbsSiO2[i] +0.0001));
 143     aTrCH4[i]  =TMath::Exp(-8*cm/(aAbsCH4[i]  +0.0001));
 144     aTotTr[i]    =aTrC6F14[i]*aTrSiO2[i]*aTrCH4[i]*aQeCsI[i];
 145   }
 146   TGraph *pTrC6F14=new TGraph(kNbins,aPckov,aTrC6F14);pTrC6F14->SetMarkerStyle(kC6F14Marker);pTrC6F14->SetMarkerColor(kC6F14Color);
 147   TGraph *pTrSiO2 =new TGraph(kNbins,aPckov,aTrSiO2) ;pTrSiO2 ->SetMarkerStyle(kSiO2M)      ;pTrSiO2 ->SetMarkerColor(kSiO2Color);
 148   TGraph *pTrCH4  =new TGraph(kNbins,aPckov,aTrCH4)  ;pTrCH4  ->SetMarkerStyle(kCH4Marker)  ;pTrCH4  ->SetMarkerColor(kCH4Color);
 149   TGraph *pTotTr  =new TGraph(kNbins,aPckov,aTotTr)  ;pTotTr  ->SetMarkerStyle(30)          ;pTotTr  ->SetMarkerColor(kYellow);
 150
 151   TMultiGraph *pTrMG=new TMultiGraph();  TLegend *pTrLe=new TLegend(0.2,0.4,0.2+kWidth,0.4+kHeight);
 152   pTrMG->Add(pQeCsI);       pTrLe->AddEntry(pQeCsI,   "CsI QE", "p");
 153   pTrMG->Add(pTrC6F14);     pTrLe->AddEntry(pTrC6F14, "C6F14", "p");
 154   pTrMG->Add(pTrSiO2);      pTrLe->AddEntry(pTrSiO2,  "SiO2",  "p");
 155   pTrMG->Add(pTrCH4);       pTrLe->AddEntry(pTrCH4,   "CH4",   "p");
 156   pTrMG->Add(pTotTr);       pTrLe->AddEntry(pTotTr,   "total", "p");
 157   pTrMG->Draw("APL");
 158   pTrMG->SetTitle("Transmission");  pTrMG->GetXaxis()->SetTitle("energy, GeV");
 159   pTrMG->Draw("APL");
 160   pTrLe->Draw();
 161 //comparison of new and old staff
 162   pC->cd(6);
 163   TGraph *pQeCsIold=new TGraph(kNbins,aPckov,aQeCsIold); pQeCsIold->SetMarkerStyle(kC6F14Marker);pQeCsIold->SetMarkerColor(kC6F14Color);
 164
 165   TMultiGraph *pCompMG=new TMultiGraph;  TLegend *pCompLe=new TLegend(0.2,0.6,0.2+kWidth,0.6+kHeight);
 166   pCompMG->Add(pQeCsI);       pCompLe->AddEntry(pQeCsI,    "QE new 30.10.03", "p");
 167   pCompMG->Add(pQeCsIold);    pCompLe->AddEntry(pQeCsIold, "QE old 01.01.02", "p");
 168   pCompMG->Draw("APL");
 169   pCompMG->SetTitle("Comparison of new and old staff");
 170   pCompMG->GetXaxis()->SetTitle("energy, GeV");
 171   pCompMG->Draw("APL");
 172   pCompLe->Draw();
 173
 174   return;
 175   TCanvas *pQeC=new TCanvas("pQeC","CsI QE currently all the same",800,900);
 176   pQeC->Divide(2,4);
 177   for(int i=1;i<=7;i++){
 178     pQeC->cd(i);
 179     switch(i){
 180       case 1: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI);pQeCsIGr->SetTitle("Module 1");break;
 181       case 2: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI);pQeCsIGr->SetTitle("Module 2");break;
 182       case 3: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI);pQeCsIGr->SetTitle("Module 3");break;
 183       case 4: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI);pQeCsIGr->SetTitle("Module 4");break;
 184       case 5: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI);pQeCsIGr->SetTitle("Module 5");break;
 185       case 6: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI);pQeCsIGr->SetTitle("Module 6");break;
 186       case 7: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI);pQeCsIGr->SetTitle("Module 7");break;
 187     }
 188     pQeCsIGr->SetMarkerStyle(kCsIMarker); pQeCsIGr->SetMarkerColor(kCsIColor);
 189     pQeCsIGr->Draw("APL");
 190     pQeCsIGr->GetXaxis()->SetTitle("energy, GeV");
 191     pQeCsIGr->Draw("APL");
 192   }
 193 }//main
 194 //__________________________________________________________________________________________________
 195 Float_t Fresnel(Float_t ene,Float_t pdoti, Bool_t pola)
 196 {//ENE(EV), PDOTI=COS(INC.ANG.), PDOTR=COS(POL.PLANE ROT.ANG.)
 197
 198     Float_t en[36] = {5.0,5.1,5.2,5.3,5.4,5.5,5.6,5.7,5.8,5.9,6.0,6.1,6.2,
 199                       6.3,6.4,6.5,6.6,6.7,6.8,6.9,7.0,7.1,7.2,7.3,7.4,7.5,7.6,7.7,
 200                       7.8,7.9,8.0,8.1,8.2,8.3,8.4,8.5};
 201
 202
 203     Float_t csin[36] = {2.14,2.21,2.33,2.48,2.76,2.97,2.99,2.59,2.81,3.05,
 204                         2.86,2.53,2.55,2.66,2.79,2.96,3.18,3.05,2.84,2.81,2.38,2.11,
 205                         2.01,2.13,2.39,2.73,3.08,3.15,2.95,2.73,2.56,2.41,2.12,1.95,
 206                         1.72,1.53};
 207
 208     Float_t csik[36] = {0.,0.,0.,0.,0.,0.196,0.408,0.208,0.118,0.49,0.784,0.543,
 209                         0.424,0.404,0.371,0.514,0.922,1.102,1.139,1.376,1.461,1.253,0.878,
 210                         0.69,0.612,0.649,0.824,1.347,1.571,1.678,1.763,1.857,1.824,1.824,
 211                         1.714,1.498};
 212     Float_t xe=ene;
 213     Int_t  j=Int_t(xe*10)-49;
 214     Float_t cn=csin[j]+((csin[j+1]-csin[j])/0.1)*(xe-en[j]);
 215     Float_t ck=csik[j]+((csik[j+1]-csik[j])/0.1)*(xe-en[j]);
 216
 217     //FORMULAE FROM HANDBOOK OF OPTICS, 33.23 OR
 218     //W.R. HUNTER, J.O.S.A. 54 (1964),15 , J.O.S.A. 55(1965),1197
 219
 220     Float_t sinin=TMath::Sqrt(1-pdoti*pdoti);
 221     Float_t tanin=sinin/pdoti;
 222
 223     Float_t c1=cn*cn-ck*ck-sinin*sinin;
 224     Float_t c2=4*cn*cn*ck*ck;
 225     Float_t aO=TMath::Sqrt(0.5*(TMath::Sqrt(c1*c1+c2)+c1));
 226     Float_t b2=0.5*(TMath::Sqrt(c1*c1+c2)-c1);
 227
 228     Float_t rs=((aO-pdoti)*(aO-pdoti)+b2)/((aO+pdoti)*(aO+pdoti)+b2);
 229     Float_t rp=rs*((aO-sinin*tanin)*(aO-sinin*tanin)+b2)/((aO+sinin*tanin)*(aO+sinin*tanin)+b2);
 230
 231
 232     //CORRECTION FACTOR FOR SURFACE ROUGHNESS
 233     //B.J. STAGG  APPLIED OPTICS, 30(1991),4113
 234
 235     Float_t sigraf=18.;
 236     Float_t lamb=1240/ene;
 237     Float_t fresn;
 238
 239     Float_t  rO=TMath::Exp(-(4*TMath::Pi()*pdoti*sigraf/lamb)*(4*TMath::Pi()*pdoti*sigraf/lamb));
 240
 241     if(pola)
 242     {
 243         Float_t pdotr=0.8;                                 //DEGREE OF POLARIZATION : 1->P , -1->S
 244         fresn=0.5*(rp*(1+pdotr)+rs*(1-pdotr));
 245     }
 246     else
 247         fresn=0.5*(rp+rs);
 248
 249     fresn = fresn*rO;
 250     return(fresn);
 251 }//Fresnel()
 252 #endif