source: trunk/examples/advanced/Rich/src/AerogelRefData.cc@ 1197

Last change on this file since 1197 was 807, checked in by garnier, 17 years ago

update

File size: 6.7 KB
Line 
1//
2// ********************************************************************
3// * License and Disclaimer *
4// * *
5// * The Geant4 software is copyright of the Copyright Holders of *
6// * the Geant4 Collaboration. It is provided under the terms and *
7// * conditions of the Geant4 Software License, included in the file *
8// * LICENSE and available at http://cern.ch/geant4/license . These *
9// * include a list of copyright holders. *
10// * *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work make any representation or warranty, express or implied, *
14// * regarding this software system or assume any liability for its *
15// * use. Please see the license in the file LICENSE and URL above *
16// * for the full disclaimer and the limitation of liability. *
17// * *
18// * This code implementation is the result of the scientific and *
19// * technical work of the GEANT4 collaboration. *
20// * By using, copying, modifying or distributing the software (or *
21// * any work based on the software) you agree to acknowledge its *
22// * use in resulting scientific publications, and indicate your *
23// * acceptance of all terms of the Geant4 Software license. *
24// ********************************************************************
25//
26// Rich advanced example for Geant4
27// AerogelRefData.cc for Rich of LHCb
28// History:
29// Created: Sajan Easo (Sajan.Easo@cern.ch)
30// Revision and changes: Patricia Mendez (Patricia.Mendez@cern.ch)
31/////////////////////////////////////////////////////////////////////////////
32#include <iostream>
33#include <fstream>
34#include "AerogelRefData.hh"
35#include "RichTbGeometryParameters.hh"
36#include "RichTbMaterialParameters.hh"
37AerogelRefData::AerogelRefData(G4String AgelRefIndFile)
38 :NumberOfRefIndBins(NumAerogelRefIndexPhotonEnergyBins),
39 StdAerogelRefphotE(std::vector<G4double>(NumAerogelRefIndexPhotonEnergyBins)),
40 StdAerogelRefIndexValue(std::vector<G4double>(NumAerogelRefIndexPhotonEnergyBins)),
41 AerogelRefIndShift(std::vector<G4double>( MaxNumberOfAerogelTypes)),
42 AerogelWavelengthRef(std::vector<G4double>( MaxNumberOfAerogelTypes)) {
43
44 AerogelRefIndexFileName=AgelRefIndFile;
45
46 StdAerogelNominalRefractiveIndex= StdAerogelNominalRefIndex;
47
48 ReadStdAerogelRefIndex();
49 G4double AerogelRefEnergy,AerogelReferenceRefInd;
50 for(G4int ityp=0; ityp< MaxNumberOfAerogelTypes; ityp++ ){
51 AerogelWavelengthRef[ityp]=AerogelReferenceRefIndWavelength[ityp];
52 AerogelRefEnergy=
53 (PhotMomWaveConv/(AerogelWavelengthRef[ityp]/nanometer))*eV;
54
55 AerogelReferenceRefInd = GetStdAgelRefIndexE(AerogelRefEnergy);
56 AerogelRefIndShift[ityp] = AerogelReferenceRefInd -
57 StdAerogelNominalRefractiveIndex ;
58 }
59
60}
61G4double AerogelRefData::GetStdAgelRefIndexE(G4double StdAgelRefEnergy) {
62
63 G4double refint= StdAerogelNominalRefractiveIndex;
64 if(StdAgelRefEnergy <= StdAerogelRefphotE[0] ) {
65 refint= StdAerogelRefphotE[0];
66 G4cout<<"Ref index reference value outside bounds "<<G4endl;
67
68 }
69 if(StdAgelRefEnergy > StdAerogelRefphotE[ NumberOfRefIndBins-1 ] ) {
70 refint= StdAerogelRefphotE[ NumberOfRefIndBins-1 ];
71 G4cout<<"Ref index reference value outside bounds "<<G4endl;
72
73 }
74
75 for(G4int iabin=0; iabin < NumberOfRefIndBins-1 ; iabin++ ){
76
77 G4double E1= StdAerogelRefphotE[iabin];
78 G4double E2= StdAerogelRefphotE[iabin+1];
79
80 if(StdAgelRefEnergy > E1 && StdAgelRefEnergy <= E2 && E1 != E2 ){
81
82 //Now interpolate
83 G4double r1= StdAerogelRefIndexValue[iabin];
84 G4double r2= StdAerogelRefIndexValue[iabin+1];
85 G4double slp=(r2-r1)/(E2-E1);
86 G4double ain= r1-slp*E1;
87 refint= slp*StdAgelRefEnergy + ain;
88 break;
89 }
90
91 }
92
93 return refint;
94
95}
96void AerogelRefData::ReadStdAerogelRefIndex() {
97
98 const char* StdAerogelRefIndFile= AerogelRefIndexFileName.c_str();
99 if( StdAerogelRefIndFile == 0 ){
100 G4cout<<" Aerogel RefIndex Input file missing"
101 << " Please provide the correct file name "
102 <<G4endl;
103 }else {
104
105 G4cout<<"Now reading Std Aerogel Ref Index from "
106 <<StdAerogelRefIndFile<<G4endl;
107 std::ifstream finpga(StdAerogelRefIndFile);
108 G4double ephot,rind;
109 for (G4int fbin=0; fbin< NumberOfRefIndBins ; fbin++){
110 finpga>>ephot;
111 finpga>>rind;
112 StdAerogelRefphotE[fbin]=ephot*eV;
113 StdAerogelRefIndexValue[fbin]=rind;
114 //for test avoid the variation in ref index. SE 18-10-01
115 // StdAerogelRefIndexValue[fbin]= StdAerogelNominalRefractiveIndex;
116
117 }
118
119 }
120}
121
122std::vector<G4double> AerogelRefData::GetCurAerogelRefIndValueVect (G4int AerogelTnum ) {
123 // no test for tdr rich. now put back to TB conditions.
124 std::vector<G4double>CAgel(NumberOfRefIndBins);
125 // AerogelType Ctype;
126 G4double CRn = GetRefnominal( AerogelTnum );
127 for (G4int ib=0; ib < NumberOfRefIndBins; ib++ ){
128
129 CAgel[ib]= StdAerogelRefIndexValue[ib]
130 - StdAerogelNominalRefractiveIndex-AerogelRefIndShift[ AerogelTnum] + CRn;
131
132 }
133 return CAgel;
134}
135
136G4double AerogelRefData:: GetCurAerogelRefIndValue( G4int rbinw , G4int AerogelTnum ){
137 // no test for tdr. now in TB conditions.
138 G4double crn= GetRefnominal( AerogelTnum );
139 return StdAerogelRefIndexValue[rbinw]
140 - StdAerogelNominalRefractiveIndex-AerogelRefIndShift[ AerogelTnum] +crn;
141 // - StdAerogelNominalRefractiveIndex-AerogelRefIndShift[ AerogelTnum] +crn;
142 // return StdAerogelRefIndexValue[rbinw]-1.03269+crn;
143}
144
145G4double AerogelRefData::GetRefnominal(G4int AgelTnum ) {
146
147 G4double CRnominal = StdAerogelNominalRefractiveIndex;
148 if(AgelTnum == 0 ){
149
150 CRnominal= AerogelTypeANominalRefIndex;
151 }else if (AgelTnum == 1 ) {
152
153 CRnominal= AerogelTypeBNominalRefIndex;
154 }else if (AgelTnum == 2 ) {
155
156 CRnominal= AerogelTypeCNominalRefIndex;
157 }else if (AgelTnum == 3 ) {
158
159 CRnominal= AerogelTypeDNominalRefIndex;
160 }else if (AgelTnum == 4 ) {
161
162 CRnominal= AerogelTypeENominalRefIndex;
163
164 } else {G4cout<<"Unknown Aerogel Type "<<G4endl; }
165
166 return CRnominal;
167}
168
169AerogelType AerogelRefData::GetAerogelType(G4int AgelTnum ) {
170 AerogelType Ctype=AerogelTypeA;
171
172 if(AgelTnum == 0 ){
173 Ctype=AerogelTypeA;
174
175 }else if (AgelTnum == 1 ) {
176 Ctype=AerogelTypeB;
177
178 }else if (AgelTnum == 2 ) {
179 Ctype=AerogelTypeC;
180
181 }else if (AgelTnum == 3 ) {
182 Ctype=AerogelTypeD;
183
184 }else if (AgelTnum == 4 ) {
185 Ctype=AerogelTypeE;
186
187
188 } else {G4cout<<"Unknown Aerogel Type "<<G4endl; }
189
190 return Ctype;
191}
192
193AerogelRefData::~AerogelRefData(){ ; }
194
195
196
197
198
199
200
201
202
Note: See TracBrowser for help on using the repository browser.