Geant4  10.02.p02
G4GMocrenIO.cc
Go to the documentation of this file.
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 //
27 // $Id: G4GMocrenIO.cc 75179 2013-10-29 10:09:31Z gcosmo $
28 //
29 //
30 // File I/O manager class for writing or reading calcuated dose
31 // distribution and some event information
32 //
33 // Created: Mar. 31, 2009 Akinori Kimura : release for the gMocrenFile driver
34 //
35 // Akinori Kimura
36 // gMocren home page:
37 // http://geant4.kek.jp/gMocren/
38 //
39 //
40 #include "G4GMocrenIO.hh"
41 #include <iostream>
42 #include <ctime>
43 #include <sstream>
44 #include <iomanip>
45 #include <cstdlib>
46 #include <cstring>
47 
48 #include "globals.hh"
49 #include "G4VisManager.hh"
50 
51 #if defined(_WIN32)
52 #define LITTLE_ENDIAN 1234
53 #define BYTE_ORDER LITTLE_ENDIAN
54 #endif
55 
56 const int DOSERANGE = 25000;
57 
58 //----- GMocrenDataPrimitive class in the GMocrenDataIO class-----//
59 template <typename T>
61  clear();
62 }
63 template <typename T>
65  /*
66  std::vector<short *>::iterator itr = image.begin();
67  for(; itr != image.end(); itr++) {
68  delete [] *itr;
69  }
70  */
71 }
72 
73 template <typename T> GMocrenDataPrimitive<T> &
75  if (this == &_right) return *this;
76  for(int i = 0; i < 3; i++) {
77  kSize[i] = _right.kSize[i];
78  kCenter[i] = _right.kCenter[i];
79  }
80  kScale = _right.kScale;
81  for(int i = 0; i < 2; i++) kMinmax[i] = _right.kMinmax[i];
82  int num = kSize[0]*kSize[1];
83  kImage.clear();
84  for(int z = 0; z < kSize[2]; z++) {
85  T * img = new T[num];
86  for(int i = 0; i < num; i++) img[i] =_right.kImage[z][i];
87  kImage.push_back(img);
88  }
89  return *this;
90 }
91 
92 template <typename T> GMocrenDataPrimitive<T> &
94 
96  bool stat = true;
97  for(int i = 0; i < 3; i++) {
98  if(kSize[i] != _right.kSize[i]) stat = false;
99  if(kCenter[i] != _right.kCenter[i]) stat = false;
100  }
101  if(!stat) {
103  G4cout << "Warning: operator + "
104  << " Cannot do the operator +"
105  << G4endl;
106  return *this;
107  }
108 
109  rprim.setSize(kSize);
110  rprim.setCenterPosition(kCenter);
111 
112  T mms[2] = {9e100,-9e100};
113  //if(mms[0] > _right.minmax[0]) mms[0] = _right.minmax[0];
114  //if(mms[1] < _right.minmax[1]) mms[1] = _right.minmax[1];
115 
116  int num = kSize[0]*kSize[1];
117  for(int z = 0; z < kSize[2]; z++) {
118  T * img = new T[num];
119  for(int xy = 0; xy < num; xy++) {
120  img[xy] = kImage[z][xy] + _right.kImage[z][xy];
121  if(mms[0] > img[xy]) mms[0] = img[xy];
122  if(mms[1] < img[xy]) mms[1] = img[xy];
123  }
124  rprim.addImage(img);
125  }
126  rprim.setMinMax(mms);
127 
128  T scl = mms[1]/DOSERANGE;
129  rprim.setScale(scl);
130 
131  return rprim;
132 }
133 
134 template <typename T> GMocrenDataPrimitive<T> &
136 
137  bool stat = true;
138  for(int i = 0; i < 3; i++) {
139  if(kSize[i] != _right.kSize[i]) stat = false;
140  if(kCenter[i] != _right.kCenter[i]) stat = false;
141  }
142  if(!stat) {
144  G4cout << "Warning: operator += " << G4endl
145  << " Cannot do the operator +="
146  << G4endl;
147  return *this;
148  }
149 
150  if(kMinmax[0] > _right.kMinmax[0]) kMinmax[0] = _right.kMinmax[0];
151  if(kMinmax[1] < _right.kMinmax[1]) kMinmax[1] = _right.kMinmax[1];
152 
153  int num = kSize[0]*kSize[1];
154  for(int z = 0; z < kSize[2]; z++) {
155  for(int xy = 0; xy < num; xy++) {
156  kImage[z][xy] += _right.kImage[z][xy];
157  if(kMinmax[0] > kImage[z][xy]) kMinmax[0] = kImage[z][xy];
158  if(kMinmax[1] < kImage[z][xy]) kMinmax[1] = kImage[z][xy];
159  }
160  }
161 
162  kScale = kMinmax[1]/DOSERANGE;
163 
164  return *this;
165 }
166 
167 template <typename T>
169  for(int i = 0; i < 3; i++) {
170  kSize[i] = 0;
171  kCenter[i] = 0.;
172  }
173  kScale = 1.;
174  kMinmax[0] = (T)32109;
175  kMinmax[1] = (T)-32109;
176 
177  clearImage();
178 }
179 template <typename T>
181  typename std::vector<T *>::iterator itr;
182  for(itr = kImage.begin(); itr != kImage.end(); itr++) {
183  delete [] *itr;
184  }
185  kImage.clear();
186 }
187 template <typename T>
189  for(int i = 0; i < 3; i++) kSize[i] = _size[i];
190 }
191 template <typename T>
193  for(int i = 0; i < 3; i++) _size[i] = kSize[i];
194 }
195 template <typename T>
196 void GMocrenDataPrimitive<T>::setScale(double & _scale) {
197  kScale = _scale;
198 }
199 template <typename T>
201  return kScale;
202 }
203 template <typename T>
205  for(int i = 0; i < 2; i++) kMinmax[i] = _minmax[i];
206 }
207 template <typename T>
209  for(int i = 0; i < 2; i++) _minmax[i] = kMinmax[i];
210 
211 }
212 template <typename T>
213 void GMocrenDataPrimitive<T>::setImage(std::vector<T *> & _image) {
214  kImage = _image;
215 }
216 template <typename T>
218  kImage.push_back(_image);
219 }
220 template <typename T>
221 std::vector<T *> & GMocrenDataPrimitive<T>::getImage() {
222  return kImage;
223 }
224 template <typename T>
226  if(_z >= (int)kImage.size()) return 0;
227  return kImage[_z];
228 }
229 template <typename T>
231  for(int i = 0; i < 3; i++) kCenter[i] = _center[i];
232 }
233 template <typename T>
235  for(int i = 0; i < 3; i++) _center[i] = kCenter[i];
236 }
237 template <typename T>
238 void GMocrenDataPrimitive<T>::setName(std::string & _name) {
239  kDataName = _name;
240 }
241 template <typename T>
243  return kDataName;
244 }
245 
246 
247 
248 
249 
251  kTrack.clear();
252  for(int i = 0; i < 3; i++) kColor[i] = 0;
253 }
254 
255 void GMocrenTrack::addStep(float _startx, float _starty, float _startz,
256  float _endx, float _endy, float _endz) {
257  struct Step step;
258  step.startPoint[0] = _startx;
259  step.startPoint[1] = _starty;
260  step.startPoint[2] = _startz;
261  step.endPoint[0] = _endx;
262  step.endPoint[1] = _endy;
263  step.endPoint[2] = _endz;
264  kTrack.push_back(step);
265 }
266 void GMocrenTrack::getStep(float & _startx, float & _starty, float & _startz,
267  float & _endx, float & _endy, float & _endz,
268  int _num) {
269  if(_num >= (int)kTrack.size()) {
271  G4cout << "GMocrenTrack::getStep(...) Error: "
272  << "invalid step # : " << _num << G4endl;
273  return;
274  }
275 
276  _startx = kTrack[_num].startPoint[0];
277  _starty = kTrack[_num].startPoint[1];
278  _startz = kTrack[_num].startPoint[2];
279  _endx = kTrack[_num].endPoint[0];
280  _endy = kTrack[_num].endPoint[1];
281  _endz = kTrack[_num].endPoint[2];
282 }
283 void GMocrenTrack::translate(std::vector<float> & _translate) {
284  std::vector<struct Step>::iterator itr = kTrack.begin();
285  for(; itr != kTrack.end(); itr++) {
286  for(int i = 0; i < 3; i++ ) {
287  itr->startPoint[i] += _translate[i];
288  itr->endPoint[i] += _translate[i];
289  }
290  }
291 }
292 
293 
294 
295 
296 
297 
298 
299 
300 
302  kDetector.clear();
303  for(int i = 0; i < 3; i++) kColor[i] = 0;
304 }
305 
306 void GMocrenDetector::addEdge(float _startx, float _starty, float _startz,
307  float _endx, float _endy, float _endz) {
308  struct Edge edge;
309  edge.startPoint[0] = _startx;
310  edge.startPoint[1] = _starty;
311  edge.startPoint[2] = _startz;
312  edge.endPoint[0] = _endx;
313  edge.endPoint[1] = _endy;
314  edge.endPoint[2] = _endz;
315  kDetector.push_back(edge);
316 }
317 void GMocrenDetector::getEdge(float & _startx, float & _starty, float & _startz,
318  float & _endx, float & _endy, float & _endz,
319  int _num) {
320  if(_num >= (int)kDetector.size()) {
322  G4cout << "GMocrenDetector::getEdge(...) Error: "
323  << "invalid edge # : " << _num << G4endl;
324  return;
325  }
326 
327  _startx = kDetector[_num].startPoint[0];
328  _starty = kDetector[_num].startPoint[1];
329  _startz = kDetector[_num].startPoint[2];
330  _endx = kDetector[_num].endPoint[0];
331  _endy = kDetector[_num].endPoint[1];
332  _endz = kDetector[_num].endPoint[2];
333 }
334 void GMocrenDetector::translate(std::vector<float> & _translate) {
335  std::vector<struct Edge>::iterator itr = kDetector.begin();
336  for(; itr != kDetector.end(); itr++) {
337  for(int i = 0; i < 3; i++) {
338  itr->startPoint[i] += _translate[i];
339  itr->endPoint[i] += _translate[i];
340  }
341  }
342 }
343 
344 
345 
346 
347 
348 
349 
350 
351 
352 // file information
353 std::string G4GMocrenIO::kId;
354 std::string G4GMocrenIO::kVersion = "2.0.0";
357 
358 #if BYTE_ORDER == LITTLE_ENDIAN
360 #else
361 char G4GMocrenIO::kLittleEndianOutput = false; // Big endian
362 #endif
363 std::string G4GMocrenIO::kComment;
364 //
365 std::string G4GMocrenIO::kFileName = "dose.gdd";
366 
367 //
368 unsigned int G4GMocrenIO::kPointerToModalityData = 0;
369 std::vector<unsigned int> G4GMocrenIO::kPointerToDoseDistData;
370 unsigned int G4GMocrenIO::kPointerToROIData = 0;
371 unsigned int G4GMocrenIO::kPointerToTrackData = 0;
372 unsigned int G4GMocrenIO::kPointerToDetectorData = 0;
373 
374 // modality
375 float G4GMocrenIO::kVoxelSpacing[3] = {0., 0., 0.};
376 class GMocrenDataPrimitive<short> G4GMocrenIO::kModality;
377 std::vector<float> G4GMocrenIO::kModalityImageDensityMap;
378 std::string G4GMocrenIO::kModalityUnit = "g/cm3 "; // 12 Bytes
379 
380 // dose
381 std::vector<class GMocrenDataPrimitive<double> > G4GMocrenIO::kDose;
382 std::string G4GMocrenIO::kDoseUnit = "keV "; // 12 Bytes
383 
384 // ROI
385 std::vector<class GMocrenDataPrimitive<short> > G4GMocrenIO::kRoi;
386 
387 // track
388 std::vector<float *> G4GMocrenIO::kSteps;
389 std::vector<unsigned char *> G4GMocrenIO::kStepColors;
390 std::vector<class GMocrenTrack> G4GMocrenIO::kTracks;
391 
392 // detector
393 std::vector<class GMocrenDetector> G4GMocrenIO::kDetectors;
394 
395 // verbose
396 int G4GMocrenIO::kVerbose = 0;
397 
398 const int IDLENGTH = 21;
399 const int VERLENGTH = 6;
400 
401 // constructor
403  : kTracksWillBeStored(true) {
404  ;
405 }
406 
407 // destructor
409  ;
410 }
411 
412 // initialize
414 
415  kId.clear();
416  kVersion = "2.0.0";
417  kNumberOfEvents = 0;
418  kLittleEndianInput = true;
419 #if BYTE_ORDER == LITTLE_ENDIAN
420  kLittleEndianOutput = true;
421 #else // Big endian
422  kLittleEndianOutput = false;
423 #endif
424  kComment.clear();
425  kFileName = "dose.gdd";
427  kPointerToDoseDistData.clear();
428  kPointerToROIData = 0;
430  // modality
431  for(int i = 0; i < 3; i++) kVoxelSpacing[i] = 0.;
432  kModality.clear();
433  kModalityImageDensityMap.clear();
434  kModalityUnit = "g/cm3 "; // 12 Bytes
435  // dose
436  kDose.clear();
437  kDoseUnit = "keV "; // 12 Bytes
438  // ROI
439  kRoi.clear();
440  // track
441  std::vector<float *>::iterator itr;
442  for(itr = kSteps.begin(); itr != kSteps.end(); itr++) delete [] *itr;
443  kSteps.clear();
444  std::vector<unsigned char *>::iterator citr;
445  for(citr = kStepColors.begin(); citr != kStepColors.end(); citr++)
446  delete [] *citr;
447  kStepColors.clear();
448  kTracksWillBeStored = true;
449 
450  // verbose
451  kVerbose = 0;
452 }
453 
455  return storeData4();
456 }
457 //
458 bool G4GMocrenIO::storeData(char * _filename) {
459  return storeData4(_filename);
460 }
461 
463 
464  bool DEBUG = false;//
465 
466  if(DEBUG || kVerbose > 0)
467  G4cout << ">>>>>>> store data (ver.4) <<<<<<<" << G4endl;
468  if(DEBUG || kVerbose > 0)
469  G4cout << " " << kFileName << G4endl;
470 
471  // output file open
472  std::ofstream ofile(kFileName.c_str(),
473  std::ios_base::out|std::ios_base::binary);
474  if(DEBUG || kVerbose > 0)
475  G4cout << " file open status: " << ofile.rdbuf() << G4endl;
476 
477  // file identifier
478  ofile.write("gMocren ", 8);
479 
480  // file version
481  unsigned char ver = 0x04;
482  ofile.write((char *)&ver, 1);
483 
484  // endian
485  //ofile.write((char *)&kLittleEndianOutput, sizeof(char));
486  char littleEndian = 0x01;
487  ofile.write((char *)&littleEndian, sizeof(char));
488  if(DEBUG || kVerbose > 0) {
489  //G4cout << "Endian: " << (int)kLittleEndianOutput << G4endl;
490  G4cout << "Endian: " << (int)littleEndian << G4endl;
491  }
492 
493  // for inverting the byte order
494  float ftmp[6];
495  int itmp[6];
496  short stmp[6];
497 
498  // comment length (fixed size)
499  int commentLength = 1024;
500  if(kLittleEndianOutput) {
501  ofile.write((char *)&commentLength, 4);
502  } else {
503  invertByteOrder((char *)&commentLength, itmp[0]);
504  ofile.write((char *)itmp, 4);
505  }
506 
507  // comment
508  char cmt[1025];
509  for(int i = 0; i < 1025; i++) cmt[i] = '\0';
510  //std::strncpy(cmt, kComment.c_str(), 1024);
511  const char * cmnt = kComment.c_str();
512  size_t lcm = std::strlen(cmnt);
513  if(lcm > 1024) lcm = 1024;
514  std::strncpy(cmt, cmnt, lcm);
515  ofile.write((char *)cmt, 1024);
516  if(DEBUG || kVerbose > 0) {
517  G4cout << "Data comment : "
518  << kComment << G4endl;
519  }
520 
521  // voxel spacings for all images
522  if(kLittleEndianOutput) {
523  ofile.write((char *)kVoxelSpacing, 12);
524  } else {
525  for(int j = 0; j < 3; j++)
526  invertByteOrder((char *)&kVoxelSpacing[j], ftmp[j]);
527  ofile.write((char *)ftmp, 12);
528  }
529  if(DEBUG || kVerbose > 0) {
530  G4cout << "Voxel spacing : ("
531  << kVoxelSpacing[0] << ", "
532  << kVoxelSpacing[1] << ", "
533  << kVoxelSpacing[2]
534  << ") mm " << G4endl;
535  }
536 
537  calcPointers4();
539 
540  // offset from file starting point to the modality image data
541  if(kLittleEndianOutput) {
542  ofile.write((char *)&kPointerToModalityData, 4);
543  } else {
544  invertByteOrder((char *)&kPointerToModalityData, itmp[0]);
545  ofile.write((char *)itmp, 4);
546  }
547 
548  // # of dose distributions
549  //int nDoseDist = (int)pointerToDoseDistData.size();
550  int nDoseDist = getNumDoseDist();
551  if(kLittleEndianOutput) {
552  ofile.write((char *)&nDoseDist, 4);
553  } else {
554  invertByteOrder((char *)&nDoseDist, itmp[0]);
555  ofile.write((char *)itmp, 4);
556  }
557 
558  // offset from file starting point to the dose image data
559  if(kLittleEndianOutput) {
560  for(int i = 0; i < nDoseDist; i++) {
561  ofile.write((char *)&kPointerToDoseDistData[i], 4);
562  }
563  } else {
564  for(int i = 0; i < nDoseDist; i++) {
565  invertByteOrder((char *)&kPointerToDoseDistData[i], itmp[0]);
566  ofile.write((char *)itmp, 4);
567  }
568  }
569 
570  // offset from file starting point to the ROI image data
571  if(kLittleEndianOutput) {
572  ofile.write((char *)&kPointerToROIData, 4);
573  } else {
574  invertByteOrder((char *)&kPointerToROIData, itmp[0]);
575  ofile.write((char *)itmp, 4);
576  }
577 
578  // offset from file starting point to the track data
579  if(kLittleEndianOutput) {
580  ofile.write((char *)&kPointerToTrackData, 4);
581  } else {
582  invertByteOrder((char *)&kPointerToTrackData, itmp[0]);
583  ofile.write((char *)itmp, 4);
584  }
585 
586  // offset from file starting point to the detector data
587  if(kLittleEndianOutput) {
588  ofile.write((char *)&kPointerToDetectorData, 4);
589  } else {
590  invertByteOrder((char *)&kPointerToDetectorData, itmp[0]);
591  ofile.write((char *)itmp, 4);
592  }
593 
594  if(DEBUG || kVerbose > 0) {
595  G4cout << "Each pointer to data : "
596  << kPointerToModalityData << ", ";
597  for(int i = 0; i < nDoseDist; i++) {
598  G4cout << kPointerToDoseDistData[i] << ", ";
599  }
600  G4cout << kPointerToROIData << ", "
601  << kPointerToTrackData << ", "
603  << G4endl;
604  }
605 
606  //----- modality image -----//
607 
608  int size[3];
609  float scale;
610  short minmax[2];
611  float fCenter[3];
612  int iCenter[3];
613  // modality image size
614  kModality.getSize(size);
615 
616  if(kLittleEndianOutput) {
617  ofile.write((char *)size, 3*sizeof(int));
618  } else {
619  for(int j = 0; j < 3; j++)
620  invertByteOrder((char *)&size[j], itmp[j]);
621  ofile.write((char *)itmp, 12);
622  }
623 
624  if(DEBUG || kVerbose > 0) {
625  G4cout << "Modality image size : ("
626  << size[0] << ", "
627  << size[1] << ", "
628  << size[2] << ")"
629  << G4endl;
630  }
631 
632  // modality image max. & min.
633  kModality.getMinMax(minmax);
634  if(kLittleEndianOutput) {
635  ofile.write((char *)minmax, 4);
636  } else {
637  for(int j = 0; j < 2; j++)
638  invertByteOrder((char *)&minmax[j], stmp[j]);
639  ofile.write((char *)stmp, 4);
640  }
641 
642  // modality image unit
643  char munit[13] = "g/cm3\0";
644  ofile.write((char *)munit, 12);
645 
646  // modality image scale
647  scale = (float)kModality.getScale();
648  if(kLittleEndianOutput) {
649  ofile.write((char *)&scale, 4);
650  } else {
651  invertByteOrder((char *)&scale, ftmp[0]);
652  ofile.write((char *)ftmp, 4);
653  }
654  if(DEBUG || kVerbose > 0) {
655  G4cout << "Modality image min., max., scale : "
656  << minmax[0] << ", "
657  << minmax[1] << ", "
658  << scale << G4endl;
659  }
660 
661  // modality image
662  int psize = size[0]*size[1];
663  if(DEBUG || kVerbose > 0) G4cout << "Modality image : ";
664  for(int i = 0; i < size[2]; i++) {
665  short * image = kModality.getImage(i);
666  if(kLittleEndianOutput) {
667  ofile.write((char *)image, psize*sizeof(short));
668  } else {
669  for(int j = 0; j < psize; j++) {
670  invertByteOrder((char *)&image[j], stmp[0]);
671  ofile.write((char *)stmp, 2);
672  }
673  }
674 
675  if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << image[(size_t)(psize*0.55)] << ", ";
676  }
677  if(DEBUG || kVerbose > 0) G4cout << G4endl;
678 
679  // modality desity map for CT value
680  size_t msize = minmax[1] - minmax[0]+1;
681  if(DEBUG || kVerbose > 0)
682  G4cout << "modality image : " << minmax[0] << ", " << minmax[1] << G4endl;
683  float * pdmap = new float[msize];
684  for(int i = 0; i < (int)msize; i++) pdmap[i] =kModalityImageDensityMap[i];
685 
686  if(kLittleEndianOutput) {
687  ofile.write((char *)pdmap, msize*sizeof(float));
688  } else {
689  for(int j = 0; j < (int)msize; j++) {
690  invertByteOrder((char *)&pdmap[j], ftmp[0]);
691  ofile.write((char *)ftmp, 4);
692  }
693  }
694 
695  if(DEBUG || kVerbose > 0) {
696  G4cout << "density map : " << std::ends;
697  for(int i = 0; i < (int)msize; i+=50)
698  G4cout <<kModalityImageDensityMap[i] << ", ";
699  G4cout << G4endl;
700  }
701  delete [] pdmap;
702 
703 
704  //----- dose distribution image -----//
705 
706  if(!isDoseEmpty()) {
707 
709 
710  for(int ndose = 0; ndose < nDoseDist; ndose++) {
711  // dose distrbution image size
712  kDose[ndose].getSize(size);
713  if(kLittleEndianOutput) {
714  ofile.write((char *)size, 3*sizeof(int));
715  } else {
716  for(int j = 0; j < 3; j++)
717  invertByteOrder((char *)&size[j], itmp[j]);
718  ofile.write((char *)itmp, 12);
719  }
720  if(DEBUG || kVerbose > 0) {
721  G4cout << "Dose dist. [" << ndose << "] image size : ("
722  << size[0] << ", "
723  << size[1] << ", "
724  << size[2] << ")"
725  << G4endl;
726  }
727 
728  // dose distribution max. & min.
729  getShortDoseDistMinMax(minmax, ndose);
730  if(kLittleEndianOutput) {
731  ofile.write((char *)minmax, 2*2); // sizeof(shorft)*2
732  } else {
733  for(int j = 0; j < 2; j++)
734  invertByteOrder((char *)&minmax[j], stmp[j]);
735  ofile.write((char *)stmp, 4);
736  }
737 
738  // dose distribution unit
739  char cdunit[13];
740  for(int i = 0; i < 13; i++) cdunit[i] = '\0';
741  const char * cu = kDoseUnit.c_str();
742  size_t lcu = std::strlen(cu);
743  if(lcu > 12) lcu = 12;
744  std::strncpy(cdunit, cu, lcu);
745  ofile.write((char *)cdunit, 12);
746  if(DEBUG || kVerbose > 0) {
747  G4cout << "Dose dist. unit : " << kDoseUnit << G4endl;
748  }
749 
750  // dose distribution scaling
751  double dscale;
752  dscale = getDoseDistScale(ndose);
753  scale = float(dscale);
754  if(kLittleEndianOutput) {
755  ofile.write((char *)&scale, 4);
756  } else {
757  invertByteOrder((char *)&scale, ftmp[0]);
758  ofile.write((char *)ftmp, 4);
759  }
760  if(DEBUG || kVerbose > 0) {
761  G4cout << "Dose dist. [" << ndose
762  << "] image min., max., scale : "
763  << minmax[0] << ", "
764  << minmax[1] << ", "
765  << scale << G4endl;
766  }
767 
768  // dose distribution image
769  int dsize = size[0]*size[1];
770  short * dimage = new short[dsize];
771  for(int z = 0; z < size[2]; z++) {
772  getShortDoseDist(dimage, z, ndose);
773  if(kLittleEndianOutput) {
774  ofile.write((char *)dimage, dsize*2); //sizeof(short)
775  } else {
776  for(int j = 0; j < dsize; j++) {
777  invertByteOrder((char *)&dimage[j], stmp[0]);
778  ofile.write((char *)stmp, 2);
779  }
780  }
781 
782  if(DEBUG || kVerbose > 0) {
783  for(int j = 0; j < dsize; j++) {
784  if(dimage[j] < 0)
785  G4cout << "[" << j << "," << z << "]"
786  << dimage[j] << ", ";
787  }
788  }
789  }
790  if(DEBUG || kVerbose > 0) G4cout << G4endl;
791  delete [] dimage;
792 
793  // relative location of the dose distribution image for
794  // the modality image
795  getDoseDistCenterPosition(fCenter, ndose);
796  for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
797  if(kLittleEndianOutput) {
798  ofile.write((char *)iCenter, 3*4); // 3*sizeof(int)
799  } else {
800  for(int j = 0; j < 3; j++)
801  invertByteOrder((char *)&iCenter[j], itmp[j]);
802  ofile.write((char *)itmp, 12);
803  }
804  if(DEBUG || kVerbose > 0) {
805  G4cout << "Dose dist. [" << ndose
806  << "]image relative location : ("
807  << iCenter[0] << ", "
808  << iCenter[1] << ", "
809  << iCenter[2] << ")" << G4endl;
810  }
811 
812  // dose distribution name
813  std::string name = getDoseDistName(ndose);
814  if(name.size() == 0) name = "dose";
815  name.resize(80);
816  ofile.write((char *)name.c_str(), 80);
817  if(DEBUG || kVerbose > 0) {
818  G4cout << "Dose dist. name : " << name << G4endl;
819  }
820 
821  }
822  }
823 
824  //----- ROI image -----//
825  if(!isROIEmpty()) {
826  // ROI image size
827  kRoi[0].getSize(size);
828  if(kLittleEndianOutput) {
829  ofile.write((char *)size, 3*sizeof(int));
830  } else {
831  for(int j = 0; j < 3; j++)
832  invertByteOrder((char *)&size[j], itmp[j]);
833  ofile.write((char *)itmp, 12);
834  }
835  if(DEBUG || kVerbose > 0) {
836  G4cout << "ROI image size : ("
837  << size[0] << ", "
838  << size[1] << ", "
839  << size[2] << ")"
840  << G4endl;
841  }
842 
843  // ROI max. & min.
844  kRoi[0].getMinMax(minmax);
845  if(kLittleEndianOutput) {
846  ofile.write((char *)minmax, sizeof(short)*2);
847  } else {
848  for(int j = 0; j < 2; j++)
849  invertByteOrder((char *)&minmax[j], stmp[j]);
850  ofile.write((char *)stmp, 4);
851  }
852 
853  // ROI distribution scaling
854  scale = (float)kRoi[0].getScale();
855  if(kLittleEndianOutput) {
856  ofile.write((char *)&scale, sizeof(float));
857  } else {
858  invertByteOrder((char *)&scale, ftmp[0]);
859  ofile.write((char *)ftmp, 4);
860  }
861  if(DEBUG || kVerbose > 0) {
862  G4cout << "ROI image min., max., scale : "
863  << minmax[0] << ", "
864  << minmax[1] << ", "
865  << scale << G4endl;
866  }
867 
868  // ROI image
869  int rsize = size[0]*size[1];
870  for(int i = 0; i < size[2]; i++) {
871  short * rimage = kRoi[0].getImage(i);
872  if(kLittleEndianOutput) {
873  ofile.write((char *)rimage, rsize*sizeof(short));
874  } else {
875  for(int j = 0; j < rsize; j++) {
876  invertByteOrder((char *)&rimage[j], stmp[0]);
877  ofile.write((char *)stmp, 2);
878  }
879  }
880 
881  }
882 
883  // ROI relative location
884  kRoi[0].getCenterPosition(fCenter);
885  for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
886  if(kLittleEndianOutput) {
887  ofile.write((char *)iCenter, 3*sizeof(int));
888  } else {
889  for(int j = 0; j < 3; j++)
890  invertByteOrder((char *)&iCenter[j], itmp[j]);
891  ofile.write((char *)itmp, 12);
892  }
893  if(DEBUG || kVerbose > 0) {
894  G4cout << "ROI image relative location : ("
895  << iCenter[0] << ", "
896  << iCenter[1] << ", "
897  << iCenter[2] << ")" << G4endl;
898  }
899  }
900 
901  //----- track information -----//
902  // number of track
903  if(kPointerToTrackData > 0) {
904 
905  int ntrk = kTracks.size();
906  if(kLittleEndianOutput) {
907  ofile.write((char *)&ntrk, sizeof(int));
908  } else {
909  invertByteOrder((char *)&ntrk, itmp[0]);
910  ofile.write((char *)itmp, 4);
911  }
912  if(DEBUG || kVerbose > 0) {
913  G4cout << "# of tracks : "
914  << ntrk << G4endl;
915  }
916 
917  for(int nt = 0; nt < ntrk; nt++) {
918 
919  // # of steps in a track
920  int nsteps = kTracks[nt].getNumberOfSteps();
921  if(kLittleEndianOutput) {
922  ofile.write((char *)&nsteps, sizeof(int));
923  } else {
924  invertByteOrder((char *)&nsteps, itmp[0]);
925  ofile.write((char *)itmp, 4);
926  }
927  if(DEBUG || kVerbose > 0) {
928  G4cout << "# of steps : " << nsteps << G4endl;
929  }
930 
931  // track color
932  unsigned char tcolor[3];
933  kTracks[nt].getColor(tcolor);
934  ofile.write((char *)tcolor, 3);
935 
936  // steps
937  float stepPoints[6];
938  for(int isteps = 0; isteps < nsteps; isteps++) {
939  kTracks[nt].getStep(stepPoints[0], stepPoints[1], stepPoints[2],
940  stepPoints[3], stepPoints[4], stepPoints[5],
941  isteps);
942 
943  if(kLittleEndianOutput) {
944  ofile.write((char *)stepPoints, sizeof(float)*6);
945  } else {
946  for(int j = 0; j < 6; j++)
947  invertByteOrder((char *)&stepPoints[j], ftmp[j]);
948  ofile.write((char *)ftmp, 24);
949  }
950  }
951  }
952  }
953 
954  //----- detector information -----//
955  // number of detectors
956  if(kPointerToDetectorData > 0) {
957  int ndet = kDetectors.size();
958  if(kLittleEndianOutput) {
959  ofile.write((char *)&ndet, sizeof(int));
960  } else {
961  invertByteOrder((char *)&ndet, itmp[0]);
962  ofile.write((char *)itmp, 4);
963  }
964  if(DEBUG || kVerbose > 0) {
965  G4cout << "# of detectors : "
966  << ndet << G4endl;
967  }
968 
969  for(int nd = 0; nd < ndet; nd++) {
970 
971  // # of edges of a detector
972  int nedges = kDetectors[nd].getNumberOfEdges();
973  if(kLittleEndianOutput) {
974  ofile.write((char *)&nedges, sizeof(int));
975  } else {
976  invertByteOrder((char *)&nedges, itmp[0]);
977  ofile.write((char *)itmp, 4);
978  }
979  if(DEBUG || kVerbose > 0) {
980  G4cout << "# of edges in a detector : " << nedges << G4endl;
981  }
982 
983  // edges
984  float edgePoints[6];
985  for(int ne = 0; ne < nedges; ne++) {
986  kDetectors[nd].getEdge(edgePoints[0], edgePoints[1], edgePoints[2],
987  edgePoints[3], edgePoints[4], edgePoints[5],
988  ne);
989 
990  if(kLittleEndianOutput) {
991  ofile.write((char *)edgePoints, sizeof(float)*6);
992  } else {
993  for(int j = 0; j < 6; j++)
994  invertByteOrder((char *)&edgePoints[j], ftmp[j]);
995  ofile.write((char *)ftmp, 24);
996  }
997 
998  if(DEBUG || kVerbose > 0) {
999  if(ne < 1) {
1000  G4cout << " edge : (" << edgePoints[0] << ", "
1001  << edgePoints[1] << ", "
1002  << edgePoints[2] << ") - ("
1003  << edgePoints[3] << ", "
1004  << edgePoints[4] << ", "
1005  << edgePoints[5] << ")" << G4endl;
1006  }
1007  }
1008  }
1009 
1010  // detector color
1011  unsigned char dcolor[3];
1012  kDetectors[nd].getColor(dcolor);
1013  ofile.write((char *)dcolor, 3);
1014  if(DEBUG || kVerbose > 0) {
1015  G4cout << " rgb : (" << (int)dcolor[0] << ", "
1016  << (int)dcolor[1] << ", "
1017  << (int)dcolor[2] << ")" << G4endl;
1018  }
1019 
1020  // detector name
1021  std::string dname = kDetectors[nd].getName();
1022  dname.resize(80);
1023  ofile.write((char *)dname.c_str(), 80);
1024  if(DEBUG || kVerbose > 0) {
1025  G4cout << " detector name : " << dname << G4endl;
1026 
1027  }
1028  }
1029  }
1030 
1031  // file end mark
1032  ofile.write("END", 3);
1033 
1034  ofile.close();
1035  if(DEBUG || kVerbose > 0)
1036  G4cout << ">>>> closed gdd file: " << kFileName << G4endl;
1037 
1038  return true;
1039 }
1041 
1042  if(kVerbose > 0) G4cout << ">>>>>>> store data (ver.3) <<<<<<<" << G4endl;
1043  if(kVerbose > 0) G4cout << " " << kFileName << G4endl;
1044 
1045  bool DEBUG = false;//
1046 
1047  // output file open
1048  std::ofstream ofile(kFileName.c_str(),
1049  std::ios_base::out|std::ios_base::binary);
1050 
1051  // file identifier
1052  ofile.write("gMocren ", 8);
1053 
1054  // file version
1055  unsigned char ver = 0x03;
1056  ofile.write((char *)&ver, 1);
1057 
1058  // endian
1059  ofile.write((char *)&kLittleEndianOutput, sizeof(char));
1060 
1061  // comment length (fixed size)
1062  int commentLength = 1024;
1063  ofile.write((char *)&commentLength, 4);
1064 
1065  // comment
1066  char cmt[1025];
1067  std::strncpy(cmt, kComment.c_str(), 1024);
1068  ofile.write((char *)cmt, 1024);
1069  if(DEBUG || kVerbose > 0) {
1070  G4cout << "Data comment : "
1071  << kComment << G4endl;
1072  }
1073 
1074  // voxel spacings for all images
1075  ofile.write((char *)kVoxelSpacing, 12);
1076  if(DEBUG || kVerbose > 0) {
1077  G4cout << "Voxel spacing : ("
1078  << kVoxelSpacing[0] << ", "
1079  << kVoxelSpacing[1] << ", "
1080  << kVoxelSpacing[2]
1081  << ") mm " << G4endl;
1082  }
1083 
1084  calcPointers3();
1085 
1086  // offset from file starting point to the modality image data
1087  ofile.write((char *)&kPointerToModalityData, 4);
1088 
1089  // # of dose distributions
1090  //int nDoseDist = (int)pointerToDoseDistData.size();
1091  int nDoseDist = getNumDoseDist();
1092  ofile.write((char *)&nDoseDist, 4);
1093 
1094  // offset from file starting point to the dose image data
1095  for(int i = 0; i < nDoseDist; i++) {
1096  ofile.write((char *)&kPointerToDoseDistData[i], 4);
1097  }
1098 
1099  // offset from file starting point to the ROI image data
1100  ofile.write((char *)&kPointerToROIData, 4);
1101 
1102  // offset from file starting point to the track data
1103  ofile.write((char *)&kPointerToTrackData, 4);
1104  if(DEBUG || kVerbose > 0) {
1105  G4cout << "Each pointer to data : "
1106  << kPointerToModalityData << ", ";
1107  for(int i = 0; i < nDoseDist; i++) {
1108  G4cout << kPointerToDoseDistData[i] << ", ";
1109  }
1110  G4cout << kPointerToROIData << ", "
1112  }
1113 
1114  //----- modality image -----//
1115 
1116  int size[3];
1117  float scale;
1118  short minmax[2];
1119  float fCenter[3];
1120  int iCenter[3];
1121  // modality image size
1122  kModality.getSize(size);
1123  ofile.write((char *)size, 3*sizeof(int));
1124  if(DEBUG || kVerbose > 0) {
1125  G4cout << "Modality image size : ("
1126  << size[0] << ", "
1127  << size[1] << ", "
1128  << size[2] << ")"
1129  << G4endl;
1130  }
1131 
1132  // modality image max. & min.
1133  kModality.getMinMax(minmax);
1134  ofile.write((char *)minmax, 4);
1135 
1136  // modality image unit
1137  char munit[13] = "g/cm3 ";
1138  ofile.write((char *)munit, 12);
1139 
1140  // modality image scale
1141  scale = (float)kModality.getScale();
1142  ofile.write((char *)&scale, 4);
1143  if(DEBUG || kVerbose > 0) {
1144  G4cout << "Modality image min., max., scale : "
1145  << minmax[0] << ", "
1146  << minmax[1] << ", "
1147  << scale << G4endl;
1148  }
1149 
1150  // modality image
1151  int psize = size[0]*size[1];
1152  if(DEBUG || kVerbose > 0) G4cout << "Modality image : ";
1153  for(int i = 0; i < size[2]; i++) {
1154  short * image = kModality.getImage(i);
1155  ofile.write((char *)image, psize*sizeof(short));
1156 
1157  if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << image[(size_t)(psize*0.55)] << ", ";
1158  }
1159  if(DEBUG || kVerbose > 0) G4cout << G4endl;
1160 
1161  // modality desity map for CT value
1162  size_t msize = minmax[1] - minmax[0]+1;
1163  float * pdmap = new float[msize];
1164  for(int i = 0; i < (int)msize; i++) pdmap[i] =kModalityImageDensityMap[i];
1165  ofile.write((char *)pdmap, msize*sizeof(float));
1166  if(DEBUG || kVerbose > 0) {
1167  G4cout << "density map : " << std::ends;
1168  for(int i = 0; i < (int)msize; i+=50)
1169  G4cout <<kModalityImageDensityMap[i] << ", ";
1170  G4cout << G4endl;
1171  }
1172  delete [] pdmap;
1173 
1174 
1175  //----- dose distribution image -----//
1176 
1177  if(!isDoseEmpty()) {
1178 
1180 
1181  for(int ndose = 0; ndose < nDoseDist; ndose++) {
1182  // dose distrbution image size
1183  kDose[ndose].getSize(size);
1184  ofile.write((char *)size, 3*sizeof(int));
1185  if(DEBUG || kVerbose > 0) {
1186  G4cout << "Dose dist. [" << ndose << "] image size : ("
1187  << size[0] << ", "
1188  << size[1] << ", "
1189  << size[2] << ")"
1190  << G4endl;
1191  }
1192 
1193  // dose distribution max. & min.
1194  getShortDoseDistMinMax(minmax, ndose);
1195  ofile.write((char *)minmax, 2*2); // sizeof(shorft)*2
1196 
1197  // dose distribution unit
1198  ofile.write((char *)kDoseUnit.c_str(), 12);
1199  if(DEBUG || kVerbose > 0) {
1200  G4cout << "Dose dist. unit : " << kDoseUnit << G4endl;
1201  }
1202 
1203  // dose distribution scaling
1204  double dscale;
1205  dscale = getDoseDistScale(ndose);
1206  scale = float(dscale);
1207  ofile.write((char *)&scale, 4);
1208  if(DEBUG || kVerbose > 0) {
1209  G4cout << "Dose dist. [" << ndose
1210  << "] image min., max., scale : "
1211  << minmax[0] << ", "
1212  << minmax[1] << ", "
1213  << scale << G4endl;
1214  }
1215 
1216  // dose distribution image
1217  int dsize = size[0]*size[1];
1218  short * dimage = new short[dsize];
1219  for(int z = 0; z < size[2]; z++) {
1220  getShortDoseDist(dimage, z, ndose);
1221  ofile.write((char *)dimage, dsize*2); //sizeof(short)
1222 
1223  if(DEBUG || kVerbose > 0) {
1224  for(int j = 0; j < dsize; j++) {
1225  if(dimage[j] < 0)
1226  G4cout << "[" << j << "," << z << "]"
1227  << dimage[j] << ", ";
1228  }
1229  }
1230  }
1231  if(DEBUG || kVerbose > 0) G4cout << G4endl;
1232  delete [] dimage;
1233 
1234  // relative location of the dose distribution image for
1235  // the modality image
1236  getDoseDistCenterPosition(fCenter, ndose);
1237  for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
1238  ofile.write((char *)iCenter, 3*4); // 3*sizeof(int)
1239  if(DEBUG || kVerbose > 0) {
1240  G4cout << "Dose dist. [" << ndose
1241  << "]image relative location : ("
1242  << iCenter[0] << ", "
1243  << iCenter[1] << ", "
1244  << iCenter[2] << ")" << G4endl;
1245  }
1246  }
1247  }
1248 
1249  //----- ROI image -----//
1250  if(!isROIEmpty()) {
1251  // ROI image size
1252  kRoi[0].getSize(size);
1253  ofile.write((char *)size, 3*sizeof(int));
1254  if(DEBUG || kVerbose > 0) {
1255  G4cout << "ROI image size : ("
1256  << size[0] << ", "
1257  << size[1] << ", "
1258  << size[2] << ")"
1259  << G4endl;
1260  }
1261 
1262  // ROI max. & min.
1263  kRoi[0].getMinMax(minmax);
1264  ofile.write((char *)minmax, sizeof(short)*2);
1265 
1266  // ROI distribution scaling
1267  scale = (float)kRoi[0].getScale();
1268  ofile.write((char *)&scale, sizeof(float));
1269  if(DEBUG || kVerbose > 0) {
1270  G4cout << "ROI image min., max., scale : "
1271  << minmax[0] << ", "
1272  << minmax[1] << ", "
1273  << scale << G4endl;
1274  }
1275 
1276  // ROI image
1277  int rsize = size[0]*size[1];
1278  for(int i = 0; i < size[2]; i++) {
1279  short * rimage = kRoi[0].getImage(i);
1280  ofile.write((char *)rimage, rsize*sizeof(short));
1281 
1282  }
1283 
1284  // ROI relative location
1285  kRoi[0].getCenterPosition(fCenter);
1286  for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
1287  ofile.write((char *)iCenter, 3*sizeof(int));
1288  if(DEBUG || kVerbose > 0) {
1289  G4cout << "ROI image relative location : ("
1290  << iCenter[0] << ", "
1291  << iCenter[1] << ", "
1292  << iCenter[2] << ")" << G4endl;
1293  }
1294  }
1295 
1296  //----- track information -----//
1297  // number of track
1298  int ntrk = kSteps.size();
1299  ofile.write((char *)&ntrk, sizeof(int));
1300  if(DEBUG || kVerbose > 0) {
1301  G4cout << "# of tracks : "
1302  << ntrk << G4endl;
1303  }
1304  // track position
1305  for(int i = 0; i < ntrk; i++) {
1306  float * tp = kSteps[i];
1307  ofile.write((char *)tp, sizeof(float)*6);
1308  }
1309  // track color
1310  int ntcolor = int(kStepColors.size());
1311  if(ntrk != ntcolor)
1313  G4cout << "# of track color information must be the same as # of tracks."
1314  << G4endl;
1315  unsigned char white[3] = {255,255,255}; // default color
1316  for(int i = 0; i < ntrk; i++) {
1317  if(i < ntcolor) {
1318  unsigned char * tcolor = kStepColors[i];
1319  ofile.write((char *)tcolor, 3);
1320  } else {
1321  ofile.write((char *)white, 3);
1322  }
1323  }
1324 
1325  // file end mark
1326  ofile.write("END", 3);
1327 
1328  ofile.close();
1329 
1330  return true;
1331 }
1332 //
1333 bool G4GMocrenIO::storeData4(char * _filename) {
1334  kFileName = _filename;
1335  return storeData4();
1336 }
1337 
1338 // version 2
1340 
1341  if(kVerbose > 0) G4cout << ">>>>>>> store data (ver.2) <<<<<<<" << G4endl;
1342  if(kVerbose > 0) G4cout << " " << kFileName << G4endl;
1343 
1344  bool DEBUG = false;//
1345 
1346  // output file open
1347  std::ofstream ofile(kFileName.c_str(),
1348  std::ios_base::out|std::ios_base::binary);
1349 
1350  // file identifier
1351  ofile.write("GRAPE ", 8);
1352 
1353  // file version
1354  unsigned char ver = 0x02;
1355  ofile.write((char *)&ver, 1);
1356  // file id for old file format support
1357  ofile.write(kId.c_str(), IDLENGTH);
1358  // file version for old file format support
1359  ofile.write(kVersion.c_str(), VERLENGTH);
1360  // endian
1361  ofile.write((char *)&kLittleEndianOutput, sizeof(char));
1362 
1363  /*
1364  // event number
1365  ofile.write((char *)&numberOfEvents, sizeof(int));
1366  float imageSpacing[3];
1367  imageSpacing[0] = modalityImageVoxelSpacing[0];
1368  imageSpacing[1] = modalityImageVoxelSpacing[1];
1369  imageSpacing[2] = modalityImageVoxelSpacing[2];
1370  ofile.write((char *)imageSpacing, 12);
1371  */
1372 
1373 
1374  // voxel spacings for all images
1375  ofile.write((char *)kVoxelSpacing, 12);
1376  if(DEBUG || kVerbose > 0) {
1377  G4cout << "Voxel spacing : ("
1378  << kVoxelSpacing[0] << ", "
1379  << kVoxelSpacing[1] << ", "
1380  << kVoxelSpacing[2]
1381  << ") mm " << G4endl;
1382  }
1383 
1384  calcPointers2();
1385  // offset from file starting point to the modality image data
1386  ofile.write((char *)&kPointerToModalityData, 4);
1387 
1388  // offset from file starting point to the dose image data
1389  ofile.write((char *)&kPointerToDoseDistData[0], 4);
1390 
1391  // offset from file starting point to the ROI image data
1392  ofile.write((char *)&kPointerToROIData, 4);
1393 
1394  // offset from file starting point to the track data
1395  ofile.write((char *)&kPointerToTrackData, 4);
1396  if(DEBUG || kVerbose > 0) {
1397  G4cout << "Each pointer to data : "
1398  << kPointerToModalityData << ", "
1399  << kPointerToDoseDistData[0] << ", "
1400  << kPointerToROIData << ", "
1402  }
1403 
1404  //----- modality image -----//
1405 
1406  int size[3];
1407  float scale;
1408  short minmax[2];
1409  float fCenter[3];
1410  int iCenter[3];
1411  // modality image size
1412  kModality.getSize(size);
1413  ofile.write((char *)size, 3*sizeof(int));
1414  if(DEBUG || kVerbose > 0) {
1415  G4cout << "Modality image size : ("
1416  << size[0] << ", "
1417  << size[1] << ", "
1418  << size[2] << ")"
1419  << G4endl;
1420  }
1421 
1422  // modality image max. & min.
1423  kModality.getMinMax(minmax);
1424  ofile.write((char *)minmax, 4);
1425 
1426  // modality image unit
1427  //char munit[13] = "g/cm3 ";
1428  //ofile.write((char *)&munit, 12);
1429 
1430  // modality image scale
1431  scale = (float)kModality.getScale();
1432  ofile.write((char *)&scale, 4);
1433  if(DEBUG || kVerbose > 0) {
1434  G4cout << "Modality image min., max., scale : "
1435  << minmax[0] << ", "
1436  << minmax[1] << ", "
1437  << scale << G4endl;
1438  }
1439 
1440  // modality image
1441  int psize = size[0]*size[1];
1442  if(DEBUG || kVerbose > 0) G4cout << "Modality image : ";
1443  for(int i = 0; i < size[2]; i++) {
1444  short * image =kModality.getImage(i);
1445  ofile.write((char *)image, psize*sizeof(short));
1446 
1447  if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << image[(size_t)(psize*0.55)] << ", ";
1448  }
1449  if(DEBUG || kVerbose > 0) G4cout << G4endl;
1450 
1451  // modality desity map for CT value
1452  size_t msize = minmax[1] - minmax[0]+1;
1453  float * pdmap = new float[msize];
1454  for(int i = 0; i < (int)msize; i++) pdmap[i] =kModalityImageDensityMap[i];
1455  ofile.write((char *)pdmap, msize*sizeof(float));
1456  if(DEBUG || kVerbose > 0) {
1457  G4cout << "density map : " << std::ends;
1458  for(int i = 0; i < (int)msize; i+=50)
1459  G4cout <<kModalityImageDensityMap[i] << ", ";
1460  G4cout << G4endl;
1461  }
1462  delete [] pdmap;
1463 
1464 
1465  //----- dose distribution image -----//
1466 
1467  if(!isDoseEmpty()) {
1469 
1470  // dose distrbution image size
1471  kDose[0].getSize(size);
1472  ofile.write((char *)size, 3*sizeof(int));
1473  if(DEBUG || kVerbose > 0) {
1474  G4cout << "Dose dist. image size : ("
1475  << size[0] << ", "
1476  << size[1] << ", "
1477  << size[2] << ")"
1478  << G4endl;
1479  }
1480 
1481  // dose distribution max. & min.
1482  getShortDoseDistMinMax(minmax);
1483  ofile.write((char *)minmax, sizeof(short)*2);
1484 
1485  // dose distribution scaling
1486  scale = (float)kDose[0].getScale();
1487  ofile.write((char *)&scale, sizeof(float));
1488  if(DEBUG || kVerbose > 0) {
1489  G4cout << "Dose dist. image min., max., scale : "
1490  << minmax[0] << ", "
1491  << minmax[1] << ", "
1492  << scale << G4endl;
1493  }
1494 
1495  // dose distribution image
1496  int dsize = size[0]*size[1];
1497  short * dimage = new short[dsize];
1498  for(int z = 0; z < size[2]; z++) {
1499  getShortDoseDist(dimage, z);
1500  ofile.write((char *)dimage, dsize*sizeof(short));
1501 
1502  if(DEBUG || kVerbose > 0) {
1503  for(int j = 0; j < dsize; j++) {
1504  if(dimage[j] < 0)
1505  G4cout << "[" << j << "," << z << "]"
1506  << dimage[j] << ", ";
1507  }
1508  }
1509  }
1510  if(DEBUG || kVerbose > 0) G4cout << G4endl;
1511  delete [] dimage;
1512 
1513  // relative location of the dose distribution image for
1514  // the modality image
1515  kDose[0].getCenterPosition(fCenter);
1516  for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
1517  ofile.write((char *)iCenter, 3*sizeof(int));
1518  if(DEBUG || kVerbose > 0) {
1519  G4cout << "Dose dist. image relative location : ("
1520  << iCenter[0] << ", "
1521  << iCenter[1] << ", "
1522  << iCenter[2] << ")" << G4endl;
1523  }
1524 
1525  }
1526 
1527  //----- ROI image -----//
1528  if(!isROIEmpty()) {
1529  // ROI image size
1530  kRoi[0].getSize(size);
1531  ofile.write((char *)size, 3*sizeof(int));
1532  if(DEBUG || kVerbose > 0) {
1533  G4cout << "ROI image size : ("
1534  << size[0] << ", "
1535  << size[1] << ", "
1536  << size[2] << ")"
1537  << G4endl;
1538  }
1539 
1540  // ROI max. & min.
1541  kRoi[0].getMinMax(minmax);
1542  ofile.write((char *)minmax, sizeof(short)*2);
1543 
1544  // ROI distribution scaling
1545  scale = (float)kRoi[0].getScale();
1546  ofile.write((char *)&scale, sizeof(float));
1547  if(DEBUG || kVerbose > 0) {
1548  G4cout << "ROI image min., max., scale : "
1549  << minmax[0] << ", "
1550  << minmax[1] << ", "
1551  << scale << G4endl;
1552  }
1553 
1554  // ROI image
1555  int rsize = size[0]*size[1];
1556  for(int i = 0; i < size[2]; i++) {
1557  short * rimage = kRoi[0].getImage(i);
1558  ofile.write((char *)rimage, rsize*sizeof(short));
1559 
1560  }
1561 
1562  // ROI relative location
1563  kRoi[0].getCenterPosition(fCenter);
1564  for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
1565  ofile.write((char *)iCenter, 3*sizeof(int));
1566  if(DEBUG || kVerbose > 0) {
1567  G4cout << "ROI image relative location : ("
1568  << iCenter[0] << ", "
1569  << iCenter[1] << ", "
1570  << iCenter[2] << ")" << G4endl;
1571  }
1572  }
1573 
1574 
1575  //----- track information -----//
1576  // track
1577  int ntrk = kSteps.size();
1578  ofile.write((char *)&ntrk, sizeof(int));
1579  if(DEBUG || kVerbose > 0) {
1580  G4cout << "# of tracks : "
1581  << ntrk << G4endl;
1582  }
1583  for(int i = 0; i < ntrk; i++) {
1584  float * tp = kSteps[i];
1585  ofile.write((char *)tp, sizeof(float)*6);
1586  }
1587 
1588 
1589  // file end mark
1590  ofile.write("END", 3);
1591 
1592  ofile.close();
1593 
1594  return true;
1595 }
1596 //
1597 bool G4GMocrenIO::storeData2(char * _filename) {
1598  kFileName = _filename;
1599  return storeData();
1600 }
1601 
1603 
1604  // input file open
1605  std::ifstream ifile(kFileName.c_str(), std::ios_base::in|std::ios_base::binary);
1606  if(!ifile) {
1608  G4cout << "Cannot open file: " << kFileName
1609  << " in G4GMocrenIO::retrieveData()." << G4endl;
1610  return false;
1611  }
1612 
1613  // file identifier
1614  char verid[9];
1615  ifile.read((char *)verid, 8);
1616  // file version
1617  unsigned char ver;
1618  ifile.read((char *)&ver, 1);
1619  ifile.close();
1620 
1621  if(std::strncmp(verid, "gMocren", 7) == 0) {
1622  if(ver == 0x03) {
1623  G4cout << ">>>>>>> retrieve data (ver.3) <<<<<<<" << G4endl;
1624  G4cout << " " << kFileName << G4endl;
1625  retrieveData3();
1626  } else if (ver == 0x04) {
1627  G4cout << ">>>>>>> retrieve data (ver.4) <<<<<<<" << G4endl;
1628  G4cout << " " << kFileName << G4endl;
1629  retrieveData4();
1630  } else {
1632  G4cout << "Error -- invalid file version : " << (int)ver
1633  << G4endl;
1634  G4cout << " " << kFileName << G4endl;
1635  }
1636  G4Exception("G4GMocrenIO::retrieveDadta()",
1637  "gMocren2001", FatalException,
1638  "Error.");
1639 
1640  }
1641  } else if(std::strncmp(verid, "GRAPE", 5) == 0) {
1642  G4cout << ">>>>>>> retrieve data (ver.2) <<<<<<<" << G4endl;
1643  G4cout << " " << kFileName << G4endl;
1644  retrieveData2();
1645  } else {
1647  G4cout << kFileName << " was not gdd file." << G4endl;
1648  return false;
1649  }
1650 
1651  return true;
1652 }
1653 
1654 bool G4GMocrenIO::retrieveData(char * _filename) {
1655  kFileName = _filename;
1656  return retrieveData();
1657 }
1658 
1659 //
1661 
1662  bool DEBUG = false;//
1663 
1664  // input file open
1665  std::ifstream ifile(kFileName.c_str(), std::ios_base::in|std::ios_base::binary);
1666  if(!ifile) {
1668  G4cout << "Cannot open file: " << kFileName
1669  << " in G4GMocrenIO::retrieveData3()." << G4endl;
1670  return false;
1671  }
1672 
1673  // data buffer
1674  char ctmp[24];
1675 
1676  // file identifier
1677  char verid[9];
1678  ifile.read((char *)verid, 8);
1679 
1680  // file version
1681  unsigned char ver;
1682  ifile.read((char *)&ver, 1);
1683  std::stringstream ss;
1684  ss << (int)ver;
1685  kVersion = ss.str();
1686  if(DEBUG || kVerbose > 0) G4cout << "File version : " << kVersion << G4endl;
1687 
1688  // endian
1689  ifile.read((char *)&kLittleEndianInput, sizeof(char));
1690  if(DEBUG || kVerbose > 0) {
1691  G4cout << "Endian : ";
1692  if(kLittleEndianInput == 1)
1693  G4cout << " little" << G4endl;
1694  else {
1695  G4cout << " big" << G4endl;
1696  }
1697  }
1698 
1699  // comment length (fixed size)
1700  int clength;
1701  ifile.read((char *)ctmp, 4);
1702  convertEndian(ctmp, clength);
1703  // comment
1704  char cmt[1025];
1705  ifile.read((char *)cmt, clength);
1706  std::string scmt = cmt;
1707  scmt += '\0';
1708  setComment(scmt);
1709  if(DEBUG || kVerbose > 0) {
1710  G4cout << "Data comment : "
1711  << kComment << G4endl;
1712  }
1713 
1714  // voxel spacings for all images
1715  ifile.read((char *)ctmp, 12);
1716  convertEndian(ctmp, kVoxelSpacing[0]);
1717  convertEndian(ctmp+4, kVoxelSpacing[1]);
1718  convertEndian(ctmp+8, kVoxelSpacing[2]);
1719  if(DEBUG || kVerbose > 0) {
1720  G4cout << "Voxel spacing : ("
1721  << kVoxelSpacing[0] << ", "
1722  << kVoxelSpacing[1] << ", "
1723  << kVoxelSpacing[2]
1724  << ") mm " << G4endl;
1725  }
1726 
1727 
1728  // offset from file starting point to the modality image data
1729  ifile.read((char *)ctmp, 4);
1731 
1732  // # of dose distributions
1733  ifile.read((char *)ctmp, 4);
1734  int nDoseDist;
1735  convertEndian(ctmp, nDoseDist);
1736 
1737  // offset from file starting point to the dose image data
1738  for(int i = 0; i < nDoseDist; i++) {
1739  ifile.read((char *)ctmp, 4);
1740  unsigned int dptr;
1741  convertEndian(ctmp, dptr);
1743  }
1744 
1745  // offset from file starting point to the ROI image data
1746  ifile.read((char *)ctmp, 4);
1748 
1749  // offset from file starting point to the track data
1750  ifile.read((char *)ctmp, 4);
1752 
1753  // offset from file starting point to the detector data
1754  ifile.read((char *)ctmp, 4);
1756 
1757  if(DEBUG || kVerbose > 0) {
1758  G4cout << "Each pointer to data : "
1759  << kPointerToModalityData << ", ";
1760  for(int i = 0; i < nDoseDist; i++)
1761  G4cout << kPointerToDoseDistData[i] << ", ";
1762  G4cout << kPointerToROIData << ", "
1763  << kPointerToTrackData << ", "
1765  << G4endl;
1766  }
1767 
1768 
1769 
1770  if(kPointerToModalityData == 0 && kPointerToDoseDistData.size() == 0 &&
1771  kPointerToROIData == 0 && kPointerToTrackData == 0) {
1772  if(DEBUG || kVerbose > 0) {
1773  G4cout << "No data." << G4endl;
1774  }
1775  return false;
1776  }
1777 
1778  // event number
1779  /* ver 1
1780  ifile.read(ctmp, sizeof(int));
1781  convertEndian(ctmp, numberOfEvents);
1782  */
1783 
1784  int size[3];
1785  float scale;
1786  double dscale;
1787  short minmax[2];
1788  float fCenter[3];
1789  int iCenter[3];
1790 
1791  //----- Modality image -----//
1792  // modality image size
1793  ifile.read(ctmp, 3*sizeof(int));
1794  convertEndian(ctmp, size[0]);
1795  convertEndian(ctmp+sizeof(int), size[1]);
1796  convertEndian(ctmp+2*sizeof(int), size[2]);
1797  if(DEBUG || kVerbose > 0) {
1798  G4cout << "Modality image size : ("
1799  << size[0] << ", "
1800  << size[1] << ", "
1801  << size[2] << ")"
1802  << G4endl;
1803  }
1804  kModality.setSize(size);
1805 
1806  // modality image voxel spacing
1807  /*
1808  ifile.read(ctmp, 3*sizeof(float));
1809  convertEndian(ctmp, modalityImageVoxelSpacing[0]);
1810  convertEndian(ctmp+sizeof(float), modalityImageVoxelSpacing[1]);
1811  convertEndian(ctmp+2*sizeof(float), modalityImageVoxelSpacing[2]);
1812  */
1813 
1814  if(kPointerToModalityData != 0) {
1815 
1816  // modality density max. & min.
1817  ifile.read((char *)ctmp, 4);
1818  convertEndian(ctmp, minmax[0]);
1819  convertEndian(ctmp+2, minmax[1]);
1820  kModality.setMinMax(minmax);
1821 
1822  // modality image unit
1823  char munit[13];
1824  munit[12] = '\0';
1825  ifile.read((char *)munit, 12);
1826  std::string smunit = munit;
1827  setModalityImageUnit(smunit);
1828 
1829  // modality density scale
1830  ifile.read((char *)ctmp, 4);
1831  convertEndian(ctmp, scale);
1832  kModality.setScale(dscale = scale);
1833  if(DEBUG || kVerbose > 0) {
1834  G4cout << "Modality image min., max., scale : "
1835  << minmax[0] << ", "
1836  << minmax[1] << ", "
1837  << scale << G4endl;
1838  }
1839 
1840  // modality density
1841  int psize = size[0]*size[1];
1842  if(DEBUG || kVerbose > 0) G4cout << "Modality image (" << psize << "): ";
1843  char * cimage = new char[psize*sizeof(short)];
1844  for(int i = 0; i < size[2]; i++) {
1845  ifile.read((char *)cimage, psize*sizeof(short));
1846  short * mimage = new short[psize];
1847  for(int j = 0; j < psize; j++) {
1848  convertEndian(cimage+j*sizeof(short), mimage[j]);
1849  }
1850  kModality.addImage(mimage);
1851 
1852  if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << mimage[(size_t)(psize*0.55)] << ", ";
1853  }
1854  if(DEBUG || kVerbose > 0) G4cout << G4endl;
1855  delete [] cimage;
1856 
1857  // modality desity map for CT value
1858  size_t msize = minmax[1]-minmax[0]+1;
1859  if(DEBUG || kVerbose > 0) G4cout << "msize: " << msize << G4endl;
1860  char * pdmap = new char[msize*sizeof(float)];
1861  ifile.read((char *)pdmap, msize*sizeof(float));
1862  float ftmp;
1863  for(int i = 0; i < (int)msize; i++) {
1864  convertEndian(pdmap+i*sizeof(float), ftmp);
1865  kModalityImageDensityMap.push_back(ftmp);
1866  }
1867  delete [] pdmap;
1868 
1869  if(DEBUG || kVerbose > 0) {
1870  G4cout << "density map : " << std::ends;
1871  for(int i = 0; i < 10; i++)
1872  G4cout <<kModalityImageDensityMap[i] << ", ";
1873  G4cout << G4endl;
1874  for(int i = 0; i < 10; i++) G4cout << "..";
1875  G4cout << G4endl;
1876  for(size_t i =kModalityImageDensityMap.size() - 10; i <kModalityImageDensityMap.size(); i++)
1877  G4cout <<kModalityImageDensityMap[i] << ", ";
1878  G4cout << G4endl;
1879  }
1880 
1881  }
1882 
1883 
1884  //----- dose distribution image -----//
1885  for(int ndose = 0; ndose < nDoseDist; ndose++) {
1886 
1887  newDoseDist();
1888 
1889  // dose distrbution image size
1890  ifile.read((char *)ctmp, 3*sizeof(int));
1891  convertEndian(ctmp, size[0]);
1892  convertEndian(ctmp+sizeof(int), size[1]);
1893  convertEndian(ctmp+2*sizeof(int), size[2]);
1894  if(DEBUG || kVerbose > 0) {
1895  G4cout << "Dose dist. image size : ("
1896  << size[0] << ", "
1897  << size[1] << ", "
1898  << size[2] << ")"
1899  << G4endl;
1900  }
1901  kDose[ndose].setSize(size);
1902 
1903  // dose distribution max. & min.
1904  ifile.read((char *)ctmp, sizeof(short)*2);
1905  convertEndian(ctmp, minmax[0]);
1906  convertEndian(ctmp+2, minmax[1]);
1907 
1908  // dose distribution unit
1909  char dunit[13];
1910  dunit[12] = '\0';
1911  ifile.read((char *)dunit, 12);
1912  std::string sdunit = dunit;
1913  setDoseDistUnit(sdunit, ndose);
1914  if(DEBUG || kVerbose > 0) {
1915  G4cout << "Dose dist. unit : " << kDoseUnit << G4endl;
1916  }
1917 
1918  // dose distribution scaling
1919  ifile.read((char *)ctmp, 4); // sizeof(float)
1920  convertEndian(ctmp, scale);
1921  kDose[ndose].setScale(dscale = scale);
1922 
1923  double dminmax[2];
1924  for(int i = 0; i < 2; i++) dminmax[i] = minmax[i]*dscale;
1925  kDose[ndose].setMinMax(dminmax);
1926 
1927  if(DEBUG || kVerbose > 0) {
1928  G4cout << "Dose dist. image min., max., scale : "
1929  << dminmax[0] << ", "
1930  << dminmax[1] << ", "
1931  << scale << G4endl;
1932  }
1933 
1934  // dose distribution image
1935  int dsize = size[0]*size[1];
1936  if(DEBUG || kVerbose > 0) G4cout << "Dose dist. (" << dsize << "): ";
1937  char * di = new char[dsize*sizeof(short)];
1938  short * shimage = new short[dsize];
1939  for(int z = 0; z < size[2]; z++) {
1940  ifile.read((char *)di, dsize*sizeof(short));
1941  double * dimage = new double[dsize];
1942  for(int xy = 0; xy < dsize; xy++) {
1943  convertEndian(di+xy*sizeof(short), shimage[xy]);
1944  dimage[xy] = shimage[xy]*dscale;
1945  }
1946  kDose[ndose].addImage(dimage);
1947 
1948  if(DEBUG || kVerbose > 0) G4cout << "[" << z << "]" << dimage[(size_t)(dsize*0.55)] << ", ";
1949 
1950  if(DEBUG || kVerbose > 0) {
1951  for(int j = 0; j < dsize; j++) {
1952  if(dimage[j] < 0)
1953  G4cout << "[" << j << "," << z << "]"
1954  << dimage[j] << ", ";
1955  }
1956  }
1957  }
1958  delete [] shimage;
1959  delete [] di;
1960  if(DEBUG || kVerbose > 0) G4cout << G4endl;
1961 
1962  ifile.read((char *)ctmp, 3*4); // 3*sizeof(int)
1963  convertEndian(ctmp, iCenter[0]);
1964  convertEndian(ctmp+4, iCenter[1]);
1965  convertEndian(ctmp+8, iCenter[2]);
1966  for(int i = 0; i < 3; i++) fCenter[i] = (float)iCenter[i];
1967  kDose[ndose].setCenterPosition(fCenter);
1968 
1969  if(DEBUG || kVerbose > 0) {
1970  G4cout << "Dose dist. image relative location : ("
1971  << fCenter[0] << ", "
1972  << fCenter[1] << ", "
1973  << fCenter[2] << ")" << G4endl;
1974  }
1975 
1976 
1977  // dose distribution name
1978  char cname[81];
1979  ifile.read((char *)cname, 80);
1980  std::string dosename = cname;
1981  setDoseDistName(dosename, ndose);
1982  if(DEBUG || kVerbose > 0) {
1983  G4cout << "Dose dist. name : " << dosename << G4endl;
1984  }
1985 
1986  }
1987 
1988  //----- ROI image -----//
1989  if(kPointerToROIData != 0) {
1990 
1991  newROI();
1992 
1993  // ROI image size
1994  ifile.read((char *)ctmp, 3*sizeof(int));
1995  convertEndian(ctmp, size[0]);
1996  convertEndian(ctmp+sizeof(int), size[1]);
1997  convertEndian(ctmp+2*sizeof(int), size[2]);
1998  kRoi[0].setSize(size);
1999  if(DEBUG || kVerbose > 0) {
2000  G4cout << "ROI image size : ("
2001  << size[0] << ", "
2002  << size[1] << ", "
2003  << size[2] << ")"
2004  << G4endl;
2005  }
2006 
2007  // ROI max. & min.
2008  ifile.read((char *)ctmp, sizeof(short)*2);
2009  convertEndian(ctmp, minmax[0]);
2010  convertEndian(ctmp+sizeof(short), minmax[1]);
2011  kRoi[0].setMinMax(minmax);
2012 
2013  // ROI distribution scaling
2014  ifile.read((char *)ctmp, sizeof(float));
2015  convertEndian(ctmp, scale);
2016  kRoi[0].setScale(dscale = scale);
2017  if(DEBUG || kVerbose > 0) {
2018  G4cout << "ROI image min., max., scale : "
2019  << minmax[0] << ", "
2020  << minmax[1] << ", "
2021  << scale << G4endl;
2022  }
2023 
2024  // ROI image
2025  int rsize = size[0]*size[1];
2026  char * ri = new char[rsize*sizeof(short)];
2027  for(int i = 0; i < size[2]; i++) {
2028  ifile.read((char *)ri, rsize*sizeof(short));
2029  short * rimage = new short[rsize];
2030  for(int j = 0; j < rsize; j++) {
2031  convertEndian(ri+j*sizeof(short), rimage[j]);
2032  }
2033  kRoi[0].addImage(rimage);
2034 
2035  }
2036  delete [] ri;
2037 
2038  // ROI relative location
2039  ifile.read((char *)ctmp, 3*sizeof(int));
2040  convertEndian(ctmp, iCenter[0]);
2041  convertEndian(ctmp+sizeof(int), iCenter[1]);
2042  convertEndian(ctmp+2*sizeof(int), iCenter[2]);
2043  for(int i = 0; i < 3; i++) fCenter[i] = iCenter[i];
2044  kRoi[0].setCenterPosition(fCenter);
2045  if(DEBUG || kVerbose > 0) {
2046  G4cout << "ROI image relative location : ("
2047  << fCenter[0] << ", "
2048  << fCenter[1] << ", "
2049  << fCenter[2] << ")" << G4endl;
2050  }
2051 
2052  }
2053 
2054  //----- track information -----//
2055  if(kPointerToTrackData != 0) {
2056 
2057  // track
2058  ifile.read((char *)ctmp, sizeof(int));
2059  int ntrk;
2060  convertEndian(ctmp, ntrk);
2061  if(DEBUG || kVerbose > 0) {
2062  G4cout << "# of tracks: " << ntrk << G4endl;
2063  }
2064 
2065  // track position
2066  unsigned char rgb[3];
2067  for(int i = 0; i < ntrk; i++) {
2068 
2069 
2070  // # of steps in a track
2071  ifile.read((char *)ctmp, sizeof(int));
2072  int nsteps;
2073  convertEndian(ctmp, nsteps);
2074 
2075  // track color
2076  ifile.read((char *)rgb, 3);
2077 
2078  std::vector<float *> steps;
2079  // steps
2080  for(int j = 0; j < nsteps; j++) {
2081 
2082  float * steppoint = new float[6];
2083  ifile.read((char *)ctmp, sizeof(float)*6);
2084 
2085  for(int k = 0; k < 6; k++) {
2086  convertEndian(ctmp+k*sizeof(float), steppoint[k]);
2087  }
2088 
2089  steps.push_back(steppoint);
2090  }
2091 
2092  // add a track to the track container
2093  addTrack(steps, rgb);
2094 
2095  if(DEBUG || kVerbose > 0) {
2096  if(i < 5) {
2097  G4cout << i << ": " ;
2098  for(int j = 0; j < 3; j++) G4cout << steps[0][j] << " ";
2099  int nstp = steps.size();
2100  G4cout << "<-> ";
2101  for(int j = 3; j < 6; j++) G4cout << steps[nstp-1][j] << " ";
2102  G4cout << " rgb( ";
2103  for(int j = 0; j < 3; j++) G4cout << (int)rgb[j] << " ";
2104  G4cout << ")" << G4endl;
2105  }
2106  }
2107  }
2108 
2109 
2110  }
2111 
2112 
2113  //----- detector information -----//
2114  if(kPointerToDetectorData != 0) {
2115 
2116  // number of detectors
2117  ifile.read((char *)ctmp, sizeof(int));
2118  int ndet;
2119  convertEndian(ctmp, ndet);
2120 
2121  if(DEBUG || kVerbose > 0) {
2122  G4cout << "# of detectors : "
2123  << ndet << G4endl;
2124  }
2125 
2126  for(int nd = 0; nd < ndet; nd++) {
2127 
2128  // # of edges of a detector
2129  ifile.read((char *)ctmp, sizeof(int));
2130  int nedges;
2131  convertEndian(ctmp, nedges);
2132  if(DEBUG || kVerbose > 0) {
2133  G4cout << "# of edges in a detector : " << nedges << G4endl;
2134  }
2135 
2136  // edges
2137  std::vector<float *> detector;
2138  char cftmp[24];
2139  for(int ne = 0; ne < nedges; ne++) {
2140 
2141  ifile.read((char *)cftmp, sizeof(float)*6);
2142  float * edgePoints = new float[6];
2143  for(int j = 0; j < 6; j++) convertEndian(&cftmp[sizeof(float)*j], edgePoints[j]);
2144  detector.push_back(edgePoints);
2145 
2146  }
2147 
2148  if(DEBUG || kVerbose > 0) {
2149  G4cout << " first edge : (" << detector[0][0] << ", "
2150  << detector[0][1] << ", "
2151  << detector[0][2] << ") - ("
2152  << detector[0][3] << ", "
2153  << detector[0][4] << ", "
2154  << detector[0][5] << ")" << G4endl;
2155  }
2156 
2157  // detector color
2158  unsigned char dcolor[3];
2159  ifile.read((char *)dcolor, 3);
2160  if(DEBUG || kVerbose > 0) {
2161  G4cout << " detector color : rgb("
2162  << (int)dcolor[0] << ", "
2163  << (int)dcolor[1] << ", "
2164  << (int)dcolor[2] << G4endl;
2165  }
2166 
2167 
2168  // detector name
2169  char cname[80];
2170  ifile.read((char *)cname, 80);
2171  std::string dname = cname;
2172  if(DEBUG || kVerbose > 0) {
2173  G4cout << " detector name : " << dname << G4endl;
2174  }
2175 
2176 
2177  addDetector(dname, detector, dcolor);
2178 
2179  }
2180  }
2181 
2182 
2183  ifile.close();
2184 
2185  return true;
2186 }
2187 bool G4GMocrenIO::retrieveData4(char * _filename) {
2188  kFileName = _filename;
2189  return retrieveData();
2190 }
2191 
2192 //
2194 
2195  bool DEBUG = false;//
2196 
2197  // input file open
2198  std::ifstream ifile(kFileName.c_str(), std::ios_base::in|std::ios_base::binary);
2199  if(!ifile) {
2201  G4cout << "Cannot open file: " << kFileName
2202  << " in G4GMocrenIO::retrieveData3()." << G4endl;
2203  return false;
2204  }
2205 
2206  // data buffer
2207  char ctmp[12];
2208 
2209  // file identifier
2210  char verid[9];
2211  ifile.read((char *)verid, 8);
2212 
2213  // file version
2214  unsigned char ver;
2215  ifile.read((char *)&ver, 1);
2216  std::stringstream ss;
2217  ss << (int)ver;
2218  kVersion = ss.str();
2219  if(DEBUG || kVerbose > 0) G4cout << "File version : " << kVersion << G4endl;
2220 
2221  // endian
2222  ifile.read((char *)&kLittleEndianInput, sizeof(char));
2223  if(DEBUG || kVerbose > 0) {
2224  G4cout << "Endian : ";
2225  if(kLittleEndianInput == 1)
2226  G4cout << " little" << G4endl;
2227  else {
2228  G4cout << " big" << G4endl;
2229  }
2230  }
2231 
2232  // comment length (fixed size)
2233  int clength;
2234  ifile.read((char *)ctmp, 4);
2235  convertEndian(ctmp, clength);
2236  // comment
2237  char cmt[1025];
2238  ifile.read((char *)cmt, clength);
2239  std::string scmt = cmt;
2240  setComment(scmt);
2241  if(DEBUG || kVerbose > 0) {
2242  G4cout << "Data comment : "
2243  << kComment << G4endl;
2244  }
2245 
2246  // voxel spacings for all images
2247  ifile.read((char *)ctmp, 12);
2248  convertEndian(ctmp, kVoxelSpacing[0]);
2249  convertEndian(ctmp+4, kVoxelSpacing[1]);
2250  convertEndian(ctmp+8, kVoxelSpacing[2]);
2251  if(DEBUG || kVerbose > 0) {
2252  G4cout << "Voxel spacing : ("
2253  << kVoxelSpacing[0] << ", "
2254  << kVoxelSpacing[1] << ", "
2255  << kVoxelSpacing[2]
2256  << ") mm " << G4endl;
2257  }
2258 
2259 
2260  // offset from file starting point to the modality image data
2261  ifile.read((char *)ctmp, 4);
2263 
2264  // # of dose distributions
2265  ifile.read((char *)ctmp, 4);
2266  int nDoseDist;
2267  convertEndian(ctmp, nDoseDist);
2268 
2269  // offset from file starting point to the dose image data
2270  for(int i = 0; i < nDoseDist; i++) {
2271  ifile.read((char *)ctmp, 4);
2272  unsigned int dptr;
2273  convertEndian(ctmp, dptr);
2275  }
2276 
2277  // offset from file starting point to the ROI image data
2278  ifile.read((char *)ctmp, 4);
2280 
2281  // offset from file starting point to the track data
2282  ifile.read((char *)ctmp, 4);
2284  if(DEBUG || kVerbose > 0) {
2285  G4cout << "Each pointer to data : "
2286  << kPointerToModalityData << ", ";
2287  for(int i = 0; i < nDoseDist; i++)
2288  G4cout << kPointerToDoseDistData[0] << ", ";
2289  G4cout << kPointerToROIData << ", "
2291  }
2292 
2293  if(kPointerToModalityData == 0 && kPointerToDoseDistData.size() == 0 &&
2294  kPointerToROIData == 0 && kPointerToTrackData == 0) {
2295  if(DEBUG || kVerbose > 0) {
2296  G4cout << "No data." << G4endl;
2297  }
2298  return false;
2299  }
2300 
2301  // event number
2302  /* ver 1
2303  ifile.read(ctmp, sizeof(int));
2304  convertEndian(ctmp, numberOfEvents);
2305  */
2306 
2307  int size[3];
2308  float scale;
2309  double dscale;
2310  short minmax[2];
2311  float fCenter[3];
2312  int iCenter[3];
2313 
2314  //----- Modality image -----//
2315  // modality image size
2316  ifile.read(ctmp, 3*sizeof(int));
2317  convertEndian(ctmp, size[0]);
2318  convertEndian(ctmp+sizeof(int), size[1]);
2319  convertEndian(ctmp+2*sizeof(int), size[2]);
2320  if(DEBUG || kVerbose > 0) {
2321  G4cout << "Modality image size : ("
2322  << size[0] << ", "
2323  << size[1] << ", "
2324  << size[2] << ")"
2325  << G4endl;
2326  }
2327  kModality.setSize(size);
2328 
2329  // modality image voxel spacing
2330  /*
2331  ifile.read(ctmp, 3*sizeof(float));
2332  convertEndian(ctmp, modalityImageVoxelSpacing[0]);
2333  convertEndian(ctmp+sizeof(float), modalityImageVoxelSpacing[1]);
2334  convertEndian(ctmp+2*sizeof(float), modalityImageVoxelSpacing[2]);
2335  */
2336 
2337  if(kPointerToModalityData != 0) {
2338 
2339  // modality density max. & min.
2340  ifile.read((char *)ctmp, 4);
2341  convertEndian(ctmp, minmax[0]);
2342  convertEndian(ctmp+2, minmax[1]);
2343  kModality.setMinMax(minmax);
2344 
2345  // modality image unit
2346  char munit[13];
2347  ifile.read((char *)munit, 12);
2348  std::string smunit = munit;
2349  setModalityImageUnit(smunit);
2350 
2351  // modality density scale
2352  ifile.read((char *)ctmp, 4);
2353  convertEndian(ctmp, scale);
2354  kModality.setScale(dscale = scale);
2355  if(DEBUG || kVerbose > 0) {
2356  G4cout << "Modality image min., max., scale : "
2357  << minmax[0] << ", "
2358  << minmax[1] << ", "
2359  << scale << G4endl;
2360  }
2361 
2362  // modality density
2363  int psize = size[0]*size[1];
2364  if(DEBUG || kVerbose > 0) G4cout << "Modality image (" << psize << "): ";
2365  char * cimage = new char[psize*sizeof(short)];
2366  for(int i = 0; i < size[2]; i++) {
2367  ifile.read((char *)cimage, psize*sizeof(short));
2368  short * mimage = new short[psize];
2369  for(int j = 0; j < psize; j++) {
2370  convertEndian(cimage+j*sizeof(short), mimage[j]);
2371  }
2372  kModality.addImage(mimage);
2373 
2374  if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << mimage[(size_t)(psize*0.55)] << ", ";
2375  }
2376  if(DEBUG || kVerbose > 0) G4cout << G4endl;
2377  delete [] cimage;
2378 
2379  // modality desity map for CT value
2380  size_t msize = minmax[1]-minmax[0]+1;
2381  if(DEBUG || kVerbose > 0) G4cout << "msize: " << msize << G4endl;
2382  char * pdmap = new char[msize*sizeof(float)];
2383  ifile.read((char *)pdmap, msize*sizeof(float));
2384  float ftmp;
2385  for(int i = 0; i < (int)msize; i++) {
2386  convertEndian(pdmap+i*sizeof(float), ftmp);
2387  kModalityImageDensityMap.push_back(ftmp);
2388  }
2389  delete [] pdmap;
2390  if(DEBUG || kVerbose > 0) {
2391  G4cout << "density map : " << std::ends;
2392  for(int i = 0; i < 10; i++)
2393  G4cout <<kModalityImageDensityMap[i] << ", ";
2394  G4cout << G4endl;
2395  for(int i = 0; i < 10; i++) G4cout << "..";
2396  G4cout << G4endl;
2397  for(size_t i =kModalityImageDensityMap.size() - 10; i <kModalityImageDensityMap.size(); i++)
2398  G4cout <<kModalityImageDensityMap[i] << ", ";
2399  G4cout << G4endl;
2400  }
2401 
2402  }
2403 
2404 
2405  //----- dose distribution image -----//
2406  for(int ndose = 0; ndose < nDoseDist; ndose++) {
2407 
2408  newDoseDist();
2409 
2410  // dose distrbution image size
2411  ifile.read((char *)ctmp, 3*sizeof(int));
2412  convertEndian(ctmp, size[0]);
2413  convertEndian(ctmp+sizeof(int), size[1]);
2414  convertEndian(ctmp+2*sizeof(int), size[2]);
2415  if(DEBUG || kVerbose > 0) {
2416  G4cout << "Dose dist. image size : ("
2417  << size[0] << ", "
2418  << size[1] << ", "
2419  << size[2] << ")"
2420  << G4endl;
2421  }
2422  kDose[ndose].setSize(size);
2423 
2424  // dose distribution max. & min.
2425  ifile.read((char *)ctmp, sizeof(short)*2);
2426  convertEndian(ctmp, minmax[0]);
2427  convertEndian(ctmp+2, minmax[1]);
2428 
2429  // dose distribution unit
2430  char dunit[13];
2431  ifile.read((char *)dunit, 12);
2432  std::string sdunit = dunit;
2433  setDoseDistUnit(sdunit, ndose);
2434  if(DEBUG || kVerbose > 0) {
2435  G4cout << "Dose dist. unit : " << kDoseUnit << G4endl;
2436  }
2437 
2438  // dose distribution scaling
2439  ifile.read((char *)ctmp, 4); // sizeof(float)
2440  convertEndian(ctmp, scale);
2441  kDose[ndose].setScale(dscale = scale);
2442 
2443  double dminmax[2];
2444  for(int i = 0; i < 2; i++) dminmax[i] = minmax[i]*dscale;
2445  kDose[ndose].setMinMax(dminmax);
2446 
2447  if(DEBUG || kVerbose > 0) {
2448  G4cout << "Dose dist. image min., max., scale : "
2449  << dminmax[0] << ", "
2450  << dminmax[1] << ", "
2451  << scale << G4endl;
2452  }
2453 
2454  // dose distribution image
2455  int dsize = size[0]*size[1];
2456  if(DEBUG || kVerbose > 0) G4cout << "Dose dist. (" << dsize << "): ";
2457  char * di = new char[dsize*sizeof(short)];
2458  short * shimage = new short[dsize];
2459  for(int z = 0; z < size[2]; z++) {
2460  ifile.read((char *)di, dsize*sizeof(short));
2461  double * dimage = new double[dsize];
2462  for(int xy = 0; xy < dsize; xy++) {
2463  convertEndian(di+xy*sizeof(short), shimage[xy]);
2464  dimage[xy] = shimage[xy]*dscale;
2465  }
2466  kDose[ndose].addImage(dimage);
2467 
2468  if(DEBUG || kVerbose > 0) G4cout << "[" << z << "]" << dimage[(size_t)(dsize*0.55)] << ", ";
2469 
2470  if(DEBUG || kVerbose > 0) {
2471  for(int j = 0; j < dsize; j++) {
2472  if(dimage[j] < 0)
2473  G4cout << "[" << j << "," << z << "]"
2474  << dimage[j] << ", ";
2475  }
2476  }
2477  }
2478  delete [] shimage;
2479  delete [] di;
2480  if(DEBUG || kVerbose > 0) G4cout << G4endl;
2481 
2482  ifile.read((char *)ctmp, 3*4); // 3*sizeof(int)
2483  convertEndian(ctmp, iCenter[0]);
2484  convertEndian(ctmp+4, iCenter[1]);
2485  convertEndian(ctmp+8, iCenter[2]);
2486  for(int i = 0; i < 3; i++) fCenter[i] = (float)iCenter[i];
2487  kDose[ndose].setCenterPosition(fCenter);
2488 
2489  if(DEBUG || kVerbose > 0) {
2490  G4cout << "Dose dist. image relative location : ("
2491  << fCenter[0] << ", "
2492  << fCenter[1] << ", "
2493  << fCenter[2] << ")" << G4endl;
2494  }
2495 
2496 
2497  }
2498 
2499  //----- ROI image -----//
2500  if(kPointerToROIData != 0) {
2501 
2502  newROI();
2503 
2504  // ROI image size
2505  ifile.read((char *)ctmp, 3*sizeof(int));
2506  convertEndian(ctmp, size[0]);
2507  convertEndian(ctmp+sizeof(int), size[1]);
2508  convertEndian(ctmp+2*sizeof(int), size[2]);
2509  kRoi[0].setSize(size);
2510  if(DEBUG || kVerbose > 0) {
2511  G4cout << "ROI image size : ("
2512  << size[0] << ", "
2513  << size[1] << ", "
2514  << size[2] << ")"
2515  << G4endl;
2516  }
2517 
2518  // ROI max. & min.
2519  ifile.read((char *)ctmp, sizeof(short)*2);
2520  convertEndian(ctmp, minmax[0]);
2521  convertEndian(ctmp+sizeof(short), minmax[1]);
2522  kRoi[0].setMinMax(minmax);
2523 
2524  // ROI distribution scaling
2525  ifile.read((char *)ctmp, sizeof(float));
2526  convertEndian(ctmp, scale);
2527  kRoi[0].setScale(dscale = scale);
2528  if(DEBUG || kVerbose > 0) {
2529  G4cout << "ROI image min., max., scale : "
2530  << minmax[0] << ", "
2531  << minmax[1] << ", "
2532  << scale << G4endl;
2533  }
2534 
2535  // ROI image
2536  int rsize = size[0]*size[1];
2537  char * ri = new char[rsize*sizeof(short)];
2538  for(int i = 0; i < size[2]; i++) {
2539  ifile.read((char *)ri, rsize*sizeof(short));
2540  short * rimage = new short[rsize];
2541  for(int j = 0; j < rsize; j++) {
2542  convertEndian(ri+j*sizeof(short), rimage[j]);
2543  }
2544  kRoi[0].addImage(rimage);
2545 
2546  }
2547  delete [] ri;
2548 
2549  // ROI relative location
2550  ifile.read((char *)ctmp, 3*sizeof(int));
2551  convertEndian(ctmp, iCenter[0]);
2552  convertEndian(ctmp+sizeof(int), iCenter[1]);
2553  convertEndian(ctmp+2*sizeof(int), iCenter[2]);
2554  for(int i = 0; i < 3; i++) fCenter[i] = iCenter[i];
2555  kRoi[0].setCenterPosition(fCenter);
2556  if(DEBUG || kVerbose > 0) {
2557  G4cout << "ROI image relative location : ("
2558  << fCenter[0] << ", "
2559  << fCenter[1] << ", "
2560  << fCenter[2] << ")" << G4endl;
2561  }
2562 
2563  }
2564 
2565  //----- track information -----//
2566  if(kPointerToTrackData != 0) {
2567 
2568  // track
2569  ifile.read((char *)ctmp, sizeof(int));
2570  int ntrk;
2571  convertEndian(ctmp, ntrk);
2572  if(DEBUG || kVerbose > 0) {
2573  G4cout << "# of tracks: " << ntrk << G4endl;
2574  }
2575 
2576  // v4
2577  std::vector<float *> trkv4;
2578 
2579  // track position
2580  for(int i = 0; i < ntrk; i++) {
2581  float * tp = new float[6];
2582 
2583  ifile.read((char *)ctmp, sizeof(float)*3);
2584  if(DEBUG || kVerbose > 0) if(i < 10) G4cout << i << ": " ;
2585  for(int j = 0; j < 3; j++) {
2586  convertEndian(ctmp+j*sizeof(float), tp[j]);
2587  if(DEBUG || kVerbose > 0) if(i < 10) G4cout << tp[j] << ", ";
2588  }
2589 
2590  ifile.read((char *)ctmp, sizeof(float)*3);
2591  for(int j = 0; j < 3; j++) {
2592  convertEndian(ctmp+j*sizeof(float), tp[j+3]);
2593  if(DEBUG || kVerbose > 0) if(i < 10) G4cout << tp[j+3] << ", ";
2594  }
2595  addTrack(tp);
2596  if(DEBUG || kVerbose > 0) if(i < 10) G4cout << G4endl;
2597 
2598  // v4
2599  trkv4.push_back(tp);
2600  }
2601 
2602  //v4
2603  unsigned char trkcolorv4[3];
2604 
2605  // track color
2606  for(int i = 0; i < ntrk; i++) {
2607  unsigned char * rgb = new unsigned char[3];
2608  ifile.read((char *)rgb, 3);
2609  addTrackColor(rgb);
2610 
2611  // v4
2612  for(int j = 0; j < 3; j++) trkcolorv4[j] = rgb[j];
2613  std::vector<float *> trk;
2614  trk.push_back(trkv4[i]);
2615  addTrack(trk, trkcolorv4);
2616 
2617  }
2618 
2619  }
2620 
2621  ifile.close();
2622 
2623  return true;
2624 }
2625 bool G4GMocrenIO::retrieveData3(char * _filename) {
2626  kFileName = _filename;
2627  return retrieveData();
2628 }
2629 
2630 //
2632 
2633  bool DEBUG = false;//
2634 
2635  // input file open
2636  std::ifstream ifile(kFileName.c_str(), std::ios_base::in|std::ios_base::binary);
2637  if(!ifile) {
2639  G4cout << "Cannot open file: " << kFileName
2640  << " in G4GMocrenIO::retrieveData2()." << G4endl;
2641  return false;
2642  }
2643 
2644  // data buffer
2645  char ctmp[12];
2646 
2647  // file identifier
2648  char verid[9];
2649  ifile.read((char *)verid, 8);
2650 
2651  // file version
2652  unsigned char ver;
2653  ifile.read((char *)&ver, 1);
2654  std::stringstream ss;
2655  ss << (int)ver;
2656  kVersion = ss.str();
2657  if(DEBUG || kVerbose > 0) G4cout << "File version : " << kVersion << G4endl;
2658 
2659  // id of version 1
2660  char idtmp[IDLENGTH];
2661  ifile.read((char *)idtmp, IDLENGTH);
2662  kId = idtmp;
2663  // version of version 1
2664  char vertmp[VERLENGTH];
2665  ifile.read((char *)vertmp, VERLENGTH);
2666 
2667  // endian
2668  ifile.read((char *)&kLittleEndianInput, sizeof(char));
2669  if(DEBUG || kVerbose > 0) {
2670  G4cout << "Endian : ";
2671  if(kLittleEndianInput == 1)
2672  G4cout << " little" << G4endl;
2673  else {
2674  G4cout << " big" << G4endl;
2675  }
2676  }
2677 
2678  // voxel spacings for all images
2679  ifile.read((char *)ctmp, 12);
2680  convertEndian(ctmp, kVoxelSpacing[0]);
2681  convertEndian(ctmp+4, kVoxelSpacing[1]);
2682  convertEndian(ctmp+8, kVoxelSpacing[2]);
2683  if(DEBUG || kVerbose > 0) {
2684  G4cout << "Voxel spacing : ("
2685  << kVoxelSpacing[0] << ", "
2686  << kVoxelSpacing[1] << ", "
2687  << kVoxelSpacing[2]
2688  << ") mm " << G4endl;
2689  }
2690 
2691 
2692  // offset from file starting point to the modality image data
2693  ifile.read((char *)ctmp, 4);
2695 
2696  // offset from file starting point to the dose image data
2697  unsigned int ptddd;
2698  ifile.read((char *)ctmp, 4);
2699  convertEndian(ctmp, ptddd);
2700  kPointerToDoseDistData.push_back(ptddd);
2701 
2702  // offset from file starting point to the ROI image data
2703  ifile.read((char *)ctmp, 4);
2705 
2706  // offset from file starting point to the track data
2707  ifile.read((char *)ctmp, 4);
2709  if(DEBUG || kVerbose > 0) {
2710  G4cout << "Each pointer to data : "
2711  << kPointerToModalityData << ", "
2712  << kPointerToDoseDistData[0] << ", "
2713  << kPointerToROIData << ", "
2715  }
2716 
2717  if(kPointerToModalityData == 0 && kPointerToDoseDistData.size() == 0 &&
2718  kPointerToROIData == 0 && kPointerToTrackData == 0) {
2719  if(DEBUG || kVerbose > 0) {
2720  G4cout << "No data." << G4endl;
2721  }
2722  return false;
2723  }
2724 
2725  // event number
2726  /* ver 1
2727  ifile.read(ctmp, sizeof(int));
2728  convertEndian(ctmp, numberOfEvents);
2729  */
2730 
2731  int size[3];
2732  float scale;
2733  double dscale;
2734  short minmax[2];
2735  float fCenter[3];
2736  int iCenter[3];
2737 
2738  //----- Modality image -----//
2739  // modality image size
2740  ifile.read(ctmp, 3*sizeof(int));
2741  convertEndian(ctmp, size[0]);
2742  convertEndian(ctmp+sizeof(int), size[1]);
2743  convertEndian(ctmp+2*sizeof(int), size[2]);
2744  if(DEBUG || kVerbose > 0) {
2745  G4cout << "Modality image size : ("
2746  << size[0] << ", "
2747  << size[1] << ", "
2748  << size[2] << ")"
2749  << G4endl;
2750  }
2751  kModality.setSize(size);
2752 
2753  // modality image voxel spacing
2754  /*
2755  ifile.read(ctmp, 3*sizeof(float));
2756  convertEndian(ctmp, modalityImageVoxelSpacing[0]);
2757  convertEndian(ctmp+sizeof(float), modalityImageVoxelSpacing[1]);
2758  convertEndian(ctmp+2*sizeof(float), modalityImageVoxelSpacing[2]);
2759  */
2760 
2761  if(kPointerToModalityData != 0) {
2762 
2763  // modality density max. & min.
2764  ifile.read((char *)ctmp, 4);
2765  convertEndian(ctmp, minmax[0]);
2766  convertEndian(ctmp+2, minmax[1]);
2767  kModality.setMinMax(minmax);
2768 
2769  // modality density scale
2770  ifile.read((char *)ctmp, 4);
2771  convertEndian(ctmp, scale);
2772  kModality.setScale(dscale = scale);
2773  if(DEBUG || kVerbose > 0) {
2774  G4cout << "Modality image min., max., scale : "
2775  << minmax[0] << ", "
2776  << minmax[1] << ", "
2777  << scale << G4endl;
2778  }
2779 
2780  // modality density
2781  int psize = size[0]*size[1];
2782  if(DEBUG || kVerbose > 0) G4cout << "Modality image (" << psize << "): ";
2783  char * cimage = new char[psize*sizeof(short)];
2784  for(int i = 0; i < size[2]; i++) {
2785  ifile.read((char *)cimage, psize*sizeof(short));
2786  short * mimage = new short[psize];
2787  for(int j = 0; j < psize; j++) {
2788  convertEndian(cimage+j*sizeof(short), mimage[j]);
2789  }
2790  kModality.addImage(mimage);
2791 
2792  if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << mimage[(size_t)(psize*0.55)] << ", ";
2793  }
2794  if(DEBUG || kVerbose > 0) G4cout << G4endl;
2795  delete [] cimage;
2796 
2797  // modality desity map for CT value
2798  size_t msize = minmax[1]-minmax[0]+1;
2799  if(DEBUG || kVerbose > 0) G4cout << "msize: " << msize << G4endl;
2800  char * pdmap = new char[msize*sizeof(float)];
2801  ifile.read((char *)pdmap, msize*sizeof(float));
2802  float ftmp;
2803  for(int i = 0; i < (int)msize; i++) {
2804  convertEndian(pdmap+i*sizeof(float), ftmp);
2805  kModalityImageDensityMap.push_back(ftmp);
2806  }
2807  delete [] pdmap;
2808  if(DEBUG || kVerbose > 0) {
2809  G4cout << "density map : " << std::ends;
2810  for(int i = 0; i < 10; i++)
2811  G4cout <<kModalityImageDensityMap[i] << ", ";
2812  G4cout << G4endl;
2813  for(int i = 0; i < 10; i++) G4cout << "..";
2814  G4cout << G4endl;
2815  for(size_t i =kModalityImageDensityMap.size() - 10; i <kModalityImageDensityMap.size(); i++)
2816  G4cout <<kModalityImageDensityMap[i] << ", ";
2817  G4cout << G4endl;
2818  }
2819 
2820  }
2821 
2822 
2823  //----- dose distribution image -----//
2824  if(kPointerToDoseDistData[0] != 0) {
2825 
2826  newDoseDist();
2827 
2828  // dose distrbution image size
2829  ifile.read((char *)ctmp, 3*sizeof(int));
2830  convertEndian(ctmp, size[0]);
2831  convertEndian(ctmp+sizeof(int), size[1]);
2832  convertEndian(ctmp+2*sizeof(int), size[2]);
2833  if(DEBUG || kVerbose > 0) {
2834  G4cout << "Dose dist. image size : ("
2835  << size[0] << ", "
2836  << size[1] << ", "
2837  << size[2] << ")"
2838  << G4endl;
2839  }
2840  kDose[0].setSize(size);
2841 
2842  // dose distribution max. & min.
2843  ifile.read((char *)ctmp, sizeof(short)*2);
2844  convertEndian(ctmp, minmax[0]);
2845  convertEndian(ctmp+2, minmax[1]);
2846  // dose distribution scaling
2847  ifile.read((char *)ctmp, sizeof(float));
2848  convertEndian(ctmp, scale);
2849  kDose[0].setScale(dscale = scale);
2850 
2851  double dminmax[2];
2852  for(int i = 0; i < 2; i++) dminmax[i] = minmax[i]*dscale;
2853  kDose[0].setMinMax(dminmax);
2854 
2855  if(DEBUG || kVerbose > 0) {
2856  G4cout << "Dose dist. image min., max., scale : "
2857  << dminmax[0] << ", "
2858  << dminmax[1] << ", "
2859  << scale << G4endl;
2860  }
2861 
2862  // dose distribution image
2863  int dsize = size[0]*size[1];
2864  if(DEBUG || kVerbose > 0) G4cout << "Dose dist. (" << dsize << "): ";
2865  char * di = new char[dsize*sizeof(short)];
2866  short * shimage = new short[dsize];
2867  for(int z = 0; z < size[2]; z++) {
2868  ifile.read((char *)di, dsize*sizeof(short));
2869  double * dimage = new double[dsize];
2870  for(int xy = 0; xy < dsize; xy++) {
2871  convertEndian(di+xy*sizeof(short), shimage[xy]);
2872  dimage[xy] = shimage[xy]*dscale;
2873  }
2874  kDose[0].addImage(dimage);
2875 
2876  if(DEBUG || kVerbose > 0) G4cout << "[" << z << "]" << dimage[(size_t)(dsize*0.55)] << ", ";
2877 
2878  if(DEBUG || kVerbose > 0) {
2879  for(int j = 0; j < dsize; j++) {
2880  if(dimage[j] < 0)
2881  G4cout << "[" << j << "," << z << "]"
2882  << dimage[j] << ", ";
2883  }
2884  }
2885  }
2886  delete [] shimage;
2887  delete [] di;
2888  if(DEBUG || kVerbose > 0) G4cout << G4endl;
2889 
2890  /* ver 1
2891  float doseDist;
2892  int dosePid;
2893  double * doseData = new double[numDoseImageVoxels];
2894  for(int i = 0; i < numDose; i++) {
2895  ifile.read(ctmp, sizeof(int));
2896  convertEndian(ctmp, dosePid);
2897  for(int j = 0; j < numDoseImageVoxels; j++) {
2898  ifile.read(ctmp, sizeof(float));
2899  convertEndian(ctmp, doseDist);
2900  doseData[j] = doseDist;
2901  }
2902  setDose(dosePid, doseData);
2903  }
2904  delete [] doseData;
2905  if(totalDose == NULL) totalDose = new double[numDoseImageVoxels];
2906  for(int i = 0; i < numDoseImageVoxels; i++) {
2907  ifile.read(ctmp, sizeof(float));
2908  convertEndian(ctmp, doseDist);
2909  totalDose[i] = doseDist;
2910  }
2911  */
2912 
2913  /* ver 1
2914  // relative location between the two images
2915  ifile.read(ctmp, 3*sizeof(float));
2916  convertEndian(ctmp, relativeLocation[0]);
2917  convertEndian(ctmp+sizeof(float), relativeLocation[1]);
2918  convertEndian(ctmp+2*sizeof(float), relativeLocation[2]);
2919  */
2920 
2921  // relative location of the dose distribution image for
2922  // the modality image
2923  //ofile.write((char *)relativeLocation, 3*sizeof(float));
2924  ifile.read((char *)ctmp, 3*sizeof(int));
2925  convertEndian(ctmp, iCenter[0]);
2926  convertEndian(ctmp+sizeof(int), iCenter[1]);
2927  convertEndian(ctmp+2*sizeof(int), iCenter[2]);
2928  for(int i = 0; i < 3; i++) fCenter[i] = (float)iCenter[i];
2929  kDose[0].setCenterPosition(fCenter);
2930 
2931  if(DEBUG || kVerbose > 0) {
2932  G4cout << "Dose dist. image relative location : ("
2933  << fCenter[0] << ", "
2934  << fCenter[1] << ", "
2935  << fCenter[2] << ")" << G4endl;
2936  }
2937 
2938 
2939  }
2940 
2941  //----- ROI image -----//
2942  if(kPointerToROIData != 0) {
2943 
2944  newROI();
2945 
2946  // ROI image size
2947  ifile.read((char *)ctmp, 3*sizeof(int));
2948  convertEndian(ctmp, size[0]);
2949  convertEndian(ctmp+sizeof(int), size[1]);
2950  convertEndian(ctmp+2*sizeof(int), size[2]);
2951  kRoi[0].setSize(size);
2952  if(DEBUG || kVerbose > 0) {
2953  G4cout << "ROI image size : ("
2954  << size[0] << ", "
2955  << size[1] << ", "
2956  << size[2] << ")"
2957  << G4endl;
2958  }
2959 
2960  // ROI max. & min.
2961  ifile.read((char *)ctmp, sizeof(short)*2);
2962  convertEndian(ctmp, minmax[0]);
2963  convertEndian(ctmp+sizeof(short), minmax[1]);
2964  kRoi[0].setMinMax(minmax);
2965 
2966  // ROI distribution scaling
2967  ifile.read((char *)ctmp, sizeof(float));
2968  convertEndian(ctmp, scale);
2969  kRoi[0].setScale(dscale = scale);
2970  if(DEBUG || kVerbose > 0) {
2971  G4cout << "ROI image min., max., scale : "
2972  << minmax[0] << ", "
2973  << minmax[1] << ", "
2974  << scale << G4endl;
2975  }
2976 
2977  // ROI image
2978  int rsize = size[0]*size[1];
2979  char * ri = new char[rsize*sizeof(short)];
2980  for(int i = 0; i < size[2]; i++) {
2981  ifile.read((char *)ri, rsize*sizeof(short));
2982  short * rimage = new short[rsize];
2983  for(int j = 0; j < rsize; j++) {
2984  convertEndian(ri+j*sizeof(short), rimage[j]);
2985  }
2986  kRoi[0].addImage(rimage);
2987 
2988  }
2989  delete [] ri;
2990 
2991  // ROI relative location
2992  ifile.read((char *)ctmp, 3*sizeof(int));
2993  convertEndian(ctmp, iCenter[0]);
2994  convertEndian(ctmp+sizeof(int), iCenter[1]);
2995  convertEndian(ctmp+2*sizeof(int), iCenter[2]);
2996  for(int i = 0; i < 3; i++) fCenter[i] = iCenter[i];
2997  kRoi[0].setCenterPosition(fCenter);
2998  if(DEBUG || kVerbose > 0) {
2999  G4cout << "ROI image relative location : ("
3000  << fCenter[0] << ", "
3001  << fCenter[1] << ", "
3002  << fCenter[2] << ")" << G4endl;
3003  }
3004 
3005  }
3006 
3007  //----- track information -----//
3008  if(kPointerToTrackData != 0) {
3009 
3010  // track
3011  ifile.read((char *)ctmp, sizeof(int));
3012  int ntrk;
3013  convertEndian(ctmp, ntrk);
3014  if(DEBUG || kVerbose > 0) {
3015  G4cout << "# of tracks: " << ntrk << G4endl;
3016  }
3017 
3018  //v4
3019  unsigned char trkcolorv4[3] = {255, 0, 0};
3020 
3021  for(int i = 0; i < ntrk; i++) {
3022  float * tp = new float[6];
3023  // v4
3024  std::vector<float *> trkv4;
3025 
3026  ifile.read((char *)ctmp, sizeof(float)*3);
3027  if(DEBUG || kVerbose > 0) if(i < 10) G4cout << i << ": " ;
3028  for(int j = 0; j < 3; j++) {
3029  convertEndian(ctmp+j*sizeof(float), tp[j]);
3030  if(DEBUG || kVerbose > 0) if(i < 10) G4cout << tp[j] << ", ";
3031  }
3032 
3033  ifile.read((char *)ctmp, sizeof(float)*3);
3034  for(int j = 0; j < 3; j++) {
3035  convertEndian(ctmp+j*sizeof(float), tp[j+3]);
3036  if(DEBUG || kVerbose > 0) if(i < 10) G4cout << tp[j+3] << ", ";
3037  }
3038 
3039  kSteps.push_back(tp);
3040  // v4
3041  trkv4.push_back(tp);
3042  addTrack(trkv4, trkcolorv4);
3043 
3044  if(DEBUG || kVerbose > 0) if(i < 10) G4cout << G4endl;
3045  }
3046 
3047  }
3048 
3049  /* ver 1
3050  // track
3051  int ntracks;
3052  ifile.read(ctmp, sizeof(int));
3053  convertEndian(ctmp, ntracks);
3054  // track displacement
3055  ifile.read(ctmp, 3*sizeof(float));
3056  convertEndian(ctmp, trackDisplacement[0]);
3057  convertEndian(ctmp+sizeof(float), trackDisplacement[2]); // exchanged with [1]
3058  convertEndian(ctmp+2*sizeof(float), trackDisplacement[1]);
3059  //
3060  //for(int i = 0; i < ntracks && i < 100; i++) {
3061  for(int i = 0; i < ntracks; i++) {
3062  DicomDoseTrack trk;
3063  short trackid, parentid, pid;
3064  int npoints;
3065  ifile.read(ctmp, sizeof(short));
3066  convertEndian(ctmp, trackid);
3067  trk.setID(trackid);
3068  ifile.read(ctmp, sizeof(short));
3069  convertEndian(ctmp, parentid);
3070  trk.setParentID(parentid);
3071  ifile.read(ctmp, sizeof(short));
3072  convertEndian(ctmp, pid);
3073  trk.setPID(pid);
3074  ifile.read(ctmp, sizeof(int));
3075  convertEndian(ctmp, npoints);
3076  for(int i = 0; i < npoints; i++) {
3077  ifile.read(ctmp, 3*sizeof(float));
3078  // storing only start and end points
3079  //if(i == 0 || i == npoints - 1) {
3080  float * point = new float[3];
3081  convertEndian(ctmp, point[0]);
3082  convertEndian(ctmp+sizeof(float), point[1]);
3083  convertEndian(ctmp+2*sizeof(float), point[2]);
3084  trk.addPoint(point);
3085  //}
3086  }
3087  track.push_back(trk);
3088  }
3089  */
3090 
3091  ifile.close();
3092 
3093  return true;
3094 }
3095 
3096 bool G4GMocrenIO::retrieveData2(char * _filename) {
3097  kFileName = _filename;
3098  return retrieveData();
3099 }
3100 
3102  time_t t;
3103  time(&t);
3104 
3105  tm * ti;
3106  ti = localtime(&t);
3107 
3108  char cmonth[12][4] = {"Jan", "Feb", "Mar", "Apr",
3109  "May", "Jun", "Jul", "Aug",
3110  "Sep", "Oct", "Nov", "Dec"};
3111  std::stringstream ss;
3112  ss << std::setfill('0')
3113  << std::setw(2)
3114  << ti->tm_hour << ":"
3115  << std::setw(2)
3116  << ti->tm_min << ":"
3117  << std::setw(2)
3118  << ti->tm_sec << ","
3119  << cmonth[ti->tm_mon] << "."
3120  << std::setw(2)
3121  << ti->tm_mday << ","
3122  << ti->tm_year+1900;
3123 
3124  kId = ss.str();
3125 }
3126 
3127 // get & set the file version
3128 std::string & G4GMocrenIO::getVersion() {return kVersion;}
3129 void G4GMocrenIO::setVersion(std::string & _version) {kVersion = _version;}
3130 
3131 // set endians of input/output data
3134 
3135 // voxel spacing
3136 void G4GMocrenIO::setVoxelSpacing(float _spacing[3]) {
3137  for(int i = 0; i < 3; i++) kVoxelSpacing[i] = _spacing[i];
3138 }
3139 void G4GMocrenIO::getVoxelSpacing(float _spacing[3]) {
3140  for(int i = 0; i < 3; i++) _spacing[i] = kVoxelSpacing[i];
3141 }
3142 
3143 // get & set number of events
3145  return kNumberOfEvents;
3146 }
3147 void G4GMocrenIO::setNumberOfEvents(int & _numberOfEvents) {
3148  kNumberOfEvents = _numberOfEvents;
3149 }
3151  kNumberOfEvents++;
3152 }
3153 
3154 // set/get pointer the modality image data
3155 void G4GMocrenIO::setPointerToModalityData(unsigned int & _pointer) {
3156  kPointerToModalityData = _pointer;
3157 }
3159  return kPointerToModalityData;
3160 }
3161 // set/get pointer the dose distribution image data
3162 void G4GMocrenIO::addPointerToDoseDistData(unsigned int & _pointer) {
3163  kPointerToDoseDistData.push_back(_pointer);
3164 }
3165 unsigned int G4GMocrenIO::getPointerToDoseDistData(int _elem) {
3166  if(kPointerToDoseDistData.size() == 0 ||
3167  kPointerToDoseDistData.size() < (size_t)_elem)
3168  return 0;
3169  else
3170  return kPointerToDoseDistData[_elem];
3171 }
3172 
3173 // set/get pointer the ROI image data
3174 void G4GMocrenIO::setPointerToROIData(unsigned int & _pointer) {
3175  kPointerToROIData = _pointer;
3176 }
3178  return kPointerToROIData;
3179 }
3180 // set/get pointer the track data
3181 void G4GMocrenIO::setPointerToTrackData(unsigned int & _pointer) {
3182  kPointerToTrackData = _pointer;
3183 }
3185  return kPointerToTrackData;
3186 }
3187 
3188 // calculate pointers for version 4
3190 
3191  // pointer to modality data
3192  unsigned int pointer = 1070; // up to "pointer to the detector data" except for "pointer to the dose dist data"
3193  int nDoseDist = getNumDoseDist();
3194  pointer += nDoseDist*4;
3195 
3196  setPointerToModalityData(pointer);
3197 
3198  // pointer to dose data
3199  // ct-density map for modality data
3200  int msize[3];
3201  getModalityImageSize(msize);
3202  short mminmax[2];
3203  getModalityImageMinMax(mminmax);
3204  int pmsize = 2*msize[0]*msize[1]*msize[2];
3205  int pmmap = 4*(mminmax[1] - mminmax[0] + 1);
3206  pointer += 32 + pmsize + pmmap;
3207  //
3208  kPointerToDoseDistData.clear();
3209  if(nDoseDist == 0) {
3210  unsigned int pointer0 = 0;
3211  addPointerToDoseDistData(pointer0);
3212  }
3213  for(int ndose = 0; ndose < nDoseDist; ndose++) {
3214  addPointerToDoseDistData(pointer);
3215  int dsize[3];
3216  getDoseDistSize(dsize);
3217  pointer += 44 + dsize[0]*dsize[1]*dsize[2]*2 + 80;
3218  }
3219 
3220  // pointer to roi data
3221  if(!isROIEmpty()) {
3222  setPointerToROIData(pointer);
3223 
3224  int rsize[3];
3225  getROISize(rsize);
3226  int prsize = 2*rsize[0]*rsize[1]*rsize[2];
3227  pointer += 20 + prsize + 12;
3228  } else {
3229  unsigned int pointer0 = 0;
3230  setPointerToROIData(pointer0);
3231  }
3232 
3233  // pointer to track data
3234  int ntrk = kTracks.size();
3235  if(ntrk != 0) {
3236  setPointerToTrackData(pointer);
3237 
3238  pointer += 4; // # of tracks
3239  for(int nt = 0; nt < ntrk; nt++) {
3240  int nsteps = kTracks[nt].getNumberOfSteps();
3241  pointer += 4 + 3 + nsteps*(4*6); // # of steps + color + steps(float*6)
3242  }
3243  } else {
3244  unsigned int pointer0 = 0;
3245  setPointerToTrackData(pointer0);
3246  }
3247  if(kVerbose > 0) G4cout << " pointer to the track data :"
3249 
3250  // pointer to detector data
3251  int ndet = kDetectors.size();
3252  if(ndet != 0) {
3253  kPointerToDetectorData = pointer;
3254  } else {
3256  }
3257  if(kVerbose > 0) G4cout << " pointer to the detector data :"
3259 
3260 }
3261 
3262 // calculate pointers for ver.3
3264 
3265  // pointer to modality data
3266  unsigned int pointer = 1066; // up to "pointer to the track data" except for "pointer to the dose dist data"
3267  int nDoseDist = getNumDoseDist();
3268  pointer += nDoseDist*4;
3269 
3270  setPointerToModalityData(pointer);
3271 
3272  // pointer to dose data
3273  // ct-density map for modality data
3274  int msize[3];
3275  getModalityImageSize(msize);
3276  short mminmax[2];
3277  getModalityImageMinMax(mminmax);
3278  int pmsize = 2*msize[0]*msize[1]*msize[2];
3279  int pmmap = 4*(mminmax[1] - mminmax[0] + 1);
3280  pointer += 32 + pmsize + pmmap;
3281  //
3282  kPointerToDoseDistData.clear();
3283  if(nDoseDist == 0) {
3284  unsigned int pointer0 = 0;
3285  addPointerToDoseDistData(pointer0);
3286  }
3287  for(int ndose = 0; ndose < nDoseDist; ndose++) {
3288  addPointerToDoseDistData(pointer);
3289  int dsize[3];
3290  getDoseDistSize(dsize);
3291  pointer += 44 + dsize[0]*dsize[1]*dsize[2]*2;
3292  }
3293 
3294  // pointer to roi data
3295  if(!isROIEmpty()) {
3296  setPointerToROIData(pointer);
3297 
3298  int rsize[3];
3299  getROISize(rsize);
3300  int prsize = 2*rsize[0]*rsize[1]*rsize[2];
3301  pointer += 20 + prsize + 12;
3302  } else {
3303  unsigned int pointer0 = 0;
3304  setPointerToROIData(pointer0);
3305  }
3306 
3307  //
3308  if(getNumTracks() != 0)
3309  setPointerToTrackData(pointer);
3310  else {
3311  unsigned int pointer0 = 0;
3312  setPointerToTrackData(pointer0);
3313  }
3314 
3315 }
3316 
3317 // calculate pointers for ver.2
3319 
3320  // pointer to modality data
3321  unsigned int pointer = 65;
3322  setPointerToModalityData(pointer);
3323 
3324  // pointer to dose data
3325  int msize[3];
3326  getModalityImageSize(msize);
3327  short mminmax[2];
3328  getModalityImageMinMax(mminmax);
3329  int pmsize = 2*msize[0]*msize[1]*msize[2];
3330  int pmmap = 4*(mminmax[1] - mminmax[0] + 1);
3331  pointer += 20 + pmsize + pmmap;
3332  int dsize[3];
3333  getDoseDistSize(dsize);
3334  kPointerToDoseDistData.clear();
3335  if(dsize[0] != 0) {
3336  kPointerToDoseDistData.push_back(pointer);
3337 
3338  int pdsize = 2*dsize[0]*dsize[1]*dsize[2];
3339  pointer += 20 + pdsize + 12;
3340  } else {
3341  unsigned int pointer0 = 0;
3342  kPointerToDoseDistData.push_back(pointer0);
3343  }
3344 
3345  // pointer to roi data
3346  if(!isROIEmpty()) {
3347  int rsize[3];
3348  getROISize(rsize);
3349  setPointerToROIData(pointer);
3350  int prsize = 2*rsize[0]*rsize[1]*rsize[2];
3351  pointer += 20 + prsize + 12;
3352 
3353  } else {
3354  unsigned int pointer0 = 0;
3355  setPointerToROIData(pointer0);
3356  }
3357 
3358  //
3359  if(getNumTracks() != 0)
3360  setPointerToTrackData(pointer);
3361  else {
3362  unsigned int pointer0 = 0;
3363  setPointerToTrackData(pointer0);
3364  }
3365 
3366 }
3367 
3368 
3369 //----- Modality image -----//
3371 
3372  kModality.getSize(_size);
3373 }
3375 
3376  kModality.setSize(_size);
3377 }
3378 
3379 // get & set the modality image size
3380 void G4GMocrenIO::setModalityImageScale(double & _scale) {
3381 
3382  kModality.setScale(_scale);
3383 }
3385 
3386  return kModality.getScale();
3387 }
3388 
3389 // set the modality image in CT
3390 void G4GMocrenIO::setModalityImage(short * _image) {
3391 
3392  kModality.addImage(_image);
3393 }
3395 
3396  return kModality.getImage(_z);
3397 }
3399 
3401 }
3402 // set/get the modality image density map
3403 void G4GMocrenIO::setModalityImageDensityMap(std::vector<float> & _map) {
3404  kModalityImageDensityMap = _map;
3405 }
3407  return kModalityImageDensityMap;
3408 }
3409 // set the modality image min./max.
3410 void G4GMocrenIO::setModalityImageMinMax(short _minmax[2]) {
3411 
3412  kModality.setMinMax(_minmax);
3413 }
3414 // get the modality image min./max.
3415 void G4GMocrenIO::getModalityImageMinMax(short _minmax[2]) {
3416 
3417  short minmax[2];
3418  kModality.getMinMax(minmax);
3419  for(int i = 0; i < 2; i++) _minmax[i] = minmax[i];
3420 }
3422 
3423  short minmax[2];
3424  kModality.getMinMax(minmax);
3425  return minmax[1];
3426 }
3428 
3429  short minmax[2];
3430  kModality.getMinMax(minmax);
3431  return minmax[0];
3432 }
3433 // set/get position of the modality image center
3435 
3436  kModality.setCenterPosition(_center);
3437 }
3439 
3440  if(isROIEmpty())
3441  for(int i = 0; i < 3; i++) _center[i] = 0;
3442  else
3443  kModality.getCenterPosition(_center);
3444 }
3445 // get & set the modality image unit
3447  return kModalityUnit;
3448 }
3449 void G4GMocrenIO::setModalityImageUnit(std::string & _unit) {
3450  kModalityUnit = _unit;
3451 }
3452 //
3453 short G4GMocrenIO::convertDensityToHU(float & _dens) {
3454  short rval = -1024; // default: air
3455  int nmap = (int)kModalityImageDensityMap.size();
3456  if(nmap != 0) {
3457  short minmax[2];
3458  kModality.getMinMax(minmax);
3459  rval = minmax[1];
3460  for(int i = 0; i < nmap; i++) {
3461  //G4cout << kModalityImageDensityMap[i] << G4endl;
3462  if(_dens <= kModalityImageDensityMap[i]) {
3463  rval = i + minmax[0];
3464  break;
3465  }
3466  }
3467  }
3468  return rval;
3469 }
3470 
3471 
3472 //----- Dose distribution -----//
3473 //
3476  kDose.push_back(doseData);
3477 }
3479  return (int)kDose.size();
3480 }
3481 
3482 // get & set the dose distribution unit
3483 std::string G4GMocrenIO::getDoseDistUnit(int _num) {
3484  // to avoid a warning in the compile process
3485  if(kDoseUnit.size() > static_cast<size_t>(_num)) return kDoseUnit;
3486 
3487  return kDoseUnit;
3488 }
3489 void G4GMocrenIO::setDoseDistUnit(std::string & _unit, int _num) {
3490  // to avoid a warning in the compile process
3491  if(_unit.size() > static_cast<size_t>(_num)) kDoseUnit = _unit;
3492 
3493  //char unit[13];
3494  //std::strncpy(unit, _unit.c_str(), 12);
3495  //doseUnit = unit;
3496  kDoseUnit = _unit;
3497 }
3498 //
3499 void G4GMocrenIO::getDoseDistSize(int _size[3], int _num) {
3500  if(isDoseEmpty())
3501  for(int i = 0; i < 3; i++) _size[i] = 0;
3502  else
3503  kDose[_num].getSize(_size);
3504 }
3505 void G4GMocrenIO::setDoseDistSize(int _size[3], int _num) {
3506 
3507  kDose[_num].setSize(_size);
3508 
3509  //resetDose();
3510 }
3511 
3512 void G4GMocrenIO::setDoseDistMinMax(short _minmax[2], int _num) {
3513 
3514  double minmax[2];
3515  double scale = kDose[_num].getScale();
3516  for(int i = 0; i < 2; i++) minmax[i] = (double)_minmax[i]*scale;
3517  kDose[_num].setMinMax(minmax);
3518 }
3519 void G4GMocrenIO::getDoseDistMinMax(short _minmax[2], int _num) {
3520 
3521  if(isDoseEmpty())
3522  for(int i = 0; i < 2; i++) _minmax[i] = 0;
3523  else {
3524  double minmax[2];
3525  kDose[_num].getMinMax(minmax);
3526  double scale = kDose[_num].getScale();
3527  for(int i = 0; i < 2; i++) _minmax[i] = (short)(minmax[i]/scale+0.5);
3528  }
3529 }
3530 void G4GMocrenIO::setDoseDistMinMax(double _minmax[2], int _num) {
3531 
3532  kDose[_num].setMinMax(_minmax);
3533 }
3534 void G4GMocrenIO::getDoseDistMinMax(double _minmax[2], int _num) {
3535 
3536  if(isDoseEmpty())
3537  for(int i = 0; i < 2; i++) _minmax[i] = 0.;
3538  else
3539  kDose[_num].getMinMax(_minmax);
3540 }
3541 
3542 // get & set the dose distribution image scale
3543 void G4GMocrenIO::setDoseDistScale(double & _scale, int _num) {
3544 
3545  kDose[_num].setScale(_scale);
3546 }
3548 
3549  if(isDoseEmpty())
3550  return 0.;
3551  else
3552  return kDose[_num].getScale();
3553 }
3554 
3555 /*
3556  void G4GMocrenIO::initializeShortDoseDist() {
3557  ;
3558  }
3559  void G4GMocrenIO::finalizeShortDoseDist() {
3560  ;
3561  }
3562 */
3563 // set the dose distribution image
3564 void G4GMocrenIO::setShortDoseDist(short * _image, int _num) {
3565 
3566  int size[3];
3567  kDose[_num].getSize(size);
3568  int dsize = size[0]*size[1];
3569  double * ddata = new double[dsize];
3570  double scale = kDose[_num].getScale();
3571  double minmax[2];
3572  kDose[_num].getMinMax(minmax);
3573  for(int xy = 0; xy < dsize; xy++) {
3574  ddata[xy] = _image[xy]*scale;
3575  if(ddata[xy] < minmax[0]) minmax[0] = ddata[xy];
3576  if(ddata[xy] > minmax[1]) minmax[1] = ddata[xy];
3577  }
3578  kDose[_num].addImage(ddata);
3579 
3580  // set min./max.
3581  kDose[_num].setMinMax(minmax);
3582 }
3583 void G4GMocrenIO::getShortDoseDist(short * _data, int _z, int _num) {
3584 
3585  if(_data == NULL) {
3587  G4cout << "In G4GMocrenIO::getShortDoseDist(), "
3588  << "first argument is NULL pointer. "
3589  << "The argument must be allocated array."
3590  << G4endl;
3591  G4Exception("G4GMocrenIO::getShortDoseDist()",
3592  "gMocren2002", FatalException,
3593  "Error.");
3594  return;
3595  }
3596 
3597  int size[3];
3598  kDose[_num].getSize(size);
3599  //short * shdata = new short[size[0]*size[1]];
3600  double * ddata = kDose[_num].getImage(_z);
3601  double scale = kDose[_num].getScale();
3602  for(int xy = 0; xy < size[0]*size[1]; xy++) {
3603  _data[xy] = (short)(ddata[xy]/scale+0.5); //there is never negative value
3604  }
3605 }
3606 void G4GMocrenIO::getShortDoseDistMinMax(short _minmax[2], int _num) {
3607  double scale = kDose[_num].getScale();
3608  double minmax[2];
3609  kDose[_num].getMinMax(minmax);
3610  for(int i = 0; i < 2; i++)
3611  _minmax[i] = (short)(minmax[i]/scale+0.5);
3612 }
3613 //
3614 void G4GMocrenIO::setDoseDist(double * _image, int _num) {
3615 
3616  kDose[_num].addImage(_image);
3617 }
3618 double * G4GMocrenIO::getDoseDist(int _z, int _num) {
3619 
3620  double * image;
3621  if(isDoseEmpty()) {
3622  image = 0;
3623  } else {
3624  image = kDose[_num].getImage(_z);
3625  }
3626  return image;
3627 }
3628 /*
3629  void G4GMocrenIO::getDoseDist(double * & _image, int _z, int _num) {
3630 
3631  G4cout << " <" << (void*)_image << "> ";
3632  if(isDoseEmpty()) {
3633  _image = 0;
3634  } else {
3635  _image = kDose[_num].getImage(_z);
3636  G4cout << " <" << (void*)_image << "> ";
3637  G4cout << _image[100] << " ";
3638  }
3639  }
3640 */
3641 bool G4GMocrenIO::addDoseDist(std::vector<double *> & _image, int _num) {
3642 
3643  int size[3];
3644  getDoseDistSize(size, _num);
3645  std::vector<double *> dosedist = kDose[_num].getImage();
3646 
3647  int nimg = size[0]*size[1];
3648  for(int z = 0; z < size[2]; z++) {
3649  for(int xy = 0; xy < nimg; xy++) {
3650  dosedist[z][xy] += _image[z][xy];
3651  }
3652  }
3653 
3654  return true;
3655 }
3656 //void setDoseDistDensityMap(float * _map) {doseImageDensityMap = _map;};
3657 // set the dose distribution image displacement
3658 void G4GMocrenIO::setDoseDistCenterPosition(float _center[3], int _num) {
3659 
3660  kDose[_num].setCenterPosition(_center);
3661 }
3662 void G4GMocrenIO::getDoseDistCenterPosition(float _center[3], int _num) {
3663 
3664  if(isDoseEmpty())
3665  for(int i = 0; i < 3; i++) _center[i] = 0;
3666  else
3667  kDose[_num].getCenterPosition(_center);
3668 }
3669 // set & get name of dose distribution
3670 void G4GMocrenIO::setDoseDistName(std::string _name, int _num) {
3671 
3672  kDose[_num].setName(_name);
3673 }
3674 std::string G4GMocrenIO::getDoseDistName(int _num) {
3675 
3676  std::string name;
3677  if(isDoseEmpty())
3678  return name;
3679  else
3680  return kDose[_num].getName();
3681 }
3682 // copy dose distributions
3683 void G4GMocrenIO::copyDoseDist(std::vector<class GMocrenDataPrimitive<double> > & _dose) {
3684  std::vector<class GMocrenDataPrimitive<double> >::iterator itr;
3685  for(itr = kDose.begin(); itr != kDose.end(); itr++) {
3686  _dose.push_back(*itr);
3687  }
3688 }
3689 // merge two dose distributions
3691  if(kDose.size() != _dose.size()) {
3693  G4cout << "G4GMocrenIO::mergeDoseDist() : Error" << G4endl;
3694  G4cout << " Unable to merge the dose distributions,"<< G4endl;
3695  G4cout << " because of different size of dose maps."<< G4endl;
3696  }
3697  return false;
3698  }
3699 
3700  int num = kDose.size();
3701  std::vector<class GMocrenDataPrimitive<double> >::iterator itr1 = kDose.begin();
3702  std::vector<class GMocrenDataPrimitive<double> >::iterator itr2 = _dose.begin();
3703  for(int i = 0; i < num; i++, itr1++, itr2++) {
3705  if(kVerbose > 0)
3706  G4cout << "merged dose distribution [" << i << "]" << G4endl;
3707  *itr1 += *itr2;
3708  }
3709 
3710  return true;
3711 }
3712 //
3714 
3715  if(!isDoseEmpty()) {
3716  for(int i = 0; i < getNumDoseDist(); i++) {
3717  kDose[i].clear();
3718  }
3719  kDose.clear();
3720  }
3721 }
3722 //
3724  if(kDose.empty()) {
3725  //if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
3726  // G4cout << "!!! dose distribution data is empty." << G4endl;
3727  return true;
3728  } else {
3729  return false;
3730  }
3731 }
3732 
3733 //
3735 
3736  double scale;
3737  double minmax[2];
3738 
3739  for(int i = 0; i < (int)kDose.size(); i++) {
3740  kDose[i].getMinMax(minmax);
3741  scale = minmax[1]/DOSERANGE;
3742  kDose[i].setScale(scale);
3743  }
3744 }
3745 
3746 
3747 //----- RoI -----//
3748 
3749 // add one RoI data
3752  kRoi.push_back(roiData);
3753 }
3755  return (int)kRoi.size();
3756 }
3757 
3758 // set/get the ROI image scale
3759 void G4GMocrenIO::setROIScale(double & _scale, int _num) {
3760 
3761  kRoi[_num].setScale(_scale);
3762 }
3763 double G4GMocrenIO::getROIScale(int _num) {
3764 
3765  if(isROIEmpty())
3766  return 0.;
3767  else
3768  return kRoi[_num].getScale();
3769 }
3770 // set the ROI image
3771 void G4GMocrenIO::setROI(short * _image, int _num) {
3772 
3773  kRoi[_num].addImage(_image);
3774 }
3775 short * G4GMocrenIO::getROI(int _z, int _num) {
3776 
3777  if(isROIEmpty())
3778  return 0;
3779  else
3780  return kRoi[_num].getImage(_z);
3781 }
3782 // set/get the ROI image size
3783 void G4GMocrenIO::setROISize(int _size[3], int _num) {
3784 
3785  return kRoi[_num].setSize(_size);
3786 }
3787 void G4GMocrenIO::getROISize(int _size[3], int _num) {
3788 
3789  if(isROIEmpty())
3790  for(int i = 0; i < 3; i++) _size[i] = 0;
3791  else
3792  return kRoi[_num].getSize(_size);
3793 }
3794 // set/get the ROI image min. and max.
3795 void G4GMocrenIO::setROIMinMax(short _minmax[2], int _num) {
3796 
3797  kRoi[_num].setMinMax(_minmax);
3798 }
3799 void G4GMocrenIO::getROIMinMax(short _minmax[2], int _num) {
3800 
3801  if(isROIEmpty())
3802  for(int i = 0; i < 2; i++) _minmax[i] = 0;
3803  else
3804  kRoi[_num].getMinMax(_minmax);
3805 }
3806 // set/get the ROI image displacement
3807 void G4GMocrenIO::setROICenterPosition(float _center[3], int _num) {
3808 
3809  kRoi[_num].setCenterPosition(_center);
3810 }
3811 void G4GMocrenIO::getROICenterPosition(float _center[3], int _num) {
3812 
3813  if(isROIEmpty())
3814  for(int i = 0; i < 3; i++) _center[i] = 0;
3815  else
3816  kRoi[_num].getCenterPosition(_center);
3817 }
3818 //
3820 
3821  if(!isROIEmpty()) {
3822  for(int i = 0; i < getNumROI(); i++) {
3823  kRoi[i].clear();
3824  }
3825  kRoi.clear();
3826  }
3827 }
3828 //
3830  if(kRoi.empty()) {
3831  //if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
3832  // G4cout << "!!! ROI data is empty." << G4endl;
3833  return true;
3834  } else {
3835  return false;
3836  }
3837 }
3838 
3839 
3840 
3841 //----- Track information -----//
3842 
3844  return (int)kSteps.size();
3845 }
3847  return (int)kTracks.size();
3848 }
3849 void G4GMocrenIO::addTrack(float * _tracks) {
3850  kSteps.push_back(_tracks);
3851 }
3852 void G4GMocrenIO::setTracks(std::vector<float *> & _tracks) {
3853  kSteps = _tracks;
3854 }
3855 std::vector<float *> & G4GMocrenIO::getTracks() {
3856  return kSteps;
3857 }
3858 void G4GMocrenIO::addTrackColor(unsigned char * _colors) {
3859  kStepColors.push_back(_colors);
3860 }
3861 void G4GMocrenIO::setTrackColors(std::vector<unsigned char *> & _trackColors) {
3862  kStepColors = _trackColors;
3863 }
3864 std::vector<unsigned char *> & G4GMocrenIO::getTrackColors() {
3865  return kStepColors;
3866 }
3867 void G4GMocrenIO::copyTracks(std::vector<float *> & _tracks,
3868  std::vector<unsigned char *> & _colors) {
3869  std::vector<float *>::iterator titr;
3870  for(titr = kSteps.begin(); titr != kSteps.end(); titr++) {
3871  float * pts = new float[6];
3872  for(int i = 0; i < 6; i++) {
3873  pts[i] = (*titr)[i];
3874  }
3875  _tracks.push_back(pts);
3876  }
3877 
3878  std::vector<unsigned char *>::iterator citr;
3879  for(citr = kStepColors.begin(); citr != kStepColors.end(); citr++) {
3880  unsigned char * pts = new unsigned char[3];
3881  for(int i = 0; i < 3; i++) {
3882  pts[i] = (*citr)[i];
3883  }
3884  _colors.push_back(pts);
3885  }
3886 }
3887 void G4GMocrenIO::mergeTracks(std::vector<float *> & _tracks,
3888  std::vector<unsigned char *> & _colors) {
3889  std::vector<float *>::iterator titr;
3890  for(titr = _tracks.begin(); titr != _tracks.end(); titr++) {
3891  addTrack(*titr);
3892  }
3893 
3894  std::vector<unsigned char *>::iterator citr;
3895  for(citr = _colors.begin(); citr != _colors.end(); citr++) {
3896  addTrackColor(*citr);
3897  }
3898 }
3899 void G4GMocrenIO::addTrack(std::vector<float *> & _steps, unsigned char _color[3]) {
3900 
3901  std::vector<float *>::iterator itr = _steps.begin();
3902  std::vector<struct GMocrenTrack::Step> steps;
3903  for(; itr != _steps.end(); itr++) {
3904  struct GMocrenTrack::Step step;
3905  for(int i = 0; i < 3; i++) {
3906  step.startPoint[i] = (*itr)[i];
3907  step.endPoint[i] = (*itr)[i+3];
3908  }
3909  steps.push_back(step);
3910  }
3911  GMocrenTrack track;
3912  track.setTrack(steps);
3913  track.setColor(_color);
3914  kTracks.push_back(track);
3915 
3916 }
3917 void G4GMocrenIO::getTrack(int _num, std::vector<float *> & _steps,
3918  std::vector<unsigned char *> & _color) {
3919 
3920  if(_num > (int)kTracks.size()) {
3922  G4cout << "ERROR in getTrack() : " << G4endl;
3923  G4Exception("G4GMocrenIO::getTrack()",
3924  "gMocren2003", FatalException,
3925  "Error.");
3926  }
3927  unsigned char * color = new unsigned char[3];
3928  kTracks[_num].getColor(color);
3929  _color.push_back(color);
3930 
3931  // steps
3932  int nsteps = kTracks[_num].getNumberOfSteps();
3933  for(int isteps = 0; isteps < nsteps; isteps++) {
3934  float * stepPoints = new float[6];
3935  kTracks[_num].getStep(stepPoints[0], stepPoints[1], stepPoints[2],
3936  stepPoints[3], stepPoints[4], stepPoints[5],
3937  isteps);
3938  _steps.push_back(stepPoints);
3939  }
3940 }
3941 
3942 void G4GMocrenIO::translateTracks(std::vector<float> & _translate) {
3943  std::vector<class GMocrenTrack>::iterator itr = kTracks.begin();
3944  for(; itr != kTracks.end(); itr++) {
3945  itr->translate(_translate);
3946  }
3947 }
3948 
3949 
3950 
3951 
3952 //----- Detector information -----//
3954  return (int)kDetectors.size();
3955 }
3956 void G4GMocrenIO::addDetector(std::string & _name,
3957  std::vector<float *> & _det,
3958  unsigned char _color[3]) {
3959 
3960  std::vector<float *>::iterator itr = _det.begin();
3961  std::vector<struct GMocrenDetector::Edge> edges;
3962  for(; itr != _det.end(); itr++) {
3963  struct GMocrenDetector::Edge edge;
3964  for(int i = 0; i < 3; i++) {
3965  edge.startPoint[i] = (*itr)[i];
3966  edge.endPoint[i] = (*itr)[i+3];
3967  }
3968  edges.push_back(edge);
3969  }
3970  GMocrenDetector detector;
3971  detector.setDetector(edges);
3972  detector.setColor(_color);
3973  detector.setName(_name);
3974  kDetectors.push_back(detector);
3975 
3976 }
3977 
3978 void G4GMocrenIO::getDetector(int _num, std::vector<float *> & _edges,
3979  std::vector<unsigned char *> & _color,
3980  std::string & _detName) {
3981 
3982  if(_num > (int)kDetectors.size()) {
3984  G4cout << "ERROR in getDetector() : " << G4endl;
3985 
3986  G4Exception("G4GMocrenIO::getDetector()",
3987  "gMocren2004", FatalException,
3988  "Error.");
3989  }
3990 
3991  _detName = kDetectors[_num].getName();
3992 
3993  unsigned char * color = new unsigned char[3];
3994  kDetectors[_num].getColor(color);
3995  _color.push_back(color);
3996 
3997  // edges
3998  int nedges = kDetectors[_num].getNumberOfEdges();
3999  for(int ne = 0; ne < nedges; ne++) {
4000  float * edgePoints = new float[6];
4001  kDetectors[_num].getEdge(edgePoints[0], edgePoints[1], edgePoints[2],
4002  edgePoints[3], edgePoints[4], edgePoints[5],
4003  ne);
4004  _edges.push_back(edgePoints);
4005  }
4006 }
4007 
4008 void G4GMocrenIO::translateDetector(std::vector<float> & _translate) {
4009  std::vector<class GMocrenDetector>::iterator itr = kDetectors.begin();
4010  for(; itr != kDetectors.end(); itr++) {
4011  itr->translate(_translate);
4012  }
4013 }
4014 
4015 // endian conversion
4016 template <typename T>
4017 void G4GMocrenIO::convertEndian(char * _val, T & _rval) {
4018 
4019  if((kLittleEndianOutput && !kLittleEndianInput) || // big endian
4020  (!kLittleEndianOutput && kLittleEndianInput)) { // little endian
4021 
4022  const int SIZE = sizeof(_rval);
4023  char ctemp;
4024  for(int i = 0; i < SIZE/2; i++) {
4025  ctemp = _val[i];
4026  _val[i] = _val[SIZE - 1 - i];
4027  _val[SIZE - 1 - i] = ctemp;
4028  }
4029  }
4030  _rval = *(T *)_val;
4031 }
4032 
4033 // inversion of byte order
4034 template <typename T>
4035 void G4GMocrenIO::invertByteOrder(char * _val, T & _rval) {
4036 
4037  const int SIZE = sizeof(_rval);
4038  //char * cval = new char[SIZE];
4039  union {
4040  char cu[16];
4041  T tu;
4042  } uni;
4043  for(int i = 0; i < SIZE; i++) {
4044  uni.cu[i] = _val[SIZE-1-i];
4045  //cval[i] = _val[SIZE-i-1];
4046  }
4047  //_rval = *(T *)cval;
4048  _rval = uni.tu;
4049  //delete [] cval;
4050 }
4051 
4052 //----- kVerbose information -----//
4054  kVerbose = _level;
4055 }
4056 
void getEdge(float &_startx, float &_starty, float &_startz, float &_endx, float &_endy, float &_endz, int _num)
Definition: G4GMocrenIO.cc:317
void getModalityImageSize(int _size[3])
std::vector< T * > kImage
Definition: G4GMocrenIO.hh:56
void convertEndian(char *, Type &)
unsigned int getPointerToTrackData()
void initialize()
Definition: G4GMocrenIO.cc:413
bool storeData4()
Definition: G4GMocrenIO.cc:462
static std::vector< float > kModalityImageDensityMap
Definition: G4GMocrenIO.hh:203
void getModalityCenterPosition(float _center[3])
void setVerboseLevel(int _level)
void newDoseDist()
short convertDensityToHU(float &_dens)
void setDoseDistUnit(std::string &_unit, int _num=0)
int getNumTracks()
double getModalityImageScale()
unsigned int getPointerToROIData()
bool isROIEmpty()
void addImage(T *_image)
Definition: G4GMocrenIO.cc:217
void setVersion(std::string &_version)
double getDoseDistScale(int _num=0)
bool storeData3()
void clearModalityImage()
static char kLittleEndianOutput
Definition: G4GMocrenIO.hh:179
double getROIScale(int _num=0)
void setROI(short *_image, int _num=0)
void addTrack(float *_tracks)
std::vector< float > & getModalityImageDensityMap()
bool kTracksWillBeStored
Definition: G4GMocrenIO.hh:219
G4double z
Definition: TRTMaterials.hh:39
const int DOSERANGE
Definition: G4GMocrenIO.cc:56
std::vector< float * > & getTracks()
unsigned int getPointerToModalityData()
void invertByteOrder(char *_val, T &_rval)
G4String name
Definition: TRTMaterials.hh:40
short * getROI(int _z, int _num=0)
void setDoseDistMinMax(short _minmax[2], int _num=0)
void setModalityImageUnit(std::string &_unit)
static std::vector< class GMocrenDataPrimitive< short > > kRoi
Definition: G4GMocrenIO.hh:212
void setTrack(std::vector< struct Step > &_aTrack)
Definition: G4GMocrenIO.hh:112
void setModalityImageSize(int _size[3])
void getShortDoseDist(short *_data, int _z, int _num=0)
bool storeData2()
bool addDoseDist(std::vector< double * > &_image, int _num=0)
void setComment(std::string &_comment)
Definition: G4GMocrenIO.hh:275
short getModalityImageMax()
std::vector< struct Step > kTrack
Definition: G4GMocrenIO.hh:98
std::vector< unsigned char * > & getTrackColors()
GMocrenDataPrimitive< T > & operator=(const GMocrenDataPrimitive< T > &_right)
Definition: G4GMocrenIO.cc:74
void getROICenterPosition(float _center[3], int _num=0)
GMocrenDataPrimitive< T > & operator+(const GMocrenDataPrimitive< T > &_right)
Definition: G4GMocrenIO.cc:93
static std::vector< class GMocrenDataPrimitive< double > > kDose
Definition: G4GMocrenIO.hh:207
void translate(std::vector< float > &_tranlate)
Definition: G4GMocrenIO.cc:334
void copyTracks(std::vector< float * > &_tracks, std::vector< unsigned char * > &_colors)
std::string getDoseDistName(int _num=0)
unsigned char kColor[3]
Definition: G4GMocrenIO.hh:135
void setROIMinMax(short _minmax[2], int _num=0)
static std::vector< unsigned char * > kStepColors
Definition: G4GMocrenIO.hh:216
std::string getName()
Definition: G4GMocrenIO.cc:242
void translateDetector(std::vector< float > &_translate)
void getDoseDistSize(int _size[3], int _num=0)
static std::vector< class GMocrenDetector > kDetectors
Definition: G4GMocrenIO.hh:222
static float kVoxelSpacing[3]
Definition: G4GMocrenIO.hh:198
void setDoseDist(double *_image, int _num=0)
double * getDoseDist(int _z, int _num=0)
void clearROIAll()
bool retrieveData()
void setColor(unsigned char _color[3])
Definition: G4GMocrenIO.hh:150
void setPointerToROIData(unsigned int &_pointer)
void addTrackColor(unsigned char *_colors)
void addOneEvent()
static int kVerbose
Definition: G4GMocrenIO.hh:225
void calcPointers3()
G4GLOB_DLL std::ostream G4cout
static char kLittleEndianInput
Definition: G4GMocrenIO.hh:178
static std::string kVersion
Definition: G4GMocrenIO.hh:172
std::vector< T * > & getImage()
Definition: G4GMocrenIO.cc:221
void getShortDoseDistMinMax(short _minmax[2], int _num=0)
void setScale(double &_scale)
Definition: G4GMocrenIO.cc:196
void mergeTracks(std::vector< float * > &_tracks, std::vector< unsigned char * > &_colors)
void getCenterPosition(float _center[3])
Definition: G4GMocrenIO.cc:234
void setMinMax(T _minmax[2])
Definition: G4GMocrenIO.cc:204
void setROICenterPosition(float _center[3], int _num=0)
static int kNumberOfEvents
Definition: G4GMocrenIO.hh:184
void clearDoseDistAll()
GMocrenDataPrimitive< T > & operator+=(const GMocrenDataPrimitive< T > &_right)
Definition: G4GMocrenIO.cc:135
void addEdge(float _startx, float _starty, float _startz, float _endx, float _endy, float _endz)
Definition: G4GMocrenIO.cc:306
static const double cm3
Definition: G4SIunits.hh:120
void setImage(std::vector< T * > &_image)
Definition: G4GMocrenIO.cc:213
int getNumTracks4()
void getROISize(int _size[3], int _num=0)
void setModalityImageMinMax(short _minmax[2])
std::string getDoseDistUnit(int _num=0)
std::string & getVersion()
void calcDoseDistScale()
void setModalityImage(short *_image)
void setCenterPosition(float _center[3])
Definition: G4GMocrenIO.cc:230
void setVoxelSpacing(float _spacing[3])
void setSize(int _size[3])
Definition: G4GMocrenIO.cc:188
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
bool retrieveData2()
void copyDoseDist(std::vector< class GMocrenDataPrimitive< double > > &_dose)
static unsigned int kPointerToModalityData
Definition: G4GMocrenIO.hh:187
int getNumDoseDist()
void getTrack(int _num, std::vector< float * > &_steps, std::vector< unsigned char * > &_color)
static std::vector< class GMocrenTrack > kTracks
Definition: G4GMocrenIO.hh:218
void setTrackColors(std::vector< unsigned char * > &_trackColors)
void getMinMax(T _minmax[2])
Definition: G4GMocrenIO.cc:208
short * getModalityImage(int _z)
void getModalityImageMinMax(short _minmax[2])
void translate(std::vector< float > &_tranlate)
Definition: G4GMocrenIO.cc:283
bool retrieveData4()
std::vector< struct Edge > kDetector
Definition: G4GMocrenIO.hh:134
static const double g
Definition: G4SIunits.hh:180
void getStep(float &_startx, float &_starty, float &_startz, float &_endx, float &_endy, float &_endz, int _num)
Definition: G4GMocrenIO.cc:266
void setModalityImageScale(double &_scale)
static std::string kFileName
Definition: G4GMocrenIO.hh:175
bool retrieveData3()
const int VERLENGTH
Definition: G4GMocrenIO.cc:399
static std::string kDoseUnit
Definition: G4GMocrenIO.hh:209
std::ofstream ofile
Definition: clparse.cc:45
bool mergeDoseDist(std::vector< class GMocrenDataPrimitive< double > > &_dose)
void getVoxelSpacing(float _spacing[3])
void setDoseDistCenterPosition(float _center[3], int _num=0)
void translateTracks(std::vector< float > &_translateo)
unsigned int getPointerToDoseDistData(int _elem=0)
void setNumberOfEvents(int &_numberOfEvents)
static Verbosity GetVerbosity()
void setPointerToModalityData(unsigned int &_pointer)
bool storeData()
Definition: G4GMocrenIO.cc:454
#define G4endl
Definition: G4ios.hh:61
void calcPointers4()
void calcPointers2()
void setName(std::string &_name)
Definition: G4GMocrenIO.hh:156
static const double keV
Definition: G4SIunits.hh:213
static std::string kId
Definition: G4GMocrenIO.hh:169
void getROIMinMax(short _minmax[2], int _num=0)
int & getNumberOfEvents()
void setShortDoseDist(short *_image, int _num=0)
void setDoseDistScale(double &_scale, int _num=0)
void addDetector(std::string &_name, std::vector< float * > &_det, unsigned char _color[3])
void setLittleEndianInput(bool _little)
static class GMocrenDataPrimitive< short > kModality
Definition: G4GMocrenIO.hh:201
void setColor(unsigned char _color[3])
Definition: G4GMocrenIO.hh:113
void setDoseDistSize(int _size[3], int _num=0)
static std::string kComment
Definition: G4GMocrenIO.hh:181
static std::vector< unsigned int > kPointerToDoseDistData
Definition: G4GMocrenIO.hh:189
static unsigned int kPointerToDetectorData
Definition: G4GMocrenIO.hh:195
void setModalityCenterPosition(float _center[3])
void getDetector(int _num, std::vector< float * > &_edges, std::vector< unsigned char * > &_color, std::string &_detectorName)
void setROISize(int _size[3], int _num=0)
void setDetector(std::vector< struct Edge > &_aDetector)
Definition: G4GMocrenIO.hh:149
static std::vector< float * > kSteps
Definition: G4GMocrenIO.hh:215
unsigned char kColor[3]
Definition: G4GMocrenIO.hh:99
void setROIScale(double &_scale, int _num=0)
void getDoseDistCenterPosition(float _center[3], int _num=0)
void getSize(int _size[3])
Definition: G4GMocrenIO.cc:192
void setName(std::string &_name)
Definition: G4GMocrenIO.cc:238
bool isDoseEmpty()
void setPointerToTrackData(unsigned int &_pointer)
void addPointerToDoseDistData(unsigned int &_pointer)
void setTracks(std::vector< float * > &_tracks)
void getDoseDistMinMax(short _minmax[2], int _num=0)
void setLittleEndianOutput(bool _little)
void setDoseDistName(std::string _name, int _num=0)
short getModalityImageMin()
int getNumberOfDetectors()
std::string getModalityImageUnit()
static std::string kModalityUnit
Definition: G4GMocrenIO.hh:204
void setModalityImageDensityMap(std::vector< float > &_map)
void addStep(float _startx, float _starty, float _startz, float _endx, float _endy, float _endz)
Definition: G4GMocrenIO.cc:255
const int IDLENGTH
Definition: G4GMocrenIO.cc:398
static unsigned int kPointerToROIData
Definition: G4GMocrenIO.hh:191
static unsigned int kPointerToTrackData
Definition: G4GMocrenIO.hh:193