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