| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 << 25 //
>> 26 // $Id$
>> 27 //
>> 28 //
26 ////////////////////////////////////////////// 29 ////////////////////////////////////////////////////////////////////////
27 // Optical Surface Class Implementation 30 // Optical Surface Class Implementation
28 ////////////////////////////////////////////// 31 ////////////////////////////////////////////////////////////////////////
29 // 32 //
30 // File: G4OpticalSurface.cc 33 // File: G4OpticalSurface.cc
31 // Description: An optical surface class for u 34 // Description: An optical surface class for use in G4OpBoundaryProcess
32 // Version: 2.0 35 // Version: 2.0
33 // Created: 1997-06-26 36 // Created: 1997-06-26
34 // Author: Peter Gumplinger 37 // Author: Peter Gumplinger
35 // updated: 2017-02-24 Mariele Stockhoff a << 38 // updated: 2017-02-24 Mariele Stockhoff add DAVIS model
>> 39 // mail: gum@triumf.ca
36 ////////////////////////////////////////////// 40 ////////////////////////////////////////////////////////////////////////
37 41
38 #include "G4OpticalSurface.hh" << 42 #include <iostream>
>> 43 #include <fstream>
39 44
>> 45 //#include "G4ios.hh"
40 #include "globals.hh" 46 #include "globals.hh"
>> 47 #include "G4OpticalSurface.hh"
41 48
42 #include <zlib.h> << 49 /////////////////////////
43 << 50 // Class Implementation
44 #include <fstream> << 51 /////////////////////////
45 #include <iostream> << 52
>> 53 //////////////
>> 54 // Operators
>> 55 //////////////
46 56
47 G4OpticalSurface& G4OpticalSurface::operator=( 57 G4OpticalSurface& G4OpticalSurface::operator=(const G4OpticalSurface& right)
48 { 58 {
49 if (this != &right) { << 59 if (this != &right)
50 theName = right.theName; << 60 {
51 theType = right.theType; << 61 theName = right.theName;
52 theModel = right.theModel; << 62 theType = right.theType;
53 theFinish = right.theFinish; << 63 theModel = right.theModel;
54 sigma_alpha = right.sigma_alpha; << 64 theFinish = right.theFinish;
55 polish = right.polish; << 65 sigma_alpha = right.sigma_alpha;
56 theMaterialPropertiesTable = right.theMate << 66 polish = right.polish;
57 << 67 theMaterialPropertiesTable = right.theMaterialPropertiesTable;
58 delete[] AngularDistribution; << 68
59 AngularDistribution = new G4float[incident << 69 if (AngularDistribution) delete [] AngularDistribution;
60 *(AngularDistribution) = *(right.AngularDi << 70 AngularDistribution =
61 << 71 new G4float[incidentIndexMax*thetaIndexMax*phiIndexMax];
62 delete[] AngularDistributionLUT; << 72 *(AngularDistribution) = *(right.AngularDistribution);
63 AngularDistributionLUT = new G4float[index << 73
64 *(AngularDistributionLUT) = *(right.Angula << 74 if (AngularDistributionLUT) delete [] AngularDistributionLUT;
65 << 75 AngularDistributionLUT = new G4float[indexmax];
66 delete[] Reflectivity; << 76 *(AngularDistributionLUT) = *(right.AngularDistributionLUT);
67 Reflectivity = new G4float[RefMax]; << 77
68 *(Reflectivity) = *(right.Reflectivity); << 78 if (Reflectivity) delete [] Reflectivity;
69 << 79 Reflectivity = new G4float[RefMax];
70 delete DichroicVector; << 80 *(Reflectivity) = *(right.Reflectivity);
71 DichroicVector = new G4Physics2DVector(); << 81
72 *DichroicVector = *(right.DichroicVector); << 82 if (DichroicVector) delete DichroicVector;
73 } << 83 DichroicVector = new G4Physics2DVector();
>> 84 *DichroicVector = *(right.DichroicVector);
>> 85 }
74 return *this; 86 return *this;
75 } 87 }
76 88
77 G4OpticalSurface::G4OpticalSurface(const G4Str << 89 /////////////////
78 G4OpticalSurfaceFinish finish, G4SurfaceType << 90 // Constructors
79 : G4SurfaceProperty(name, type), theModel(mo << 91 /////////////////
80 { << 92
81 AngularDistribution = nullptr; << 93 G4OpticalSurface::G4OpticalSurface(const G4String& name,
82 << 94 G4OpticalSurfaceModel model,
83 AngularDistributionLUT = nullptr; << 95 G4OpticalSurfaceFinish finish,
84 Reflectivity = nullptr; << 96 G4SurfaceType type,
85 << 97 G4double value)
86 DichroicVector = nullptr; << 98 : G4SurfaceProperty(name,type),
87 << 99 theModel(model),
88 switch (theModel) { << 100 theFinish(finish),
89 case glisur: << 101 theMaterialPropertiesTable(0)
90 polish = value; << 102 {
91 sigma_alpha = 0.0; << 103 if ( model == glisur ){
92 break; << 104 polish = value;
93 case LUT: << 105 sigma_alpha = 0.0;
94 case dichroic: << 106 }
95 case DAVIS: << 107 else if ( model == unified ) {
96 ReadDataFile(); << 108 sigma_alpha = value;
97 // fall through << 109 polish = 0.0;
98 case unified: << 110 }
99 sigma_alpha = value; << 111 else if ( model == LUT ) {
100 polish = 0.0; << 112 sigma_alpha = value;
101 break; << 113 polish = 0.0;
102 default: << 114 }
103 G4Exception("G4OpticalSurface::G4Optical << 115 else if ( model == dichroic ) {
104 "Constructor called with INVALID model << 116 sigma_alpha = value;
105 } << 117 polish = 0.0;
>> 118 }
>> 119 else if ( model == DAVIS ) {
>> 120 sigma_alpha = value;
>> 121 polish = 0.0;
>> 122 }
>> 123 else {
>> 124 G4Exception("G4OpticalSurface::G4OpticalSurface()", "mat309",
>> 125 FatalException,
>> 126 "Constructor called with INVALID model.");
>> 127 }
>> 128
>> 129 AngularDistribution = NULL;
>> 130
>> 131 AngularDistributionLUT = NULL;
>> 132 Reflectivity = NULL;
>> 133
>> 134 DichroicVector = NULL;
>> 135
>> 136 if (type == dielectric_LUT) {
>> 137 AngularDistribution =
>> 138 new G4float[incidentIndexMax*thetaIndexMax*phiIndexMax];
>> 139 ReadLUTFile();
>> 140 }
>> 141
>> 142 if (type == dielectric_LUTDAVIS) {
>> 143 AngularDistributionLUT = new G4float[indexmax];
>> 144 ReadLUTDAVISFile();
>> 145
>> 146 Reflectivity = new G4float[RefMax];
>> 147 ReadReflectivityLUTFile();
>> 148 }
>> 149
>> 150 if (type == dielectric_dichroic) {
>> 151 DichroicVector = new G4Physics2DVector();
>> 152 ReadDichroicFile();
>> 153 }
106 } 154 }
107 155
108 G4OpticalSurface::~G4OpticalSurface() 156 G4OpticalSurface::~G4OpticalSurface()
109 { 157 {
110 delete[] AngularDistribution; << 158 if (AngularDistribution) delete [] AngularDistribution;
111 159
112 delete[] AngularDistributionLUT; << 160 if (AngularDistributionLUT) delete [] AngularDistributionLUT;
>> 161 if (Reflectivity) delete [] Reflectivity;
113 162
114 delete[] Reflectivity; << 163 if (DichroicVector) delete DichroicVector;
115 <<
116 delete DichroicVector; <<
117 } 164 }
118 165
119 G4OpticalSurface::G4OpticalSurface(const G4Opt << 166 G4OpticalSurface::G4OpticalSurface(const G4OpticalSurface &right)
120 : G4SurfaceProperty(right.theName, right.the << 167 : G4SurfaceProperty(right.theName,right.theType)
121 { 168 {
122 *this = right; << 169 *this = right;
123 this->theName = right.theName; << 170 this->theName = right.theName;
124 this->theType = right.theType; << 171 this->theType = right.theType;
125 this->theModel = right.theModel; << 172 this->theModel = right.theModel;
126 this->theFinish = right.theFinish; << 173 this->theFinish = right.theFinish;
127 this->sigma_alpha = right.sigma_alpha; << 174 this->sigma_alpha = right.sigma_alpha;
128 this->polish = right.polish; << 175 this->polish = right.polish;
129 this->theMaterialPropertiesTable = right.the << 176 this->theMaterialPropertiesTable = right.theMaterialPropertiesTable;
130 << 177
131 delete[] AngularDistribution; << 178 if (this->AngularDistribution) delete [] AngularDistribution;
132 this->AngularDistribution = new G4float[inci << 179 this->AngularDistribution =
133 *(this->AngularDistribution) = *(right.Angul << 180 new G4float[incidentIndexMax*thetaIndexMax*phiIndexMax];
134 << 181 *(this->AngularDistribution) = *(right.AngularDistribution);
135 delete[] AngularDistributionLUT; << 182
136 this->AngularDistributionLUT = new G4float[i << 183 if (this->AngularDistributionLUT) delete [] AngularDistributionLUT;
137 *(this->AngularDistributionLUT) = *(right.An << 184 this->AngularDistributionLUT = new G4float[indexmax];
>> 185 *(this->AngularDistributionLUT) = *(right.AngularDistributionLUT);
>> 186 if (this->Reflectivity) delete [] Reflectivity;
>> 187 this->Reflectivity = new G4float[RefMax];
>> 188 *(this->Reflectivity) = *(right.Reflectivity);
138 189
139 delete[] Reflectivity; << 190 if (this->DichroicVector) delete DichroicVector;
140 this->Reflectivity = new G4float[RefMax]; << 191 this->DichroicVector = new G4Physics2DVector();
141 *(this->Reflectivity) = *(right.Reflectivity << 192 *(this->DichroicVector) = *(right.DichroicVector);
142 <<
143 delete DichroicVector; <<
144 this->DichroicVector = new G4Physics2DVector <<
145 *(this->DichroicVector) = *(right.DichroicVe <<
146 } 193 }
147 194
148 G4bool G4OpticalSurface::operator==(const G4Op << 195 G4int G4OpticalSurface::operator==(const G4OpticalSurface &right) const
149 { 196 {
150 return (this == (G4OpticalSurface*)&right); << 197 return (this == (G4OpticalSurface *) &right);
151 } 198 }
152 199
153 G4bool G4OpticalSurface::operator!=(const G4Op << 200 G4int G4OpticalSurface::operator!=(const G4OpticalSurface &right) const
154 { 201 {
155 return (this != (G4OpticalSurface*)&right); << 202 return (this != (G4OpticalSurface *) &right);
156 } 203 }
>> 204 ////////////
>> 205 // Methods
>> 206 ////////////
157 207
158 G4int G4OpticalSurface::GetInmax() const { ret << 208 G4int G4OpticalSurface::GetInmax(void) const
159 << 209 {
160 G4int G4OpticalSurface::GetLUTbins() const { r << 210 return indexmax;
>> 211 }
161 212
162 G4int G4OpticalSurface::GetRefMax() const { re << 213 G4int G4OpticalSurface::GetLUTbins(void) const
>> 214 {
>> 215 return LUTbins;
>> 216 }
163 217
164 G4int G4OpticalSurface::GetThetaIndexMax() con << 218 G4int G4OpticalSurface::GetRefMax(void) const
>> 219 {
>> 220 return RefMax;
>> 221 }
>> 222
>> 223 G4int G4OpticalSurface::GetThetaIndexMax(void) const
>> 224 {
>> 225 return thetaIndexMax;
>> 226 }
165 227
166 G4int G4OpticalSurface::GetPhiIndexMax() const << 228 G4int G4OpticalSurface::GetPhiIndexMax(void) const
>> 229 {
>> 230 return phiIndexMax;
>> 231 }
167 232
168 void G4OpticalSurface::DumpInfo() const << 233 void G4OpticalSurface::DumpInfo() const
169 { 234 {
170 // Dump info for surface << 235 // Dump info for surface
171 236
172 G4cout << " Surface type = " << G4int(the << 237 G4cout <<
173 << " Surface finish = " << G4int(the << 238 " Surface type = " << G4int(theType) << G4endl <<
174 << " Surface model = " << G4int(the << 239 " Surface finish = " << G4int(theFinish) << G4endl <<
175 << G4endl << " ----------------- " < << 240 " Surface model = " << G4int(theModel) << G4endl;
176 241
177 if (theModel == glisur) { << 242 G4cout << G4endl;
178 G4cout << " polish: " << polish << G4endl; << 243
179 } << 244 G4cout << " Surface parameter " << G4endl;
180 else { << 245 G4cout << " ----------------- " << G4endl;
181 G4cout << " sigma_alpha: " << sigma_alpha << 246 if (theModel == glisur ){
182 } << 247 G4cout << polish << G4endl;
183 G4cout << G4endl; << 248 }
>> 249 else if (theModel == LUT ){
>> 250 G4cout << sigma_alpha << G4endl;
>> 251 }
>> 252 else {
>> 253 G4cout << sigma_alpha << G4endl;
>> 254 }
>> 255 G4cout << G4endl;
184 } 256 }
185 257
186 void G4OpticalSurface::SetType(const G4Surface 258 void G4OpticalSurface::SetType(const G4SurfaceType& type)
187 { 259 {
188 theType = type; 260 theType = type;
189 ReadDataFile(); << 261
>> 262 if (type == dielectric_LUTDAVIS) {
>> 263 if (!AngularDistributionLUT) AngularDistributionLUT =
>> 264 new G4float[indexmax];
>> 265 ReadLUTDAVISFile();
>> 266 }
>> 267
>> 268 if (type == dielectric_LUT) {
>> 269 if (!AngularDistribution) AngularDistribution =
>> 270 new G4float[incidentIndexMax*thetaIndexMax*phiIndexMax];
>> 271
>> 272 ReadLUTFile();
>> 273 }
>> 274 if (type == dielectric_dichroic) {
>> 275 if (!DichroicVector) DichroicVector = new G4Physics2DVector();
>> 276 ReadDichroicFile();
>> 277 }
190 } 278 }
191 279
192 void G4OpticalSurface::SetFinish(const G4Optic 280 void G4OpticalSurface::SetFinish(const G4OpticalSurfaceFinish finish)
193 { 281 {
194 theFinish = finish; 282 theFinish = finish;
195 ReadDataFile(); << 283 if (theType == dielectric_LUT) {
196 } << 284 if (!AngularDistribution) AngularDistribution =
197 << 285 new G4float[incidentIndexMax*thetaIndexMax*phiIndexMax];
198 void G4OpticalSurface::ReadDataFile() << 286 ReadLUTFile();
199 { << 287 }
200 // type and finish can be set in either orde << 288
201 // for consistency. Need to read file on set << 289 if (theType == dielectric_LUTDAVIS) {
202 switch (theType) { << 290 if (!AngularDistributionLUT) AngularDistributionLUT =
203 case dielectric_LUT: << 291 new G4float[indexmax];
204 if (AngularDistribution == nullptr) { << 292 ReadLUTDAVISFile();
205 AngularDistribution = new G4float[inci << 293
206 } << 294 if (!Reflectivity) Reflectivity = new G4float[RefMax];
207 ReadLUTFile(); << 295 ReadReflectivityLUTFile();
208 break; << 296 }
209 case dielectric_LUTDAVIS: << 297
210 if (AngularDistributionLUT == nullptr) { << 298 if (theType == dielectric_dichroic) {
211 AngularDistributionLUT = new G4float[i << 299 if (!DichroicVector) DichroicVector = new G4Physics2DVector();
212 } << 300 ReadDichroicFile();
213 ReadLUTDAVISFile(); <<
214 <<
215 if (Reflectivity == nullptr) { <<
216 Reflectivity = new G4float[RefMax]; <<
217 } <<
218 ReadReflectivityLUTFile(); <<
219 break; <<
220 case dielectric_dichroic: <<
221 if (DichroicVector == nullptr) { <<
222 DichroicVector = new G4Physics2DVector <<
223 } <<
224 ReadDichroicFile(); <<
225 break; <<
226 default: <<
227 break; <<
228 } 301 }
229 } 302 }
230 303
231 void G4OpticalSurface::ReadLUTFile() 304 void G4OpticalSurface::ReadLUTFile()
232 { 305 {
233 G4String readLUTFileName; << 306 G4String readLUTFileName = " ";
>> 307
>> 308 if (theFinish == polishedlumirrorglue) {
>> 309 readLUTFileName = "PolishedLumirrorGlue.dat";
>> 310 }
>> 311 else if (theFinish == polishedlumirrorair) {
>> 312 readLUTFileName = "PolishedLumirror.dat";
>> 313 }
>> 314 else if (theFinish == polishedteflonair) {
>> 315 readLUTFileName = "PolishedTeflon.dat";
>> 316 }
>> 317 else if (theFinish == polishedtioair) {
>> 318 readLUTFileName = "PolishedTiO.dat";
>> 319 }
>> 320 else if (theFinish == polishedtyvekair) {
>> 321 readLUTFileName = "PolishedTyvek.dat";
>> 322 }
>> 323 else if (theFinish == polishedvm2000glue) {
>> 324 readLUTFileName = "PolishedVM2000Glue.dat";
>> 325 }
>> 326 else if (theFinish == polishedvm2000air) {
>> 327 readLUTFileName = "PolishedVM2000.dat";
>> 328 }
>> 329 else if (theFinish == etchedlumirrorglue) {
>> 330 readLUTFileName = "EtchedLumirrorGlue.dat";
>> 331 }
>> 332 else if (theFinish == etchedlumirrorair) {
>> 333 readLUTFileName = "EtchedLumirror.dat";
>> 334 }
>> 335 else if (theFinish == etchedteflonair) {
>> 336 readLUTFileName = "EtchedTeflon.dat";
>> 337 }
>> 338 else if (theFinish == etchedtioair) {
>> 339 readLUTFileName = "EtchedTiO.dat";
>> 340 }
>> 341 else if (theFinish == etchedtyvekair) {
>> 342 readLUTFileName = "EtchedTyvek.dat";
>> 343 }
>> 344 else if (theFinish == etchedvm2000glue) {
>> 345 readLUTFileName = "EtchedVM2000Glue.dat";
>> 346 }
>> 347 else if (theFinish == etchedvm2000air) {
>> 348 readLUTFileName = "EtchedVM2000.dat";
>> 349 }
>> 350 else if (theFinish == groundlumirrorglue) {
>> 351 readLUTFileName = "GroundLumirrorGlue.dat";
>> 352 }
>> 353 else if (theFinish == groundlumirrorair) {
>> 354 readLUTFileName = "GroundLumirror.dat";
>> 355 }
>> 356 else if (theFinish == groundteflonair) {
>> 357 readLUTFileName = "GroundTeflon.dat";
>> 358 }
>> 359 else if (theFinish == groundtioair) {
>> 360 readLUTFileName = "GroundTiO.dat";
>> 361 }
>> 362 else if (theFinish == groundtyvekair) {
>> 363 readLUTFileName = "GroundTyvek.dat";
>> 364 }
>> 365 else if (theFinish == groundvm2000glue) {
>> 366 readLUTFileName = "GroundVM2000Glue.dat";
>> 367 }
>> 368 else if (theFinish == groundvm2000air) {
>> 369 readLUTFileName = "GroundVM2000.dat";
>> 370 }
234 371
235 switch (theFinish) { << 372 if (readLUTFileName == " ") return;
236 case polishedlumirrorglue: << 373
237 readLUTFileName = "PolishedLumirrorGlue. << 374 char* path = getenv("G4REALSURFACEDATA");
238 break; << 375 if (!path) {
239 case polishedlumirrorair: << 376 G4String excep =
240 readLUTFileName = "PolishedLumirror.z"; << 377 "G4OpBoundaryProcess - G4REALSURFACEDATA environment variable not set";
241 break; << 378 G4Exception("G4OpticalSurface::ReadLUTFile()", "mat310",
242 case polishedteflonair: << 379 FatalException, excep);
243 readLUTFileName = "PolishedTeflon.z"; << 380 return;
244 break; << 381 }
245 case polishedtioair: << 382 G4String pathString(path);
246 readLUTFileName = "PolishedTiO.z"; << 383
247 break; << 384 readLUTFileName = pathString + "/" + readLUTFileName;
248 case polishedtyvekair: << 385
249 readLUTFileName = "PolishedTyvek.z"; << 386 std::ifstream readLUTFileHandle(readLUTFileName, std::ios::in);
250 break; << 387
251 case polishedvm2000glue: << 388 if (readLUTFileHandle) {
252 readLUTFileName = "PolishedVM2000Glue.z" << 389 G4int idxmax = incidentIndexMax*thetaIndexMax*phiIndexMax;
253 break; << 390 for (G4int i = 0; i<idxmax; i++) {
254 case polishedvm2000air: << 391 if (readLUTFileHandle.eof()) break;
255 readLUTFileName = "PolishedVM2000.z"; << 392 readLUTFileHandle >> AngularDistribution[i];
256 break; << 393 }
257 case etchedlumirrorglue: << 394 if (!readLUTFileHandle.bad()) {
258 readLUTFileName = "EtchedLumirrorGlue.z" << 395 G4cout <<"LUT - data file: "<< readLUTFileName <<" read in! "<< G4endl;
259 break; << 396 }
260 case etchedlumirrorair: << 397 else {
261 readLUTFileName = "EtchedLumirror.z"; << 398 G4String excep="LUT - data file: "+readLUTFileName+" not read propery";
262 break; << 399 G4Exception("G4OpticalSurface::ReadLUTFile()", "mat312",
263 case etchedteflonair: << 400 FatalException, excep);
264 readLUTFileName = "EtchedTeflon.z"; << 401 return;
265 break; << 402 }
266 case etchedtioair: <<
267 readLUTFileName = "EtchedTiO.z"; <<
268 break; <<
269 case etchedtyvekair: <<
270 readLUTFileName = "EtchedTyvek.z"; <<
271 break; <<
272 case etchedvm2000glue: <<
273 readLUTFileName = "EtchedVM2000Glue.z"; <<
274 break; <<
275 case etchedvm2000air: <<
276 readLUTFileName = "EtchedVM2000.z"; <<
277 break; <<
278 case groundlumirrorglue: <<
279 readLUTFileName = "GroundLumirrorGlue.z" <<
280 break; <<
281 case groundlumirrorair: <<
282 readLUTFileName = "GroundLumirror.z"; <<
283 break; <<
284 case groundteflonair: <<
285 readLUTFileName = "GroundTeflon.z"; <<
286 break; <<
287 case groundtioair: <<
288 readLUTFileName = "GroundTiO.z"; <<
289 break; <<
290 case groundtyvekair: <<
291 readLUTFileName = "GroundTyvek.z"; <<
292 break; <<
293 case groundvm2000glue: <<
294 readLUTFileName = "GroundVM2000Glue.z"; <<
295 break; <<
296 case groundvm2000air: <<
297 readLUTFileName = "GroundVM2000.z"; <<
298 break; <<
299 default: <<
300 return; <<
301 } <<
302 <<
303 std::istringstream iss; <<
304 ReadCompressedFile(readLUTFileName, iss); <<
305 <<
306 size_t idxmax = incidentIndexMax * thetaInde <<
307 for (size_t i = 0; i < idxmax; ++i) { <<
308 iss >> AngularDistribution[i]; <<
309 } 403 }
310 G4cout << "LUT - data file: " << readLUTFile << 404 else {
>> 405 G4String excep ="LUT - data file: "+readLUTFileName+" not found";
>> 406 G4Exception("G4OpticalSurface::ReadLUTFile()", "mat311",
>> 407 FatalException, excep);
>> 408 return;
>> 409 }
>> 410 readLUTFileHandle.close();
311 } 411 }
312 412
313 void G4OpticalSurface::ReadLUTDAVISFile() 413 void G4OpticalSurface::ReadLUTDAVISFile()
314 { 414 {
315 G4String readLUTDAVISFileName; << 415 G4String readLUTDAVISFileName = " ";
316 416
317 switch (theFinish) { << 417 if (theFinish == Rough_LUT) {
318 case Rough_LUT: << 418 readLUTDAVISFileName = "Rough_LUT.dat";
319 readLUTDAVISFileName = "Rough_LUT.z"; << 419 }else if (theFinish == RoughTeflon_LUT) {
320 break; << 420 readLUTDAVISFileName = "RoughTeflon_LUT.dat";
321 case RoughTeflon_LUT: << 421 }else if (theFinish == RoughESR_LUT) {
322 readLUTDAVISFileName = "RoughTeflon_LUT. << 422 readLUTDAVISFileName = "RoughESR_LUT.dat";
323 break; << 423 }else if (theFinish == RoughESRGrease_LUT) {
324 case RoughESR_LUT: << 424 readLUTDAVISFileName = "RoughESRGrease_LUT.dat";
325 readLUTDAVISFileName = "RoughESR_LUT.z"; << 425 }else if (theFinish == Polished_LUT) {
326 break; << 426 readLUTDAVISFileName = "Polished_LUT.dat";
327 case RoughESRGrease_LUT: << 427 }else if (theFinish == PolishedTeflon_LUT) {
328 readLUTDAVISFileName = "RoughESRGrease_L << 428 readLUTDAVISFileName = "PolishedTeflon_LUT.dat";
329 break; << 429 }else if (theFinish == PolishedESR_LUT) {
330 case Polished_LUT: << 430 readLUTDAVISFileName = "PolishedESR_LUT.dat";
331 readLUTDAVISFileName = "Polished_LUT.z"; << 431 }else if (theFinish == PolishedESRGrease_LUT) {
332 break; << 432 readLUTDAVISFileName = "PolishedESRGrease_LUT.dat";
333 case PolishedTeflon_LUT: << 433 }else if (theFinish == Detector_LUT) {
334 readLUTDAVISFileName = "PolishedTeflon_L << 434 readLUTDAVISFileName = "Detector_LUT.dat";
335 break; << 435 }
336 case PolishedESR_LUT: << 436
337 readLUTDAVISFileName = "PolishedESR_LUT. << 437 if (readLUTDAVISFileName == " ") return;
338 break; << 438
339 case PolishedESRGrease_LUT: << 439 char* path = getenv("G4REALSURFACEDATA");
340 readLUTDAVISFileName = "PolishedESRGreas << 440 if (!path) {
341 break; << 441 G4String excep =
342 case Detector_LUT: << 442 "G4OpBoundaryProcess - G4REALSURFACEDATA environment variable not set";
343 readLUTDAVISFileName = "Detector_LUT.z"; << 443 G4Exception("G4OpticalSurface::ReadLUTDAVISFile()", "mat310",
344 break; << 444 FatalException, excep);
345 default: << 445 return;
346 return; << 446 }
>> 447 G4String pathString(path);
>> 448
>> 449 readLUTDAVISFileName = pathString + "/" + readLUTDAVISFileName;
>> 450
>> 451 std::ifstream readLUTDAVISFileHandle(readLUTDAVISFileName, std::ios::in);
>> 452
>> 453 if (readLUTDAVISFileHandle) {
>> 454
>> 455 for (G4int i = 0; i<indexmax; i++) {
>> 456 if (readLUTDAVISFileHandle.eof()) break;
>> 457 readLUTDAVISFileHandle >> AngularDistributionLUT[i];
>> 458 }
>> 459 if (!readLUTDAVISFileHandle.bad()) {
>> 460 G4cout <<"LUT DAVIS - data file: " <<
>> 461 readLUTDAVISFileName << " read in! " << G4endl;
>> 462 }
>> 463 else {
>> 464 G4String excep="LUT DAVIS- data file: "+readLUTDAVISFileName+" not read propery";
>> 465 G4Exception("G4OpticalSurface::ReadLUTFile()", "mat312",
>> 466 FatalException, excep);
>> 467 return;
>> 468 }
347 } 469 }
348 << 470 else {
349 std::istringstream iss; << 471 G4String excep ="LUT DAVIS- data file: "+readLUTDAVISFileName+" not found";
350 ReadCompressedFile(readLUTDAVISFileName, iss << 472 G4Exception("G4OpticalSurface::ReadLUTFile()", "mat311",
351 << 473 FatalException, excep);
352 for (size_t i = 0; i < indexmax; ++i) { << 474 return;
353 iss >> AngularDistributionLUT[i]; <<
354 } 475 }
355 G4cout << "LUT DAVIS - data file: " << readL << 476 readLUTDAVISFileHandle.close();
356 } 477 }
357 478
358 void G4OpticalSurface::ReadReflectivityLUTFile 479 void G4OpticalSurface::ReadReflectivityLUTFile()
359 { 480 {
360 G4String readReflectivityLUTFileName; << 481 G4String readReflectivityLUTFileName = " ";
361 482
362 switch (theFinish) { << 483 if (theFinish == Rough_LUT) {
363 case Rough_LUT: << 484 readReflectivityLUTFileName = "Rough_LUTR.dat";
364 readReflectivityLUTFileName = "Rough_LUT <<
365 break; <<
366 case RoughTeflon_LUT: <<
367 readReflectivityLUTFileName = "RoughTefl <<
368 break; <<
369 case RoughESR_LUT: <<
370 readReflectivityLUTFileName = "RoughESR_ <<
371 break; <<
372 case RoughESRGrease_LUT: <<
373 readReflectivityLUTFileName = "RoughESRG <<
374 break; <<
375 case Polished_LUT: <<
376 readReflectivityLUTFileName = "Polished_ <<
377 break; <<
378 case PolishedTeflon_LUT: <<
379 readReflectivityLUTFileName = "PolishedT <<
380 break; <<
381 case PolishedESR_LUT: <<
382 readReflectivityLUTFileName = "PolishedE <<
383 break; <<
384 case PolishedESRGrease_LUT: <<
385 readReflectivityLUTFileName = "PolishedE <<
386 break; <<
387 case Detector_LUT: <<
388 readReflectivityLUTFileName = "Detector_ <<
389 break; <<
390 default: <<
391 return; <<
392 } <<
393 <<
394 std::istringstream iss; <<
395 ReadCompressedFile(readReflectivityLUTFileNa <<
396 <<
397 for (size_t i = 0; i < RefMax; ++i) { <<
398 iss >> Reflectivity[i]; <<
399 } <<
400 G4cout << "LUT DAVIS - reflectivity data fil <<
401 << G4endl; <<
402 } <<
403 <<
404 // uncompress one data file into the input str <<
405 void G4OpticalSurface::ReadCompressedFile(cons <<
406 { <<
407 G4String* dataString = nullptr; <<
408 G4String path = G4FindDataDir("G4REALSURFACE <<
409 G4String compfilename = path + "/" + filenam <<
410 // create input stream with binary mode oper <<
411 std::ifstream in(compfilename, std::ios::bin <<
412 if (in.good()) { <<
413 // get current position in the stream (was <<
414 G4int fileSize = (G4int)in.tellg(); <<
415 // set current position being the beginnin <<
416 in.seekg(0, std::ios::beg); <<
417 // create (zlib) byte buffer for the data <<
418 auto compdata = new Bytef[fileSize]; <<
419 while (in) { <<
420 in.read((char*)compdata, fileSize); <<
421 } <<
422 // create (zlib) byte buffer for the uncom <<
423 auto complen = (uLongf)(fileSize * 4); <<
424 auto uncompdata = new Bytef[complen]; <<
425 while (Z_OK != uncompress(uncompdata, &com <<
426 // increase uncompressed byte buffer <<
427 delete[] uncompdata; <<
428 complen *= 2; <<
429 uncompdata = new Bytef[complen]; <<
430 } <<
431 // delete the compressed data buffer <<
432 delete[] compdata; <<
433 // create a string from uncompressed data <<
434 dataString = new G4String((char*)uncompdat <<
435 // delete the uncompressed data buffer <<
436 delete[] uncompdata; <<
437 } 485 }
438 else { << 486 if (theFinish == RoughTeflon_LUT) {
439 G4ExceptionDescription ed; << 487 readReflectivityLUTFileName = "RoughTeflon_LUTR.dat";
440 ed << "Problem while trying to read " + co << 488 }
441 G4Exception("G4OpticalSurface::ReadCompres << 489 if (theFinish== RoughESR_LUT) {
442 return; << 490 readReflectivityLUTFileName = "RoughESR_LUTR.dat";
>> 491 }
>> 492 if (theFinish == RoughESRGrease_LUT) {
>> 493 readReflectivityLUTFileName = "RoughESRGrease_LUTR.dat";
>> 494 }
>> 495 if (theFinish== Polished_LUT) {
>> 496 readReflectivityLUTFileName = "Polished_LUTR.dat";
>> 497 }
>> 498 if (theFinish == PolishedTeflon_LUT) {
>> 499 readReflectivityLUTFileName = "PolishedTeflon_LUTR.dat";
>> 500 }
>> 501 if (theFinish== PolishedESR_LUT) {
>> 502 readReflectivityLUTFileName = "PolishedESR_LUTR.dat";
>> 503 }
>> 504 if (theFinish== PolishedESRGrease_LUT) {
>> 505 readReflectivityLUTFileName = "PolishedESRGrease_LUTR.dat";
>> 506 }
>> 507 if (theFinish== Detector_LUT) {
>> 508 readReflectivityLUTFileName = "Detector_LUTR.dat";
>> 509 }
>> 510
>> 511 if (readReflectivityLUTFileName == " ") return;
>> 512
>> 513 char* path = getenv("G4REALSURFACEDATA");
>> 514 if (!path) {
>> 515 G4String excep =
>> 516 "G4OpBoundaryProcess - G4REALSURFACEDATA environment variable not set";
>> 517 G4Exception("G4OpticalSurface::ReadReflectivityLUTFile()", "mat310",
>> 518 FatalException, excep);
>> 519 return;
>> 520 }
>> 521 G4String pathString(path);
>> 522
>> 523 readReflectivityLUTFileName = pathString + "/" + readReflectivityLUTFileName;
>> 524
>> 525 std::ifstream readReflectivityLUTFileHandle(readReflectivityLUTFileName,
>> 526 std::ios::in);
>> 527
>> 528 if (readReflectivityLUTFileHandle) {
>> 529 for (G4int i = 0; i<RefMax; i++) {
>> 530 if (readReflectivityLUTFileHandle.eof()) break;
>> 531 readReflectivityLUTFileHandle >> Reflectivity[i];
>> 532 }
>> 533 if (!readReflectivityLUTFileHandle.bad()) {
>> 534 G4cout << "Reflectivity LUT DAVIS - data file: " <<
>> 535 readReflectivityLUTFileName << " read in! " << G4endl;
>> 536 }
>> 537 else {
>> 538 G4String excep="Reflectivity LUT DAVIS- data file: "+readReflectivityLUTFileName+" not read propery";
>> 539 G4Exception("G4OpticalSurface::ReadReflectivityLUTFile()", "mat312",
>> 540 FatalException, excep);
>> 541 return;
>> 542 }
443 } 543 }
444 // create the input string stream from the d << 544 else {
445 if (dataString != nullptr) { << 545 G4String excep ="Reflectivity LUT DAVIS - data file: "+readReflectivityLUTFileName+" not found";
446 iss.str(*dataString); << 546 G4Exception("G4OpticalSurface::ReadReflectivityLUTFile()", "mat311",
447 in.close(); << 547 FatalException, excep);
448 delete dataString; << 548 return;
449 G4cout << "G4OpticalSurface: data file " < <<
450 } 549 }
>> 550 readReflectivityLUTFileHandle.close();
451 } 551 }
452 552
453 void G4OpticalSurface::ReadDichroicFile() 553 void G4OpticalSurface::ReadDichroicFile()
454 { 554 {
455 const char* datadir = G4FindDataDir("G4DICHR << 555 const char* datadir = getenv("G4DICHROICDATA");
456 556
457 if (datadir == nullptr) { << 557 if(!datadir) {
458 G4Exception("G4OpticalSurface::ReadDichroi << 558 G4Exception("G4OpticalSurface::ReadDichroicFile()","mat313",
459 "Environment variable G4DICHROICDATA not << 559 FatalException,"Environment variable G4DICHROICDATA not defined");
460 return; 560 return;
461 } 561 }
462 562
463 std::ostringstream ost; 563 std::ostringstream ost;
464 ost << datadir; 564 ost << datadir;
465 std::ifstream fin(ost.str().c_str()); 565 std::ifstream fin(ost.str().c_str());
466 if (! fin.is_open()) { << 566 if( !fin.is_open()) {
467 G4ExceptionDescription ed; 567 G4ExceptionDescription ed;
468 ed << "Dichroic surface data file <" << os << 568 ed << "Dichroic surface data file <" << ost.str().c_str()
469 G4Exception("G4OpticalSurface::ReadDichroi << 569 << "> is not opened!" << G4endl;
>> 570 G4Exception("G4OpticalSurface::ReadDichroicFile()","mat314",
>> 571 FatalException,ed," ");
470 return; 572 return;
471 } 573 }
472 574
473 if (! (DichroicVector->Retrieve(fin))) { << 575 if( !(DichroicVector->Retrieve(fin)) ) {
474 G4ExceptionDescription ed; 576 G4ExceptionDescription ed;
475 ed << "Dichroic surface data file <" << os << 577 ed << "Dichroic surface data file <" << ost.str().c_str()
476 G4Exception("G4OpticalSurface::ReadDichroi << 578 << "> is not opened!" << G4endl;
>> 579 G4Exception("G4OpticalSurface::ReadDichroicFile()","mat315",
>> 580 FatalException,ed," ");
477 return; 581 return;
478 } 582 }
479 583
480 // DichroicVector->SetBicubicInterpolation( << 584 // DichroicVector->SetBicubicInterpolation(true);
481 585
482 G4cout << " *** Dichroic surface data file * 586 G4cout << " *** Dichroic surface data file *** " << G4endl;
483 587
484 auto numberOfXNodes = (G4int)DichroicVector- << 588 G4int numberOfXNodes = DichroicVector->GetLengthX();
485 auto numberOfYNodes = (G4int)DichroicVector- << 589 G4int numberOfYNodes = DichroicVector->GetLengthY();
486 590
487 G4cout << "numberOfXNodes: " << numberOfXNod 591 G4cout << "numberOfXNodes: " << numberOfXNodes << G4endl;
488 G4cout << "numberOfYNodes: " << numberOfYNod 592 G4cout << "numberOfYNodes: " << numberOfYNodes << G4endl;
489 593
490 if (0 > numberOfXNodes || numberOfXNodes >= << 594 if (0 > numberOfXNodes || numberOfXNodes >= INT_MAX) numberOfXNodes = 0;
491 numberOfXNodes = 0; << 595 if (0 > numberOfYNodes || numberOfYNodes >= INT_MAX) numberOfYNodes = 0;
492 } <<
493 if (0 > numberOfYNodes || numberOfYNodes >= <<
494 numberOfYNodes = 0; <<
495 } <<
496 596
497 G4PV2DDataVector xVector; << 597 G4PV2DDataVector xVector;
498 G4PV2DDataVector yVector; << 598 G4PV2DDataVector yVector;
499 599
500 xVector.resize(numberOfXNodes, 0.); << 600 xVector.resize(numberOfXNodes,0.);
501 yVector.resize(numberOfYNodes, 0.); << 601 yVector.resize(numberOfYNodes,0.);
502 602
503 for (G4int i = 0; i < numberOfXNodes; ++i) { << 603 for(G4int i = 0; i<numberOfXNodes; ++i) {
504 G4cout << "i: " << DichroicVector->GetX(i) << 604 G4cout << "i: " << DichroicVector->GetX(i) << G4endl;
505 xVector[i] = DichroicVector->GetX(i); << 605 xVector[i] = DichroicVector->GetX(i);
506 } 606 }
507 for (G4int j = 0; j < numberOfYNodes; ++j) { << 607 for(G4int j = 0; j<numberOfYNodes; ++j) {
508 G4cout << "j: " << DichroicVector->GetY(j) << 608 G4cout << "j: " << DichroicVector->GetY(j) << G4endl;
509 yVector[j] = DichroicVector->GetY(j); << 609 yVector[j] = DichroicVector->GetY(j);
510 } 610 }
511 611
512 for (G4int j = 0; j < numberOfYNodes; ++j) { << 612 for(G4int j = 0; j<numberOfYNodes; ++j) {
513 for (G4int i = 0; i < numberOfXNodes; ++i) << 613 for(G4int i = 0; i<numberOfXNodes; ++i) {
514 G4cout << " i: " << i << " j: " << j << << 614 G4cout << " i: " << i << " j: " << j << " "
515 } << 615 << DichroicVector->GetValue(i,j) << G4endl;
>> 616 }
516 } 617 }
>> 618
>> 619 // G4int idx, idy;
>> 620
>> 621 // for(G4int j = 0; j<numberOfYNodes-1; ++j) {
>> 622 // G4double y = (yVector[j] + yVector[j+1])/2.;
>> 623 // for(G4int i = 0; i<numberOfXNodes-1; ++i) {
>> 624 // G4double x = (xVector[i] + xVector[i+1])/2.;
>> 625 // G4cout << " x: " << x << " y: " << y << " "
>> 626 // << DichroicVector->Value(x,y,idx,idy) << G4endl;
>> 627 // }
>> 628 // }
>> 629
517 } 630 }
518 631