ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DicomPhantomZSliceHeader.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file DicomPhantomZSliceHeader.cc
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 //
29 //
30 
31 #include "globals.hh"
32 #include "G4LogicalVolume.hh"
33 #include "G4MaterialTable.hh"
34 #include "G4Material.hh"
35 #include "G4GeometryTolerance.hh"
36 #include "G4NistManager.hh"
37 
39 
40 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..
42 : fNoVoxelX(0),fNoVoxelY(0),fNoVoxelZ(0),
43  fMinX(0),fMinY(0),fMinZ(0),
44  fMaxX(0),fMaxY(0),fMaxZ(0),
45  fFilename(fname),fSliceLocation(0)
46 {
47 
48 }
49 
50 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..
52 {
53 
54 }
55 
56 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...
58 {
59  //----- Read material indices and names
60  G4int nmate;
61  G4String mateindex;
62  G4String matename;
63  fin >> nmate;
64 #ifdef G4VERBOSE
65  G4cout << " DicomPhantomZSliceHeader reading number of materials "
66  << nmate << G4endl;
67 #endif
68 
69  for( G4int im = 0; im < nmate; im++ ){
70  fin >> mateindex >> matename;
71 #ifdef G4VERBOSE
72  //G4cout << " DicomPhantomZSliceHeader reading material "
73  // << im << " : "<< mateindex << " " << matename << G4endl;
74 #endif
75 
76  if( ! CheckMaterialExists( matename ) ) {
77  G4Exception("DicomPhantomZSliceHeader::DicomPhantomZSliceHeader",
78  "A material is found in file that is not built in the C++ code",
79  FatalErrorInArgument, matename.c_str());
80  }
81 
82  fMaterialNames.push_back(matename);
83  }
84 
85  //----- Read number of voxels
86  fin >> fNoVoxelX >> fNoVoxelY >> fNoVoxelZ;
87 #ifdef G4VERBOSE
88  G4cout << " Number of voxels " << fNoVoxelX << " " << fNoVoxelY
89  << " " << fNoVoxelZ << G4endl;
90 #endif
91 
92  //----- Read minimal and maximal extensions (= walls of phantom)
93  fin >> fMinX >> fMaxX;
94  fin >> fMinY >> fMaxY;
95  fin >> fMinZ >> fMaxZ;
96 #ifdef G4VERBOSE
97  G4cout << " Extension in X " << fMinX << " " << fMaxX << G4endl
98  << " Extension in Y " << fMinY << " " << fMaxY << G4endl
99  << " Extension in Z " << fMinZ << " " << fMaxZ << G4endl;
100 #endif
101 
102  fSliceLocation = 0.5*(fMinZ + fMaxZ);
103 }
104 
105 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
107 {
109  std::vector<G4Material*>::const_iterator matite;
110  for( matite = matTab->begin(); matite != matTab->end(); ++matite ) {
111  if( (*matite)->GetName() == mateName ) { return true; }
112  }
113 
115  if( g4mate ) {
116  return false;
117  } else {
118  return true;
119  }
120 
121 }
122 
123 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
125 {
126  *this = *this + rhs;
127 }
128 
129 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..
131  const DicomPhantomZSliceHeader& rhs )
132 {
133  //----- Check that both slices has the same dimensions
134  if( fNoVoxelX != rhs.GetNoVoxelX()
135  || fNoVoxelY != rhs.GetNoVoxelY() ) {
136  G4cerr << "DicomPhantomZSliceHeader error adding two slice headers:\
137  !!! Different number of voxels: "
138  << " X= " << fNoVoxelX << " =? " << rhs.GetNoVoxelX()
139  << " Y= " << fNoVoxelY << " =? " << rhs.GetNoVoxelY()
140  << " Z= " << fNoVoxelZ << " =? " << rhs.GetNoVoxelZ()
141  << G4endl;
142  G4Exception("DicomPhantomZSliceHeader::DicomPhantomZSliceHeader",
143  "",FatalErrorInArgument,"");
144  }
145  //----- Check that both slices has the same extensions
146  if( fMinX != rhs.GetMinX() || fMaxX != rhs.GetMaxX()
147  || fMinY != rhs.GetMinY() || fMaxY != rhs.GetMaxY() ) {
148  G4cerr << "DicomPhantomZSliceHeader error adding two slice headers:\
149  !!! Different extensions: "
150  << " Xmin= " << fMinX << " =? " << rhs.GetMinX()
151  << " Xmax= " << fMaxX << " =? " << rhs.GetMaxX()
152  << " Ymin= " << fMinY << " =? " << rhs.GetMinY()
153  << " Ymax= " << fMaxY << " =? " << rhs.GetMaxY()
154  << G4endl;
155  G4Exception("DicomPhantomZSliceHeader::operator+","",
157  }
158 
159  //----- Check that both slices have the same materials
160  std::vector<G4String> fMaterialNames2 = rhs.GetMaterialNames();
161  if( fMaterialNames.size() != fMaterialNames2.size() ) {
162  G4cerr << "DicomPhantomZSliceHeader error adding two slice headers:\
163  !!! Different number of materials: " << fMaterialNames.size() << " =? "
164  << fMaterialNames2.size() << G4endl;
165  G4Exception("DicomPhantomZSliceHeader::operator+","",
167  }
168  for( unsigned int ii = 0; ii < fMaterialNames.size(); ii++ ) {
169  if( fMaterialNames[ii] != fMaterialNames2[ii] ) {
170  G4cerr << "DicomPhantomZSliceHeader error adding two slice headers:\
171  !!! Different material number " << ii << " : "
172  << fMaterialNames[ii] << " =? "
173  << fMaterialNames2[ii] << G4endl;
174  G4Exception("DicomPhantomZSliceHeader::operator+","",
176  }
177  }
178 
179  //----- Check that the slices are contiguous in Z
180  if( std::fabs( fMinZ - rhs.GetMaxZ() ) >
182  std::fabs( fMaxZ - rhs.GetMinZ() ) >
184  G4cerr << "DicomPhantomZSliceHeader error adding two slice headers:!!!\
185  Slices are not contiguous in Z "
186  << " Zmin= " << fMinZ << " & " << rhs.GetMinZ()
187  << " Zmax= " << fMaxZ << " & " << rhs.GetMaxZ()
188  << G4endl;
189  G4Exception("DicomPhantomZSliceHeader::operator+","",
191  }
192 
193  //----- Build slice header copying first one
194  DicomPhantomZSliceHeader temp( *this );
195 
196  //----- Add data from second slice header
197  temp.SetMinZ( std::min( fMinZ, rhs.GetMinZ() ) );
198  temp.SetMaxZ( std::max( fMaxZ, rhs.GetMaxZ() ) );
199  temp.SetNoVoxelZ( fNoVoxelZ + rhs.GetNoVoxelZ() );
200 
201  return temp;
202 }
203 
204 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..
206 {
207 
208  G4cout << "DicomPhantomZSliceHeader::Dumping Z Slice data to "
209  << fFilename << "..." << G4endl;
210  //sleep(5);
211 
212  // May seen counter-intuitive (dumping to file you are reading from), but
213  // the reason for this is modification slice spacing
214  if(fMateIDs.size() == 0 || fValues.size() == 0) { ReadDataFromFile(); }
215 
216 
217  std::ofstream out;
218  out.open(fFilename.c_str());
219 
220  if(!out) {
221  G4String descript = "DicomPhantomZSliceHeader::DumpToFile: could not open "
222  +fFilename;
223  G4Exception(descript.c_str(),"", FatalException, "");
224  }
225 
226  out << fMaterialNames.size() << std::endl;
227  for(unsigned int i = 0; i < fMaterialNames.size(); ++i) {
228  out << i << " " << fMaterialNames.at(i) << std::endl;
229  }
230 
231  out << fNoVoxelX << " " << fNoVoxelY << " " << fNoVoxelZ << std::endl;
232  out << fMinX << " " << fMaxX << std::endl;
233  out << fMinY << " " << fMaxY << std::endl;
234  out << fMinZ << " " << fMaxZ << std::endl;
235 
236  for(unsigned int i = 0; i < fMateIDs.size(); ++i)
237  { Print(out,fMateIDs.at(i)," "); }
238 
239  for(unsigned int i = 0; i < fValues.size(); ++i)
240  { Print(out,fValues.at(i)," ",6); }
241 
242  out.close();
243 
244 }
245 
246 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
248 {
249  std::ifstream in;
250  in.open(fFilename.c_str());
251 
252  if(!in) {
253  G4String descript = "DicomPhantomZSliceHeader::DumpToFile: could not open "
254  +fFilename;
255  G4Exception(descript.c_str(),"", FatalException, "");
256  }
257 
259  in >> nMaterials;
260 
261  fMaterialNames.resize(nMaterials,"");
262  for(G4int i = 0; i < nMaterials; ++i) {
263  G4String str1, str2;
264  in >> str1 >> str2;
265  if(!IsInteger(str1)) {
266  G4String descript = "String : " + str1 + " supposed to be integer";
267  G4Exception("DicomPhantomZSliceHeader::ReadDataFromFile - error in \
268  formatting: missing material index","", FatalException,descript.c_str());
269  }
270  G4int index = G4s2n<G4int>(str1);
271  if(index > nMaterials || index < 0) {
272  G4String descript = "Index : " + str1;
273  G4Exception("DicomPhantomZSliceHeader::ReadDataFromFile - error:\
274  bad material index","", FatalException,descript.c_str());
275  }
276  fMaterialNames[index] = str2;
277  }
278 
279  in >> fNoVoxelX >> fNoVoxelY >> fNoVoxelZ;
280 
281  G4double tmpMinX, tmpMinY, tmpMinZ;
282  G4double tmpMaxX, tmpMaxY, tmpMaxZ;
283 
284  in >> tmpMinX >> tmpMaxX;
285  in >> tmpMinY >> tmpMaxY;
286  in >> tmpMinZ >> tmpMaxZ;
287 
288  fMinX = (CheckConsistency(tmpMinX,fMinX,"Min X value")) ?
289  fMinX : ((fMinX == 0) ? tmpMinX : fMinX);
290  fMaxX = (CheckConsistency(tmpMaxX,fMaxX,"Max X value")) ?
291  fMaxX : ((fMaxX == 0) ? tmpMaxX : fMaxX);
292 
293  fMinY = (CheckConsistency(tmpMinY,fMinY,"Min Y value")) ?
294  fMinY : ((fMinY == 0) ? tmpMinY : fMinY);
295  fMaxY = (CheckConsistency(tmpMaxY,fMaxY,"Max Y value")) ?
296  fMaxY : ((fMaxY == 0) ? tmpMaxY : fMaxY);
297 
298  fMinZ = (CheckConsistency(tmpMinZ,fMinZ,"Min Z value")) ?
299  fMinZ : ((fMinZ == 0) ? tmpMinZ : fMinZ);
300  fMaxZ = (CheckConsistency(tmpMaxZ,fMaxZ,"Max Z value")) ?
301  fMaxZ : ((fMaxZ == 0) ? tmpMaxZ : fMaxZ);
302 
303  fMateIDs.clear();
304  fValues.clear();
305  fMateIDs.resize(fNoVoxelY*fNoVoxelZ,std::vector<G4int>(fNoVoxelX,0));
306  fValues.resize(fNoVoxelY*fNoVoxelZ,std::vector<G4double>(fNoVoxelX,0.));
307  for(G4int k = 0; k < fNoVoxelZ; ++k) {
308  for(G4int j = 0; j < fNoVoxelY; ++j) {
309  for(G4int i = 0; i < fNoVoxelX; ++i) {
310  G4int tmpMateID;
311  in >> tmpMateID;
312  G4int row = j*(k+1);
313  fMateIDs[row][i] = tmpMateID;
314  }
315  }
316  }
317 
318  for(G4int k = 0; k < fNoVoxelZ; ++k) {
319  for(G4int j = 0; j < fNoVoxelY; ++j) {
320  for(G4int i = 0; i < fNoVoxelX; ++i) {
321  G4double tmpValue;
322  in >> tmpValue;
323  G4int row = j*(k+1);
324  fValues[row][i] = tmpValue;
325  }
326  }
327  }
328 
329  in.close();
330 }
331 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....