00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040 #include "G4GMocrenIO.hh"
00041 #include <iostream>
00042 #include <ctime>
00043 #include <sstream>
00044 #include <iomanip>
00045 #include <cstdlib>
00046 #include <cstring>
00047
00048 #include "globals.hh"
00049 #include "G4VisManager.hh"
00050
00051 #if defined(_WIN32)
00052 #define LITTLE_ENDIAN 1234
00053 #define BYTE_ORDER LITTLE_ENDIAN
00054 #endif
00055
00056 const int DOSERANGE = 25000;
00057
00058
00059 template <typename T>
00060 GMocrenDataPrimitive<T>::GMocrenDataPrimitive () {
00061 clear();
00062 }
00063 template <typename T>
00064 GMocrenDataPrimitive<T>::~GMocrenDataPrimitive () {
00065
00066
00067
00068
00069
00070
00071 }
00072
00073 template <typename T> GMocrenDataPrimitive<T> &
00074 GMocrenDataPrimitive<T>::operator = (const GMocrenDataPrimitive<T> & _right) {
00075 for(int i = 0; i < 3; i++) {
00076 kSize[i] = _right.kSize[i];
00077 kCenter[i] = _right.kCenter[i];
00078 }
00079 kScale = _right.kScale;
00080 for(int i = 0; i < 2; i++) kMinmax[i] = _right.kMinmax[i];
00081 int num = kSize[0]*kSize[1];
00082 kImage.clear();
00083 for(int z = 0; z < kSize[2]; z++) {
00084 T * img = new T[num];
00085 for(int i = 0; i < num; i++) img[i] =_right.kImage[z][i];
00086 kImage.push_back(img);
00087 }
00088 return *this;
00089 }
00090
00091 template <typename T> GMocrenDataPrimitive<T> &
00092 GMocrenDataPrimitive<T>::operator + (const GMocrenDataPrimitive<T> & _right) {
00093
00094 GMocrenDataPrimitive<T> rprim;
00095 bool stat = true;
00096 for(int i = 0; i < 3; i++) {
00097 if(kSize[i] != _right.kSize[i]) stat = false;
00098 if(kCenter[i] != _right.kCenter[i]) stat = false;
00099 }
00100 if(!stat) {
00101 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
00102 G4cout << "Warning: operator + "
00103 << " Cannot do the operator +"
00104 << G4endl;
00105 return *this;
00106 }
00107
00108 rprim.setSize(kSize);
00109 rprim.setCenterPosition(kCenter);
00110
00111 T mms[2] = {9e100,-9e100};
00112
00113
00114
00115 int num = kSize[0]*kSize[1];
00116 for(int z = 0; z < kSize[2]; z++) {
00117 T * img = new T[num];
00118 for(int xy = 0; xy < num; xy++) {
00119 img[xy] = kImage[z][xy] + _right.kImage[z][xy];
00120 if(mms[0] > img[xy]) mms[0] = img[xy];
00121 if(mms[1] < img[xy]) mms[1] = img[xy];
00122 }
00123 rprim.addImage(img);
00124 }
00125 rprim.setMinMax(mms);
00126
00127 T scl = mms[1]/DOSERANGE;
00128 rprim.setScale(scl);
00129
00130 return rprim;
00131 }
00132
00133 template <typename T> GMocrenDataPrimitive<T> &
00134 GMocrenDataPrimitive<T>::operator += (const GMocrenDataPrimitive<T> & _right) {
00135
00136 bool stat = true;
00137 for(int i = 0; i < 3; i++) {
00138 if(kSize[i] != _right.kSize[i]) stat = false;
00139 if(kCenter[i] != _right.kCenter[i]) stat = false;
00140 }
00141 if(!stat) {
00142 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
00143 G4cout << "Warning: operator += " << G4endl
00144 << " Cannot do the operator +="
00145 << G4endl;
00146 return *this;
00147 }
00148
00149 if(kMinmax[0] > _right.kMinmax[0]) kMinmax[0] = _right.kMinmax[0];
00150 if(kMinmax[1] < _right.kMinmax[1]) kMinmax[1] = _right.kMinmax[1];
00151
00152 int num = kSize[0]*kSize[1];
00153 for(int z = 0; z < kSize[2]; z++) {
00154 for(int xy = 0; xy < num; xy++) {
00155 kImage[z][xy] += _right.kImage[z][xy];
00156 if(kMinmax[0] > kImage[z][xy]) kMinmax[0] = kImage[z][xy];
00157 if(kMinmax[1] < kImage[z][xy]) kMinmax[1] = kImage[z][xy];
00158 }
00159 }
00160
00161 kScale = kMinmax[1]/DOSERANGE;
00162
00163 return *this;
00164 }
00165
00166 template <typename T>
00167 void GMocrenDataPrimitive<T>::clear() {
00168 for(int i = 0; i < 3; i++) {
00169 kSize[i] = 0;
00170 kCenter[i] = 0.;
00171 }
00172 kScale = 1.;
00173 kMinmax[0] = (T)32109;
00174 kMinmax[1] = (T)-32109;
00175
00176 clearImage();
00177 }
00178 template <typename T>
00179 void GMocrenDataPrimitive<T>::clearImage() {
00180 typename std::vector<T *>::iterator itr;
00181 for(itr = kImage.begin(); itr != kImage.end(); itr++) {
00182 delete [] *itr;
00183 }
00184 kImage.clear();
00185 }
00186 template <typename T>
00187 void GMocrenDataPrimitive<T>::setSize(int _size[3]) {
00188 for(int i = 0; i < 3; i++) kSize[i] = _size[i];
00189 }
00190 template <typename T>
00191 void GMocrenDataPrimitive<T>::getSize(int _size[3]) {
00192 for(int i = 0; i < 3; i++) _size[i] = kSize[i];
00193 }
00194 template <typename T>
00195 void GMocrenDataPrimitive<T>::setScale(double & _scale) {
00196 kScale = _scale;
00197 }
00198 template <typename T>
00199 double GMocrenDataPrimitive<T>::getScale() {
00200 return kScale;
00201 }
00202 template <typename T>
00203 void GMocrenDataPrimitive<T>::setMinMax(T _minmax[2]) {
00204 for(int i = 0; i < 2; i++) kMinmax[i] = _minmax[i];
00205 }
00206 template <typename T>
00207 void GMocrenDataPrimitive<T>::getMinMax(T _minmax[2]) {
00208 for(int i = 0; i < 2; i++) _minmax[i] = kMinmax[i];
00209
00210 }
00211 template <typename T>
00212 void GMocrenDataPrimitive<T>::setImage(std::vector<T *> & _image) {
00213 kImage = _image;
00214 }
00215 template <typename T>
00216 void GMocrenDataPrimitive<T>::addImage(T * _image) {
00217 kImage.push_back(_image);
00218 }
00219 template <typename T>
00220 std::vector<T *> & GMocrenDataPrimitive<T>::getImage() {
00221 return kImage;
00222 }
00223 template <typename T>
00224 T * GMocrenDataPrimitive<T>::getImage(int _z) {
00225 if(_z >= (int)kImage.size()) return 0;
00226 return kImage[_z];
00227 }
00228 template <typename T>
00229 void GMocrenDataPrimitive<T>::setCenterPosition(float _center[3]) {
00230 for(int i = 0; i < 3; i++) kCenter[i] = _center[i];
00231 }
00232 template <typename T>
00233 void GMocrenDataPrimitive<T>::getCenterPosition(float _center[3]) {
00234 for(int i = 0; i < 3; i++) _center[i] = kCenter[i];
00235 }
00236 template <typename T>
00237 void GMocrenDataPrimitive<T>::setName(std::string & _name) {
00238 kDataName = _name;
00239 }
00240 template <typename T>
00241 std::string GMocrenDataPrimitive<T>::getName() {
00242 return kDataName;
00243 }
00244
00245
00246
00247
00248
00249 GMocrenTrack::GMocrenTrack() {
00250 kTrack.clear();
00251 for(int i = 0; i < 3; i++) kColor[i] = 0;
00252 }
00253
00254 void GMocrenTrack::addStep(float _startx, float _starty, float _startz,
00255 float _endx, float _endy, float _endz) {
00256 struct Step step;
00257 step.startPoint[0] = _startx;
00258 step.startPoint[1] = _starty;
00259 step.startPoint[2] = _startz;
00260 step.endPoint[0] = _endx;
00261 step.endPoint[1] = _endy;
00262 step.endPoint[2] = _endz;
00263 kTrack.push_back(step);
00264 }
00265 void GMocrenTrack::getStep(float & _startx, float & _starty, float & _startz,
00266 float & _endx, float & _endy, float & _endz,
00267 int _num) {
00268 if(_num >= (int)kTrack.size()) {
00269 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
00270 G4cout << "GMocrenTrack::getStep(...) Error: "
00271 << "invalid step # : " << _num << G4endl;
00272 return;
00273 }
00274
00275 _startx = kTrack[_num].startPoint[0];
00276 _starty = kTrack[_num].startPoint[1];
00277 _startz = kTrack[_num].startPoint[2];
00278 _endx = kTrack[_num].endPoint[0];
00279 _endy = kTrack[_num].endPoint[1];
00280 _endz = kTrack[_num].endPoint[2];
00281 }
00282 void GMocrenTrack::translate(std::vector<float> & _translate) {
00283 std::vector<struct Step>::iterator itr = kTrack.begin();
00284 for(; itr != kTrack.end(); itr++) {
00285 for(int i = 0; i < 3; i++ ) {
00286 itr->startPoint[i] += _translate[i];
00287 itr->endPoint[i] += _translate[i];
00288 }
00289 }
00290 }
00291
00292
00293
00294
00295
00296
00297
00298
00299
00300 GMocrenDetector::GMocrenDetector() {
00301 kDetector.clear();
00302 for(int i = 0; i < 3; i++) kColor[i] = 0;
00303 }
00304
00305 void GMocrenDetector::addEdge(float _startx, float _starty, float _startz,
00306 float _endx, float _endy, float _endz) {
00307 struct Edge edge;
00308 edge.startPoint[0] = _startx;
00309 edge.startPoint[1] = _starty;
00310 edge.startPoint[2] = _startz;
00311 edge.endPoint[0] = _endx;
00312 edge.endPoint[1] = _endy;
00313 edge.endPoint[2] = _endz;
00314 kDetector.push_back(edge);
00315 }
00316 void GMocrenDetector::getEdge(float & _startx, float & _starty, float & _startz,
00317 float & _endx, float & _endy, float & _endz,
00318 int _num) {
00319 if(_num >= (int)kDetector.size()) {
00320 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
00321 G4cout << "GMocrenDetector::getEdge(...) Error: "
00322 << "invalid edge # : " << _num << G4endl;
00323 return;
00324 }
00325
00326 _startx = kDetector[_num].startPoint[0];
00327 _starty = kDetector[_num].startPoint[1];
00328 _startz = kDetector[_num].startPoint[2];
00329 _endx = kDetector[_num].endPoint[0];
00330 _endy = kDetector[_num].endPoint[1];
00331 _endz = kDetector[_num].endPoint[2];
00332 }
00333 void GMocrenDetector::translate(std::vector<float> & _translate) {
00334 std::vector<struct Edge>::iterator itr = kDetector.begin();
00335 for(; itr != kDetector.end(); itr++) {
00336 for(int i = 0; i < 3; i++) {
00337 itr->startPoint[i] += _translate[i];
00338 itr->endPoint[i] += _translate[i];
00339 }
00340 }
00341 }
00342
00343
00344
00345
00346
00347
00348
00349
00350
00351
00352 std::string G4GMocrenIO::kId;
00353 std::string G4GMocrenIO::kVersion = "2.0.0";
00354 int G4GMocrenIO::kNumberOfEvents = 0;
00355 char G4GMocrenIO::kLittleEndianInput = true;
00356
00357 #if BYTE_ORDER == LITTLE_ENDIAN
00358 char G4GMocrenIO::kLittleEndianOutput = true;
00359 #else
00360 char G4GMocrenIO::kLittleEndianOutput = false;
00361 #endif
00362 std::string G4GMocrenIO::kComment;
00363
00364 std::string G4GMocrenIO::kFileName = "dose.gdd";
00365
00366
00367 unsigned int G4GMocrenIO::kPointerToModalityData = 0;
00368 std::vector<unsigned int> G4GMocrenIO::kPointerToDoseDistData;
00369 unsigned int G4GMocrenIO::kPointerToROIData = 0;
00370 unsigned int G4GMocrenIO::kPointerToTrackData = 0;
00371 unsigned int G4GMocrenIO::kPointerToDetectorData = 0;
00372
00373
00374 float G4GMocrenIO::kVoxelSpacing[3] = {0., 0., 0.};
00375 class GMocrenDataPrimitive<short> G4GMocrenIO::kModality;
00376 std::vector<float> G4GMocrenIO::kModalityImageDensityMap;
00377 std::string G4GMocrenIO::kModalityUnit = "g/cm3 ";
00378
00379
00380 std::vector<class GMocrenDataPrimitive<double> > G4GMocrenIO::kDose;
00381 std::string G4GMocrenIO::kDoseUnit = "keV ";
00382
00383
00384 std::vector<class GMocrenDataPrimitive<short> > G4GMocrenIO::kRoi;
00385
00386
00387 std::vector<float *> G4GMocrenIO::kSteps;
00388 std::vector<unsigned char *> G4GMocrenIO::kStepColors;
00389 std::vector<class GMocrenTrack> G4GMocrenIO::kTracks;
00390
00391
00392 std::vector<class GMocrenDetector> G4GMocrenIO::kDetectors;
00393
00394
00395 int G4GMocrenIO::kVerbose = 0;
00396
00397 const int IDLENGTH = 21;
00398 const int VERLENGTH = 6;
00399
00400
00401 G4GMocrenIO::G4GMocrenIO()
00402 : kTracksWillBeStored(true) {
00403 ;
00404 }
00405
00406
00407 G4GMocrenIO::~G4GMocrenIO() {
00408 ;
00409 }
00410
00411
00412 void G4GMocrenIO::initialize() {
00413
00414 kId.clear();
00415 kVersion = "2.0.0";
00416 kNumberOfEvents = 0;
00417 kLittleEndianInput = true;
00418 #if BYTE_ORDER == LITTLE_ENDIAN
00419 kLittleEndianOutput = true;
00420 #else // Big endian
00421 kLittleEndianOutput = false;
00422 #endif
00423 kComment.clear();
00424 kFileName = "dose.gdd";
00425 kPointerToModalityData = 0;
00426 kPointerToDoseDistData.clear();
00427 kPointerToROIData = 0;
00428 kPointerToTrackData = 0;
00429
00430 for(int i = 0; i < 3; i++) kVoxelSpacing[i] = 0.;
00431 kModality.clear();
00432 kModalityImageDensityMap.clear();
00433 kModalityUnit = "g/cm3 ";
00434
00435 kDose.clear();
00436 kDoseUnit = "keV ";
00437
00438 kRoi.clear();
00439
00440 std::vector<float *>::iterator itr;
00441 for(itr = kSteps.begin(); itr != kSteps.end(); itr++) delete [] *itr;
00442 kSteps.clear();
00443 std::vector<unsigned char *>::iterator citr;
00444 for(citr = kStepColors.begin(); citr != kStepColors.end(); citr++)
00445 delete [] *citr;
00446 kStepColors.clear();
00447 kTracksWillBeStored = true;
00448
00449
00450 kVerbose = 0;
00451 }
00452
00453 bool G4GMocrenIO::storeData() {
00454 return storeData4();
00455 }
00456
00457 bool G4GMocrenIO::storeData(char * _filename) {
00458 return storeData4(_filename);
00459 }
00460
00461 bool G4GMocrenIO::storeData4() {
00462
00463 bool DEBUG = false;
00464
00465 if(DEBUG || kVerbose > 0)
00466 G4cout << ">>>>>>> store data (ver.4) <<<<<<<" << G4endl;
00467 if(DEBUG || kVerbose > 0)
00468 G4cout << " " << kFileName << G4endl;
00469
00470
00471 std::ofstream ofile(kFileName.c_str(),
00472 std::ios_base::out|std::ios_base::binary);
00473 if(DEBUG || kVerbose > 0)
00474 G4cout << " file open status: " << ofile << G4endl;
00475
00476
00477 ofile.write("gMocren ", 8);
00478
00479
00480 unsigned char ver = 0x04;
00481 ofile.write((char *)&ver, 1);
00482
00483
00484
00485 char littleEndian = 0x01;
00486 ofile.write((char *)&littleEndian, sizeof(char));
00487 if(DEBUG || kVerbose > 0) {
00488
00489 G4cout << "Endian: " << (int)littleEndian << G4endl;
00490 }
00491
00492
00493 float ftmp[6];
00494 int itmp[6];
00495 short stmp[6];
00496
00497
00498 int commentLength = 1024;
00499 if(kLittleEndianOutput) {
00500 ofile.write((char *)&commentLength, 4);
00501 } else {
00502 invertByteOrder((char *)&commentLength, itmp[0]);
00503 ofile.write((char *)itmp, 4);
00504 }
00505
00506
00507 char cmt[1025];
00508 for(int i = 0; i < 1025; i++) cmt[i] = '\0';
00509
00510 const char * cmnt = kComment.c_str();
00511 size_t lcm = std::strlen(cmnt);
00512 if(lcm > 1024) lcm = 1024;
00513 std::strncpy(cmt, cmnt, lcm);
00514 ofile.write((char *)cmt, 1024);
00515 if(DEBUG || kVerbose > 0) {
00516 G4cout << "Data comment : "
00517 << kComment << G4endl;
00518 }
00519
00520
00521 if(kLittleEndianOutput) {
00522 ofile.write((char *)kVoxelSpacing, 12);
00523 } else {
00524 for(int j = 0; j < 3; j++)
00525 invertByteOrder((char *)&kVoxelSpacing[j], ftmp[j]);
00526 ofile.write((char *)ftmp, 12);
00527 }
00528 if(DEBUG || kVerbose > 0) {
00529 G4cout << "Voxel spacing : ("
00530 << kVoxelSpacing[0] << ", "
00531 << kVoxelSpacing[1] << ", "
00532 << kVoxelSpacing[2]
00533 << ") mm " << G4endl;
00534 }
00535
00536 calcPointers4();
00537 if(!kTracksWillBeStored) kPointerToTrackData = 0;
00538
00539
00540 if(kLittleEndianOutput) {
00541 ofile.write((char *)&kPointerToModalityData, 4);
00542 } else {
00543 invertByteOrder((char *)&kPointerToModalityData, itmp[0]);
00544 ofile.write((char *)itmp, 4);
00545 }
00546
00547
00548
00549 int nDoseDist = getNumDoseDist();
00550 if(kLittleEndianOutput) {
00551 ofile.write((char *)&nDoseDist, 4);
00552 } else {
00553 invertByteOrder((char *)&nDoseDist, itmp[0]);
00554 ofile.write((char *)itmp, 4);
00555 }
00556
00557
00558 if(kLittleEndianOutput) {
00559 for(int i = 0; i < nDoseDist; i++) {
00560 ofile.write((char *)&kPointerToDoseDistData[i], 4);
00561 }
00562 } else {
00563 for(int i = 0; i < nDoseDist; i++) {
00564 invertByteOrder((char *)&kPointerToDoseDistData[i], itmp[0]);
00565 ofile.write((char *)itmp, 4);
00566 }
00567 }
00568
00569
00570 if(kLittleEndianOutput) {
00571 ofile.write((char *)&kPointerToROIData, 4);
00572 } else {
00573 invertByteOrder((char *)&kPointerToROIData, itmp[0]);
00574 ofile.write((char *)itmp, 4);
00575 }
00576
00577
00578 if(kLittleEndianOutput) {
00579 ofile.write((char *)&kPointerToTrackData, 4);
00580 } else {
00581 invertByteOrder((char *)&kPointerToTrackData, itmp[0]);
00582 ofile.write((char *)itmp, 4);
00583 }
00584
00585
00586 if(kLittleEndianOutput) {
00587 ofile.write((char *)&kPointerToDetectorData, 4);
00588 } else {
00589 invertByteOrder((char *)&kPointerToDetectorData, itmp[0]);
00590 ofile.write((char *)itmp, 4);
00591 }
00592
00593 if(DEBUG || kVerbose > 0) {
00594 G4cout << "Each pointer to data : "
00595 << kPointerToModalityData << ", ";
00596 for(int i = 0; i < nDoseDist; i++) {
00597 G4cout << kPointerToDoseDistData[i] << ", ";
00598 }
00599 G4cout << kPointerToROIData << ", "
00600 << kPointerToTrackData << ", "
00601 << kPointerToDetectorData
00602 << G4endl;
00603 }
00604
00605
00606
00607 int size[3];
00608 float scale;
00609 short minmax[2];
00610 float fCenter[3];
00611 int iCenter[3];
00612
00613 kModality.getSize(size);
00614
00615 if(kLittleEndianOutput) {
00616 ofile.write((char *)size, 3*sizeof(int));
00617 } else {
00618 for(int j = 0; j < 3; j++)
00619 invertByteOrder((char *)&size[j], itmp[j]);
00620 ofile.write((char *)itmp, 12);
00621 }
00622
00623 if(DEBUG || kVerbose > 0) {
00624 G4cout << "Modality image size : ("
00625 << size[0] << ", "
00626 << size[1] << ", "
00627 << size[2] << ")"
00628 << G4endl;
00629 }
00630
00631
00632 kModality.getMinMax(minmax);
00633 if(kLittleEndianOutput) {
00634 ofile.write((char *)minmax, 4);
00635 } else {
00636 for(int j = 0; j < 2; j++)
00637 invertByteOrder((char *)&minmax[j], stmp[j]);
00638 ofile.write((char *)stmp, 4);
00639 }
00640
00641
00642 char munit[13] = "g/cm3\0";
00643 ofile.write((char *)munit, 12);
00644
00645
00646 scale = (float)kModality.getScale();
00647 if(kLittleEndianOutput) {
00648 ofile.write((char *)&scale, 4);
00649 } else {
00650 invertByteOrder((char *)&scale, ftmp[0]);
00651 ofile.write((char *)ftmp, 4);
00652 }
00653 if(DEBUG || kVerbose > 0) {
00654 G4cout << "Modality image min., max., scale : "
00655 << minmax[0] << ", "
00656 << minmax[1] << ", "
00657 << scale << G4endl;
00658 }
00659
00660
00661 int psize = size[0]*size[1];
00662 if(DEBUG || kVerbose > 0) G4cout << "Modality image : ";
00663 for(int i = 0; i < size[2]; i++) {
00664 short * image = kModality.getImage(i);
00665 if(kLittleEndianOutput) {
00666 ofile.write((char *)image, psize*sizeof(short));
00667 } else {
00668 for(int j = 0; j < psize; j++) {
00669 invertByteOrder((char *)&image[j], stmp[0]);
00670 ofile.write((char *)stmp, 2);
00671 }
00672 }
00673
00674 if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << image[(size_t)(psize*0.55)] << ", ";
00675 }
00676 if(DEBUG || kVerbose > 0) G4cout << G4endl;
00677
00678
00679 size_t msize = minmax[1] - minmax[0]+1;
00680 if(DEBUG || kVerbose > 0)
00681 G4cout << "modality image : " << minmax[0] << ", " << minmax[1] << G4endl;
00682 float * pdmap = new float[msize];
00683 for(int i = 0; i < (int)msize; i++) pdmap[i] =kModalityImageDensityMap[i];
00684
00685 if(kLittleEndianOutput) {
00686 ofile.write((char *)pdmap, msize*sizeof(float));
00687 } else {
00688 for(int j = 0; j < (int)msize; j++) {
00689 invertByteOrder((char *)&pdmap[j], ftmp[0]);
00690 ofile.write((char *)ftmp, 4);
00691 }
00692 }
00693
00694 if(DEBUG || kVerbose > 0) {
00695 G4cout << "density map : " << std::ends;
00696 for(int i = 0; i < (int)msize; i+=50)
00697 G4cout <<kModalityImageDensityMap[i] << ", ";
00698 G4cout << G4endl;
00699 }
00700 delete [] pdmap;
00701
00702
00703
00704
00705 if(!isDoseEmpty()) {
00706
00707 calcDoseDistScale();
00708
00709 for(int ndose = 0; ndose < nDoseDist; ndose++) {
00710
00711 kDose[ndose].getSize(size);
00712 if(kLittleEndianOutput) {
00713 ofile.write((char *)size, 3*sizeof(int));
00714 } else {
00715 for(int j = 0; j < 3; j++)
00716 invertByteOrder((char *)&size[j], itmp[j]);
00717 ofile.write((char *)itmp, 12);
00718 }
00719 if(DEBUG || kVerbose > 0) {
00720 G4cout << "Dose dist. [" << ndose << "] image size : ("
00721 << size[0] << ", "
00722 << size[1] << ", "
00723 << size[2] << ")"
00724 << G4endl;
00725 }
00726
00727
00728 getShortDoseDistMinMax(minmax, ndose);
00729 if(kLittleEndianOutput) {
00730 ofile.write((char *)minmax, 2*2);
00731 } else {
00732 for(int j = 0; j < 2; j++)
00733 invertByteOrder((char *)&minmax[j], stmp[j]);
00734 ofile.write((char *)stmp, 4);
00735 }
00736
00737
00738 char cdunit[13];
00739 for(int i = 0; i < 13; i++) cdunit[i] = '\0';
00740 const char * cu = kDoseUnit.c_str();
00741 size_t lcu = std::strlen(cu);
00742 if(lcu > 1024) lcu = 1024;
00743 std::strncpy(cdunit, cu, lcu);
00744 ofile.write((char *)cdunit, 12);
00745 if(DEBUG || kVerbose > 0) {
00746 G4cout << "Dose dist. unit : " << kDoseUnit << G4endl;
00747 }
00748
00749
00750 double dscale;
00751 dscale = getDoseDistScale(ndose);
00752 scale = float(dscale);
00753 if(kLittleEndianOutput) {
00754 ofile.write((char *)&scale, 4);
00755 } else {
00756 invertByteOrder((char *)&scale, ftmp[0]);
00757 ofile.write((char *)ftmp, 4);
00758 }
00759 if(DEBUG || kVerbose > 0) {
00760 G4cout << "Dose dist. [" << ndose
00761 << "] image min., max., scale : "
00762 << minmax[0] << ", "
00763 << minmax[1] << ", "
00764 << scale << G4endl;
00765 }
00766
00767
00768 int dsize = size[0]*size[1];
00769 short * dimage = new short[dsize];
00770 for(int z = 0; z < size[2]; z++) {
00771 getShortDoseDist(dimage, z, ndose);
00772 if(kLittleEndianOutput) {
00773 ofile.write((char *)dimage, dsize*2);
00774 } else {
00775 for(int j = 0; j < dsize; j++) {
00776 invertByteOrder((char *)&dimage[j], stmp[0]);
00777 ofile.write((char *)stmp, 2);
00778 }
00779 }
00780
00781 if(DEBUG || kVerbose > 0) {
00782 for(int j = 0; j < dsize; j++) {
00783 if(dimage[j] < 0)
00784 G4cout << "[" << j << "," << z << "]"
00785 << dimage[j] << ", ";
00786 }
00787 }
00788 }
00789 if(DEBUG || kVerbose > 0) G4cout << G4endl;
00790 delete [] dimage;
00791
00792
00793
00794 getDoseDistCenterPosition(fCenter, ndose);
00795 for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
00796 if(kLittleEndianOutput) {
00797 ofile.write((char *)iCenter, 3*4);
00798 } else {
00799 for(int j = 0; j < 3; j++)
00800 invertByteOrder((char *)&iCenter[j], itmp[j]);
00801 ofile.write((char *)itmp, 12);
00802 }
00803 if(DEBUG || kVerbose > 0) {
00804 G4cout << "Dose dist. [" << ndose
00805 << "]image relative location : ("
00806 << iCenter[0] << ", "
00807 << iCenter[1] << ", "
00808 << iCenter[2] << ")" << G4endl;
00809 }
00810
00811
00812 std::string name = getDoseDistName(ndose);
00813 if(name.size() == 0) name = "dose";
00814 name.resize(80);
00815 ofile.write((char *)name.c_str(), 80);
00816 if(DEBUG || kVerbose > 0) {
00817 G4cout << "Dose dist. name : " << name << G4endl;
00818 }
00819
00820 }
00821 }
00822
00823
00824 if(!isROIEmpty()) {
00825
00826 kRoi[0].getSize(size);
00827 if(kLittleEndianOutput) {
00828 ofile.write((char *)size, 3*sizeof(int));
00829 } else {
00830 for(int j = 0; j < 3; j++)
00831 invertByteOrder((char *)&size[j], itmp[j]);
00832 ofile.write((char *)itmp, 12);
00833 }
00834 if(DEBUG || kVerbose > 0) {
00835 G4cout << "ROI image size : ("
00836 << size[0] << ", "
00837 << size[1] << ", "
00838 << size[2] << ")"
00839 << G4endl;
00840 }
00841
00842
00843 kRoi[0].getMinMax(minmax);
00844 if(kLittleEndianOutput) {
00845 ofile.write((char *)minmax, sizeof(short)*2);
00846 } else {
00847 for(int j = 0; j < 2; j++)
00848 invertByteOrder((char *)&minmax[j], stmp[j]);
00849 ofile.write((char *)stmp, 4);
00850 }
00851
00852
00853 scale = (float)kRoi[0].getScale();
00854 if(kLittleEndianOutput) {
00855 ofile.write((char *)&scale, sizeof(float));
00856 } else {
00857 invertByteOrder((char *)&scale, ftmp[0]);
00858 ofile.write((char *)ftmp, 4);
00859 }
00860 if(DEBUG || kVerbose > 0) {
00861 G4cout << "ROI image min., max., scale : "
00862 << minmax[0] << ", "
00863 << minmax[1] << ", "
00864 << scale << G4endl;
00865 }
00866
00867
00868 int rsize = size[0]*size[1];
00869 for(int i = 0; i < size[2]; i++) {
00870 short * rimage = kRoi[0].getImage(i);
00871 if(kLittleEndianOutput) {
00872 ofile.write((char *)rimage, rsize*sizeof(short));
00873 } else {
00874 for(int j = 0; j < rsize; j++) {
00875 invertByteOrder((char *)&rimage[j], stmp[0]);
00876 ofile.write((char *)stmp, 2);
00877 }
00878 }
00879
00880 }
00881
00882
00883 kRoi[0].getCenterPosition(fCenter);
00884 for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
00885 if(kLittleEndianOutput) {
00886 ofile.write((char *)iCenter, 3*sizeof(int));
00887 } else {
00888 for(int j = 0; j < 3; j++)
00889 invertByteOrder((char *)&iCenter[j], itmp[j]);
00890 ofile.write((char *)itmp, 12);
00891 }
00892 if(DEBUG || kVerbose > 0) {
00893 G4cout << "ROI image relative location : ("
00894 << iCenter[0] << ", "
00895 << iCenter[1] << ", "
00896 << iCenter[2] << ")" << G4endl;
00897 }
00898 }
00899
00900
00901
00902 if(kPointerToTrackData > 0) {
00903
00904 int ntrk = kTracks.size();
00905 if(kLittleEndianOutput) {
00906 ofile.write((char *)&ntrk, sizeof(int));
00907 } else {
00908 invertByteOrder((char *)&ntrk, itmp[0]);
00909 ofile.write((char *)itmp, 4);
00910 }
00911 if(DEBUG || kVerbose > 0) {
00912 G4cout << "# of tracks : "
00913 << ntrk << G4endl;
00914 }
00915
00916 for(int nt = 0; nt < ntrk; nt++) {
00917
00918
00919 int nsteps = kTracks[nt].getNumberOfSteps();
00920 if(kLittleEndianOutput) {
00921 ofile.write((char *)&nsteps, sizeof(int));
00922 } else {
00923 invertByteOrder((char *)&nsteps, itmp[0]);
00924 ofile.write((char *)itmp, 4);
00925 }
00926 if(DEBUG || kVerbose > 0) {
00927 G4cout << "# of steps : " << nsteps << G4endl;
00928 }
00929
00930
00931 unsigned char tcolor[3];
00932 kTracks[nt].getColor(tcolor);
00933 ofile.write((char *)tcolor, 3);
00934
00935
00936 float stepPoints[6];
00937 for(int isteps = 0; isteps < nsteps; isteps++) {
00938 kTracks[nt].getStep(stepPoints[0], stepPoints[1], stepPoints[2],
00939 stepPoints[3], stepPoints[4], stepPoints[5],
00940 isteps);
00941
00942 if(kLittleEndianOutput) {
00943 ofile.write((char *)stepPoints, sizeof(float)*6);
00944 } else {
00945 for(int j = 0; j < 6; j++)
00946 invertByteOrder((char *)&stepPoints[j], ftmp[j]);
00947 ofile.write((char *)ftmp, 24);
00948 }
00949 }
00950 }
00951 }
00952
00953
00954
00955 if(kPointerToDetectorData > 0) {
00956 int ndet = kDetectors.size();
00957 if(kLittleEndianOutput) {
00958 ofile.write((char *)&ndet, sizeof(int));
00959 } else {
00960 invertByteOrder((char *)&ndet, itmp[0]);
00961 ofile.write((char *)itmp, 4);
00962 }
00963 if(DEBUG || kVerbose > 0) {
00964 G4cout << "# of detectors : "
00965 << ndet << G4endl;
00966 }
00967
00968 for(int nd = 0; nd < ndet; nd++) {
00969
00970
00971 int nedges = kDetectors[nd].getNumberOfEdges();
00972 if(kLittleEndianOutput) {
00973 ofile.write((char *)&nedges, sizeof(int));
00974 } else {
00975 invertByteOrder((char *)&nedges, itmp[0]);
00976 ofile.write((char *)itmp, 4);
00977 }
00978 if(DEBUG || kVerbose > 0) {
00979 G4cout << "# of edges in a detector : " << nedges << G4endl;
00980 }
00981
00982
00983 float edgePoints[6];
00984 for(int ne = 0; ne < nedges; ne++) {
00985 kDetectors[nd].getEdge(edgePoints[0], edgePoints[1], edgePoints[2],
00986 edgePoints[3], edgePoints[4], edgePoints[5],
00987 ne);
00988
00989 if(kLittleEndianOutput) {
00990 ofile.write((char *)edgePoints, sizeof(float)*6);
00991 } else {
00992 for(int j = 0; j < 6; j++)
00993 invertByteOrder((char *)&edgePoints[j], ftmp[j]);
00994 ofile.write((char *)ftmp, 24);
00995 }
00996
00997 if(DEBUG || kVerbose > 0) {
00998 if(ne < 1) {
00999 G4cout << " edge : (" << edgePoints[0] << ", "
01000 << edgePoints[1] << ", "
01001 << edgePoints[2] << ") - ("
01002 << edgePoints[3] << ", "
01003 << edgePoints[4] << ", "
01004 << edgePoints[5] << ")" << G4endl;
01005 }
01006 }
01007 }
01008
01009
01010 unsigned char dcolor[3];
01011 kDetectors[nd].getColor(dcolor);
01012 ofile.write((char *)dcolor, 3);
01013 if(DEBUG || kVerbose > 0) {
01014 G4cout << " rgb : (" << (int)dcolor[0] << ", "
01015 << (int)dcolor[1] << ", "
01016 << (int)dcolor[2] << ")" << G4endl;
01017 }
01018
01019
01020 std::string dname = kDetectors[nd].getName();
01021 dname.resize(80);
01022 ofile.write((char *)dname.c_str(), 80);
01023 if(DEBUG || kVerbose > 0) {
01024 G4cout << " detector name : " << dname << G4endl;
01025
01026 }
01027 }
01028 }
01029
01030
01031 ofile.write("END", 3);
01032
01033 ofile.close();
01034 if(DEBUG || kVerbose > 0)
01035 G4cout << ">>>> closed gdd file: " << kFileName << G4endl;
01036
01037 return true;
01038 }
01039 bool G4GMocrenIO::storeData3() {
01040
01041 if(kVerbose > 0) G4cout << ">>>>>>> store data (ver.3) <<<<<<<" << G4endl;
01042 if(kVerbose > 0) G4cout << " " << kFileName << G4endl;
01043
01044 bool DEBUG = false;
01045
01046
01047 std::ofstream ofile(kFileName.c_str(),
01048 std::ios_base::out|std::ios_base::binary);
01049
01050
01051 ofile.write("gMocren ", 8);
01052
01053
01054 unsigned char ver = 0x03;
01055 ofile.write((char *)&ver, 1);
01056
01057
01058 ofile.write((char *)&kLittleEndianOutput, sizeof(char));
01059
01060
01061 int commentLength = 1024;
01062 ofile.write((char *)&commentLength, 4);
01063
01064
01065 char cmt[1025];
01066 std::strncpy(cmt, kComment.c_str(), 1024);
01067 ofile.write((char *)cmt, 1024);
01068 if(DEBUG || kVerbose > 0) {
01069 G4cout << "Data comment : "
01070 << kComment << G4endl;
01071 }
01072
01073
01074 ofile.write((char *)kVoxelSpacing, 12);
01075 if(DEBUG || kVerbose > 0) {
01076 G4cout << "Voxel spacing : ("
01077 << kVoxelSpacing[0] << ", "
01078 << kVoxelSpacing[1] << ", "
01079 << kVoxelSpacing[2]
01080 << ") mm " << G4endl;
01081 }
01082
01083 calcPointers3();
01084
01085
01086 ofile.write((char *)&kPointerToModalityData, 4);
01087
01088
01089
01090 int nDoseDist = getNumDoseDist();
01091 ofile.write((char *)&nDoseDist, 4);
01092
01093
01094 for(int i = 0; i < nDoseDist; i++) {
01095 ofile.write((char *)&kPointerToDoseDistData[i], 4);
01096 }
01097
01098
01099 ofile.write((char *)&kPointerToROIData, 4);
01100
01101
01102 ofile.write((char *)&kPointerToTrackData, 4);
01103 if(DEBUG || kVerbose > 0) {
01104 G4cout << "Each pointer to data : "
01105 << kPointerToModalityData << ", ";
01106 for(int i = 0; i < nDoseDist; i++) {
01107 G4cout << kPointerToDoseDistData[i] << ", ";
01108 }
01109 G4cout << kPointerToROIData << ", "
01110 << kPointerToTrackData << G4endl;
01111 }
01112
01113
01114
01115 int size[3];
01116 float scale;
01117 short minmax[2];
01118 float fCenter[3];
01119 int iCenter[3];
01120
01121 kModality.getSize(size);
01122 ofile.write((char *)size, 3*sizeof(int));
01123 if(DEBUG || kVerbose > 0) {
01124 G4cout << "Modality image size : ("
01125 << size[0] << ", "
01126 << size[1] << ", "
01127 << size[2] << ")"
01128 << G4endl;
01129 }
01130
01131
01132 kModality.getMinMax(minmax);
01133 ofile.write((char *)minmax, 4);
01134
01135
01136 char munit[13] = "g/cm3 ";
01137 ofile.write((char *)munit, 12);
01138
01139
01140 scale = (float)kModality.getScale();
01141 ofile.write((char *)&scale, 4);
01142 if(DEBUG || kVerbose > 0) {
01143 G4cout << "Modality image min., max., scale : "
01144 << minmax[0] << ", "
01145 << minmax[1] << ", "
01146 << scale << G4endl;
01147 }
01148
01149
01150 int psize = size[0]*size[1];
01151 if(DEBUG || kVerbose > 0) G4cout << "Modality image : ";
01152 for(int i = 0; i < size[2]; i++) {
01153 short * image = kModality.getImage(i);
01154 ofile.write((char *)image, psize*sizeof(short));
01155
01156 if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << image[(size_t)(psize*0.55)] << ", ";
01157 }
01158 if(DEBUG || kVerbose > 0) G4cout << G4endl;
01159
01160
01161 size_t msize = minmax[1] - minmax[0]+1;
01162 float * pdmap = new float[msize];
01163 for(int i = 0; i < (int)msize; i++) pdmap[i] =kModalityImageDensityMap[i];
01164 ofile.write((char *)pdmap, msize*sizeof(float));
01165 if(DEBUG || kVerbose > 0) {
01166 G4cout << "density map : " << std::ends;
01167 for(int i = 0; i < (int)msize; i+=50)
01168 G4cout <<kModalityImageDensityMap[i] << ", ";
01169 G4cout << G4endl;
01170 }
01171 delete [] pdmap;
01172
01173
01174
01175
01176 if(!isDoseEmpty()) {
01177
01178 calcDoseDistScale();
01179
01180 for(int ndose = 0; ndose < nDoseDist; ndose++) {
01181
01182 kDose[ndose].getSize(size);
01183 ofile.write((char *)size, 3*sizeof(int));
01184 if(DEBUG || kVerbose > 0) {
01185 G4cout << "Dose dist. [" << ndose << "] image size : ("
01186 << size[0] << ", "
01187 << size[1] << ", "
01188 << size[2] << ")"
01189 << G4endl;
01190 }
01191
01192
01193 getShortDoseDistMinMax(minmax, ndose);
01194 ofile.write((char *)minmax, 2*2);
01195
01196
01197 ofile.write((char *)kDoseUnit.c_str(), 12);
01198 if(DEBUG || kVerbose > 0) {
01199 G4cout << "Dose dist. unit : " << kDoseUnit << G4endl;
01200 }
01201
01202
01203 double dscale;
01204 dscale = getDoseDistScale(ndose);
01205 scale = float(dscale);
01206 ofile.write((char *)&scale, 4);
01207 if(DEBUG || kVerbose > 0) {
01208 G4cout << "Dose dist. [" << ndose
01209 << "] image min., max., scale : "
01210 << minmax[0] << ", "
01211 << minmax[1] << ", "
01212 << scale << G4endl;
01213 }
01214
01215
01216 int dsize = size[0]*size[1];
01217 short * dimage = new short[dsize];
01218 for(int z = 0; z < size[2]; z++) {
01219 getShortDoseDist(dimage, z, ndose);
01220 ofile.write((char *)dimage, dsize*2);
01221
01222 if(DEBUG || kVerbose > 0) {
01223 for(int j = 0; j < dsize; j++) {
01224 if(dimage[j] < 0)
01225 G4cout << "[" << j << "," << z << "]"
01226 << dimage[j] << ", ";
01227 }
01228 }
01229 }
01230 if(DEBUG || kVerbose > 0) G4cout << G4endl;
01231 delete [] dimage;
01232
01233
01234
01235 getDoseDistCenterPosition(fCenter, ndose);
01236 for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
01237 ofile.write((char *)iCenter, 3*4);
01238 if(DEBUG || kVerbose > 0) {
01239 G4cout << "Dose dist. [" << ndose
01240 << "]image relative location : ("
01241 << iCenter[0] << ", "
01242 << iCenter[1] << ", "
01243 << iCenter[2] << ")" << G4endl;
01244 }
01245 }
01246 }
01247
01248
01249 if(!isROIEmpty()) {
01250
01251 kRoi[0].getSize(size);
01252 ofile.write((char *)size, 3*sizeof(int));
01253 if(DEBUG || kVerbose > 0) {
01254 G4cout << "ROI image size : ("
01255 << size[0] << ", "
01256 << size[1] << ", "
01257 << size[2] << ")"
01258 << G4endl;
01259 }
01260
01261
01262 kRoi[0].getMinMax(minmax);
01263 ofile.write((char *)minmax, sizeof(short)*2);
01264
01265
01266 scale = (float)kRoi[0].getScale();
01267 ofile.write((char *)&scale, sizeof(float));
01268 if(DEBUG || kVerbose > 0) {
01269 G4cout << "ROI image min., max., scale : "
01270 << minmax[0] << ", "
01271 << minmax[1] << ", "
01272 << scale << G4endl;
01273 }
01274
01275
01276 int rsize = size[0]*size[1];
01277 for(int i = 0; i < size[2]; i++) {
01278 short * rimage = kRoi[0].getImage(i);
01279 ofile.write((char *)rimage, rsize*sizeof(short));
01280
01281 }
01282
01283
01284 kRoi[0].getCenterPosition(fCenter);
01285 for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
01286 ofile.write((char *)iCenter, 3*sizeof(int));
01287 if(DEBUG || kVerbose > 0) {
01288 G4cout << "ROI image relative location : ("
01289 << iCenter[0] << ", "
01290 << iCenter[1] << ", "
01291 << iCenter[2] << ")" << G4endl;
01292 }
01293 }
01294
01295
01296
01297 int ntrk = kSteps.size();
01298 ofile.write((char *)&ntrk, sizeof(int));
01299 if(DEBUG || kVerbose > 0) {
01300 G4cout << "# of tracks : "
01301 << ntrk << G4endl;
01302 }
01303
01304 for(int i = 0; i < ntrk; i++) {
01305 float * tp = kSteps[i];
01306 ofile.write((char *)tp, sizeof(float)*6);
01307 }
01308
01309 int ntcolor = int(kStepColors.size());
01310 if(ntrk != ntcolor)
01311 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
01312 G4cout << "# of track color information must be the same as # of tracks."
01313 << G4endl;
01314 unsigned char white[3] = {255,255,255};
01315 for(int i = 0; i < ntrk; i++) {
01316 if(i < ntcolor) {
01317 unsigned char * tcolor = kStepColors[i];
01318 ofile.write((char *)tcolor, 3);
01319 } else {
01320 ofile.write((char *)white, 3);
01321 }
01322 }
01323
01324
01325 ofile.write("END", 3);
01326
01327 ofile.close();
01328
01329 return true;
01330 }
01331
01332 bool G4GMocrenIO::storeData4(char * _filename) {
01333 kFileName = _filename;
01334 return storeData4();
01335 }
01336
01337
01338 bool G4GMocrenIO::storeData2() {
01339
01340 if(kVerbose > 0) G4cout << ">>>>>>> store data (ver.2) <<<<<<<" << G4endl;
01341 if(kVerbose > 0) G4cout << " " << kFileName << G4endl;
01342
01343 bool DEBUG = false;
01344
01345
01346 std::ofstream ofile(kFileName.c_str(),
01347 std::ios_base::out|std::ios_base::binary);
01348
01349
01350 ofile.write("GRAPE ", 8);
01351
01352
01353 unsigned char ver = 0x02;
01354 ofile.write((char *)&ver, 1);
01355
01356 ofile.write(kId.c_str(), IDLENGTH);
01357
01358 ofile.write(kVersion.c_str(), VERLENGTH);
01359
01360 ofile.write((char *)&kLittleEndianOutput, sizeof(char));
01361
01362
01363
01364
01365
01366
01367
01368
01369
01370
01371
01372
01373
01374 ofile.write((char *)kVoxelSpacing, 12);
01375 if(DEBUG || kVerbose > 0) {
01376 G4cout << "Voxel spacing : ("
01377 << kVoxelSpacing[0] << ", "
01378 << kVoxelSpacing[1] << ", "
01379 << kVoxelSpacing[2]
01380 << ") mm " << G4endl;
01381 }
01382
01383 calcPointers2();
01384
01385 ofile.write((char *)&kPointerToModalityData, 4);
01386
01387
01388 ofile.write((char *)&kPointerToDoseDistData[0], 4);
01389
01390
01391 ofile.write((char *)&kPointerToROIData, 4);
01392
01393
01394 ofile.write((char *)&kPointerToTrackData, 4);
01395 if(DEBUG || kVerbose > 0) {
01396 G4cout << "Each pointer to data : "
01397 << kPointerToModalityData << ", "
01398 << kPointerToDoseDistData[0] << ", "
01399 << kPointerToROIData << ", "
01400 << kPointerToTrackData << G4endl;
01401 }
01402
01403
01404
01405 int size[3];
01406 float scale;
01407 short minmax[2];
01408 float fCenter[3];
01409 int iCenter[3];
01410
01411 kModality.getSize(size);
01412 ofile.write((char *)size, 3*sizeof(int));
01413 if(DEBUG || kVerbose > 0) {
01414 G4cout << "Modality image size : ("
01415 << size[0] << ", "
01416 << size[1] << ", "
01417 << size[2] << ")"
01418 << G4endl;
01419 }
01420
01421
01422 kModality.getMinMax(minmax);
01423 ofile.write((char *)minmax, 4);
01424
01425
01426
01427
01428
01429
01430 scale = (float)kModality.getScale();
01431 ofile.write((char *)&scale, 4);
01432 if(DEBUG || kVerbose > 0) {
01433 G4cout << "Modality image min., max., scale : "
01434 << minmax[0] << ", "
01435 << minmax[1] << ", "
01436 << scale << G4endl;
01437 }
01438
01439
01440 int psize = size[0]*size[1];
01441 if(DEBUG || kVerbose > 0) G4cout << "Modality image : ";
01442 for(int i = 0; i < size[2]; i++) {
01443 short * image =kModality.getImage(i);
01444 ofile.write((char *)image, psize*sizeof(short));
01445
01446 if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << image[(size_t)(psize*0.55)] << ", ";
01447 }
01448 if(DEBUG || kVerbose > 0) G4cout << G4endl;
01449
01450
01451 size_t msize = minmax[1] - minmax[0]+1;
01452 float * pdmap = new float[msize];
01453 for(int i = 0; i < (int)msize; i++) pdmap[i] =kModalityImageDensityMap[i];
01454 ofile.write((char *)pdmap, msize*sizeof(float));
01455 if(DEBUG || kVerbose > 0) {
01456 G4cout << "density map : " << std::ends;
01457 for(int i = 0; i < (int)msize; i+=50)
01458 G4cout <<kModalityImageDensityMap[i] << ", ";
01459 G4cout << G4endl;
01460 }
01461 delete [] pdmap;
01462
01463
01464
01465
01466 if(!isDoseEmpty()) {
01467 calcDoseDistScale();
01468
01469
01470 kDose[0].getSize(size);
01471 ofile.write((char *)size, 3*sizeof(int));
01472 if(DEBUG || kVerbose > 0) {
01473 G4cout << "Dose dist. image size : ("
01474 << size[0] << ", "
01475 << size[1] << ", "
01476 << size[2] << ")"
01477 << G4endl;
01478 }
01479
01480
01481 getShortDoseDistMinMax(minmax);
01482 ofile.write((char *)minmax, sizeof(short)*2);
01483
01484
01485 scale = (float)kDose[0].getScale();
01486 ofile.write((char *)&scale, sizeof(float));
01487 if(DEBUG || kVerbose > 0) {
01488 G4cout << "Dose dist. image min., max., scale : "
01489 << minmax[0] << ", "
01490 << minmax[1] << ", "
01491 << scale << G4endl;
01492 }
01493
01494
01495 int dsize = size[0]*size[1];
01496 short * dimage = new short[dsize];
01497 for(int z = 0; z < size[2]; z++) {
01498 getShortDoseDist(dimage, z);
01499 ofile.write((char *)dimage, dsize*sizeof(short));
01500
01501 if(DEBUG || kVerbose > 0) {
01502 for(int j = 0; j < dsize; j++) {
01503 if(dimage[j] < 0)
01504 G4cout << "[" << j << "," << z << "]"
01505 << dimage[j] << ", ";
01506 }
01507 }
01508 }
01509 if(DEBUG || kVerbose > 0) G4cout << G4endl;
01510 delete [] dimage;
01511
01512
01513
01514 kDose[0].getCenterPosition(fCenter);
01515 for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
01516 ofile.write((char *)iCenter, 3*sizeof(int));
01517 if(DEBUG || kVerbose > 0) {
01518 G4cout << "Dose dist. image relative location : ("
01519 << iCenter[0] << ", "
01520 << iCenter[1] << ", "
01521 << iCenter[2] << ")" << G4endl;
01522 }
01523
01524 }
01525
01526
01527 if(!isROIEmpty()) {
01528
01529 kRoi[0].getSize(size);
01530 ofile.write((char *)size, 3*sizeof(int));
01531 if(DEBUG || kVerbose > 0) {
01532 G4cout << "ROI image size : ("
01533 << size[0] << ", "
01534 << size[1] << ", "
01535 << size[2] << ")"
01536 << G4endl;
01537 }
01538
01539
01540 kRoi[0].getMinMax(minmax);
01541 ofile.write((char *)minmax, sizeof(short)*2);
01542
01543
01544 scale = (float)kRoi[0].getScale();
01545 ofile.write((char *)&scale, sizeof(float));
01546 if(DEBUG || kVerbose > 0) {
01547 G4cout << "ROI image min., max., scale : "
01548 << minmax[0] << ", "
01549 << minmax[1] << ", "
01550 << scale << G4endl;
01551 }
01552
01553
01554 int rsize = size[0]*size[1];
01555 for(int i = 0; i < size[2]; i++) {
01556 short * rimage = kRoi[0].getImage(i);
01557 ofile.write((char *)rimage, rsize*sizeof(short));
01558
01559 }
01560
01561
01562 kRoi[0].getCenterPosition(fCenter);
01563 for(int i = 0; i < 3; i++) iCenter[i] = (int)fCenter[i];
01564 ofile.write((char *)iCenter, 3*sizeof(int));
01565 if(DEBUG || kVerbose > 0) {
01566 G4cout << "ROI image relative location : ("
01567 << iCenter[0] << ", "
01568 << iCenter[1] << ", "
01569 << iCenter[2] << ")" << G4endl;
01570 }
01571 }
01572
01573
01574
01575
01576 int ntrk = kSteps.size();
01577 ofile.write((char *)&ntrk, sizeof(int));
01578 if(DEBUG || kVerbose > 0) {
01579 G4cout << "# of tracks : "
01580 << ntrk << G4endl;
01581 }
01582 for(int i = 0; i < ntrk; i++) {
01583 float * tp = kSteps[i];
01584 ofile.write((char *)tp, sizeof(float)*6);
01585 }
01586
01587
01588
01589 ofile.write("END", 3);
01590
01591 ofile.close();
01592
01593 return true;
01594 }
01595
01596 bool G4GMocrenIO::storeData2(char * _filename) {
01597 kFileName = _filename;
01598 return storeData();
01599 }
01600
01601 bool G4GMocrenIO::retrieveData() {
01602
01603
01604 std::ifstream ifile(kFileName.c_str(), std::ios_base::in|std::ios_base::binary);
01605 if(!ifile) {
01606 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
01607 G4cout << "Cannot open file: " << kFileName
01608 << " in G4GMocrenIO::retrieveData()." << G4endl;
01609 return false;
01610 }
01611
01612
01613 char verid[9];
01614 ifile.read((char *)verid, 8);
01615
01616 unsigned char ver;
01617 ifile.read((char *)&ver, 1);
01618 ifile.close();
01619
01620 if(std::strncmp(verid, "gMocren", 7) == 0) {
01621 if(ver == 0x03) {
01622 G4cout << ">>>>>>> retrieve data (ver.3) <<<<<<<" << G4endl;
01623 G4cout << " " << kFileName << G4endl;
01624 retrieveData3();
01625 } else if (ver == 0x04) {
01626 G4cout << ">>>>>>> retrieve data (ver.4) <<<<<<<" << G4endl;
01627 G4cout << " " << kFileName << G4endl;
01628 retrieveData4();
01629 } else {
01630 if (G4VisManager::GetVerbosity() >= G4VisManager::errors) {
01631 G4cout << "Error -- invalid file version : " << (int)ver
01632 << G4endl;
01633 G4cout << " " << kFileName << G4endl;
01634 }
01635 std::exit(-1);
01636 }
01637 } else if(std::strncmp(verid, "GRAPE", 5) == 0) {
01638 G4cout << ">>>>>>> retrieve data (ver.2) <<<<<<<" << G4endl;
01639 G4cout << " " << kFileName << G4endl;
01640 retrieveData2();
01641 } else {
01642 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
01643 G4cout << kFileName << " was not gdd file." << G4endl;
01644 return false;
01645 }
01646
01647 return true;
01648 }
01649
01650 bool G4GMocrenIO::retrieveData(char * _filename) {
01651 kFileName = _filename;
01652 return retrieveData();
01653 }
01654
01655
01656 bool G4GMocrenIO::retrieveData4() {
01657
01658 bool DEBUG = false;
01659
01660
01661 std::ifstream ifile(kFileName.c_str(), std::ios_base::in|std::ios_base::binary);
01662 if(!ifile) {
01663 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
01664 G4cout << "Cannot open file: " << kFileName
01665 << " in G4GMocrenIO::retrieveData3()." << G4endl;
01666 return false;
01667 }
01668
01669
01670 char ctmp[24];
01671
01672
01673 char verid[9];
01674 ifile.read((char *)verid, 8);
01675
01676
01677 unsigned char ver;
01678 ifile.read((char *)&ver, 1);
01679 std::stringstream ss;
01680 ss << (int)ver;
01681 kVersion = ss.str();
01682 if(DEBUG || kVerbose > 0) G4cout << "File version : " << kVersion << G4endl;
01683
01684
01685 ifile.read((char *)&kLittleEndianInput, sizeof(char));
01686 if(DEBUG || kVerbose > 0) {
01687 G4cout << "Endian : ";
01688 if(kLittleEndianInput == 1)
01689 G4cout << " little" << G4endl;
01690 else {
01691 G4cout << " big" << G4endl;
01692 }
01693 }
01694
01695
01696 int clength;
01697 ifile.read((char *)ctmp, 4);
01698 convertEndian(ctmp, clength);
01699
01700 char cmt[1025];
01701 ifile.read((char *)cmt, clength);
01702 std::string scmt = cmt;
01703 scmt += '\0';
01704 setComment(scmt);
01705 if(DEBUG || kVerbose > 0) {
01706 G4cout << "Data comment : "
01707 << kComment << G4endl;
01708 }
01709
01710
01711 ifile.read((char *)ctmp, 12);
01712 convertEndian(ctmp, kVoxelSpacing[0]);
01713 convertEndian(ctmp+4, kVoxelSpacing[1]);
01714 convertEndian(ctmp+8, kVoxelSpacing[2]);
01715 if(DEBUG || kVerbose > 0) {
01716 G4cout << "Voxel spacing : ("
01717 << kVoxelSpacing[0] << ", "
01718 << kVoxelSpacing[1] << ", "
01719 << kVoxelSpacing[2]
01720 << ") mm " << G4endl;
01721 }
01722
01723
01724
01725 ifile.read((char *)ctmp, 4);
01726 convertEndian(ctmp, kPointerToModalityData);
01727
01728
01729 ifile.read((char *)ctmp, 4);
01730 int nDoseDist;
01731 convertEndian(ctmp, nDoseDist);
01732
01733
01734 for(int i = 0; i < nDoseDist; i++) {
01735 ifile.read((char *)ctmp, 4);
01736 unsigned int dptr;
01737 convertEndian(ctmp, dptr);
01738 addPointerToDoseDistData(dptr);
01739 }
01740
01741
01742 ifile.read((char *)ctmp, 4);
01743 convertEndian(ctmp, kPointerToROIData);
01744
01745
01746 ifile.read((char *)ctmp, 4);
01747 convertEndian(ctmp, kPointerToTrackData);
01748
01749
01750 ifile.read((char *)ctmp, 4);
01751 convertEndian(ctmp, kPointerToDetectorData);
01752
01753 if(DEBUG || kVerbose > 0) {
01754 G4cout << "Each pointer to data : "
01755 << kPointerToModalityData << ", ";
01756 for(int i = 0; i < nDoseDist; i++)
01757 G4cout << kPointerToDoseDistData[i] << ", ";
01758 G4cout << kPointerToROIData << ", "
01759 << kPointerToTrackData << ", "
01760 << kPointerToDetectorData
01761 << G4endl;
01762 }
01763
01764
01765
01766 if(kPointerToModalityData == 0 && kPointerToDoseDistData.size() == 0 &&
01767 kPointerToROIData == 0 && kPointerToTrackData == 0) {
01768 if(DEBUG || kVerbose > 0) {
01769 G4cout << "No data." << G4endl;
01770 }
01771 return false;
01772 }
01773
01774
01775
01776
01777
01778
01779
01780 int size[3];
01781 float scale;
01782 double dscale;
01783 short minmax[2];
01784 float fCenter[3];
01785 int iCenter[3];
01786
01787
01788
01789 ifile.read(ctmp, 3*sizeof(int));
01790 convertEndian(ctmp, size[0]);
01791 convertEndian(ctmp+sizeof(int), size[1]);
01792 convertEndian(ctmp+2*sizeof(int), size[2]);
01793 if(DEBUG || kVerbose > 0) {
01794 G4cout << "Modality image size : ("
01795 << size[0] << ", "
01796 << size[1] << ", "
01797 << size[2] << ")"
01798 << G4endl;
01799 }
01800 kModality.setSize(size);
01801
01802
01803
01804
01805
01806
01807
01808
01809
01810 if(kPointerToModalityData != 0) {
01811
01812
01813 ifile.read((char *)ctmp, 4);
01814 convertEndian(ctmp, minmax[0]);
01815 convertEndian(ctmp+2, minmax[1]);
01816 kModality.setMinMax(minmax);
01817
01818
01819 char munit[13];
01820 munit[12] = '\0';
01821 ifile.read((char *)munit, 12);
01822 std::string smunit = munit;
01823 setModalityImageUnit(smunit);
01824
01825
01826 ifile.read((char *)ctmp, 4);
01827 convertEndian(ctmp, scale);
01828 kModality.setScale(dscale = scale);
01829 if(DEBUG || kVerbose > 0) {
01830 G4cout << "Modality image min., max., scale : "
01831 << minmax[0] << ", "
01832 << minmax[1] << ", "
01833 << scale << G4endl;
01834 }
01835
01836
01837 int psize = size[0]*size[1];
01838 if(DEBUG || kVerbose > 0) G4cout << "Modality image (" << psize << "): ";
01839 char * cimage = new char[psize*sizeof(short)];
01840 for(int i = 0; i < size[2]; i++) {
01841 ifile.read((char *)cimage, psize*sizeof(short));
01842 short * mimage = new short[psize];
01843 for(int j = 0; j < psize; j++) {
01844 convertEndian(cimage+j*sizeof(short), mimage[j]);
01845 }
01846 kModality.addImage(mimage);
01847
01848 if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << mimage[(size_t)(psize*0.55)] << ", ";
01849 }
01850 if(DEBUG || kVerbose > 0) G4cout << G4endl;
01851 delete [] cimage;
01852
01853
01854 size_t msize = minmax[1]-minmax[0]+1;
01855 if(DEBUG || kVerbose > 0) G4cout << "msize: " << msize << G4endl;
01856 char * pdmap = new char[msize*sizeof(float)];
01857 ifile.read((char *)pdmap, msize*sizeof(float));
01858 float ftmp;
01859 for(int i = 0; i < (int)msize; i++) {
01860 convertEndian(pdmap+i*sizeof(float), ftmp);
01861 kModalityImageDensityMap.push_back(ftmp);
01862 }
01863 delete [] pdmap;
01864
01865 if(DEBUG || kVerbose > 0) {
01866 G4cout << "density map : " << std::ends;
01867 for(int i = 0; i < 10; i++)
01868 G4cout <<kModalityImageDensityMap[i] << ", ";
01869 G4cout << G4endl;
01870 for(int i = 0; i < 10; i++) G4cout << "..";
01871 G4cout << G4endl;
01872 for(size_t i =kModalityImageDensityMap.size() - 10; i <kModalityImageDensityMap.size(); i++)
01873 G4cout <<kModalityImageDensityMap[i] << ", ";
01874 G4cout << G4endl;
01875 }
01876
01877 }
01878
01879
01880
01881 for(int ndose = 0; ndose < nDoseDist; ndose++) {
01882
01883 newDoseDist();
01884
01885
01886 ifile.read((char *)ctmp, 3*sizeof(int));
01887 convertEndian(ctmp, size[0]);
01888 convertEndian(ctmp+sizeof(int), size[1]);
01889 convertEndian(ctmp+2*sizeof(int), size[2]);
01890 if(DEBUG || kVerbose > 0) {
01891 G4cout << "Dose dist. image size : ("
01892 << size[0] << ", "
01893 << size[1] << ", "
01894 << size[2] << ")"
01895 << G4endl;
01896 }
01897 kDose[ndose].setSize(size);
01898
01899
01900 ifile.read((char *)ctmp, sizeof(short)*2);
01901 convertEndian(ctmp, minmax[0]);
01902 convertEndian(ctmp+2, minmax[1]);
01903
01904
01905 char dunit[13];
01906 dunit[12] = '\0';
01907 ifile.read((char *)dunit, 12);
01908 std::string sdunit = dunit;
01909 setDoseDistUnit(sdunit, ndose);
01910 if(DEBUG || kVerbose > 0) {
01911 G4cout << "Dose dist. unit : " << kDoseUnit << G4endl;
01912 }
01913
01914
01915 ifile.read((char *)ctmp, 4);
01916 convertEndian(ctmp, scale);
01917 kDose[ndose].setScale(dscale = scale);
01918
01919 double dminmax[2];
01920 for(int i = 0; i < 2; i++) dminmax[i] = minmax[i]*dscale;
01921 kDose[ndose].setMinMax(dminmax);
01922
01923 if(DEBUG || kVerbose > 0) {
01924 G4cout << "Dose dist. image min., max., scale : "
01925 << dminmax[0] << ", "
01926 << dminmax[1] << ", "
01927 << scale << G4endl;
01928 }
01929
01930
01931 int dsize = size[0]*size[1];
01932 if(DEBUG || kVerbose > 0) G4cout << "Dose dist. (" << dsize << "): ";
01933 char * di = new char[dsize*sizeof(short)];
01934 short * shimage = new short[dsize];
01935 for(int z = 0; z < size[2]; z++) {
01936 ifile.read((char *)di, dsize*sizeof(short));
01937 double * dimage = new double[dsize];
01938 for(int xy = 0; xy < dsize; xy++) {
01939 convertEndian(di+xy*sizeof(short), shimage[xy]);
01940 dimage[xy] = shimage[xy]*dscale;
01941 }
01942 kDose[ndose].addImage(dimage);
01943
01944 if(DEBUG || kVerbose > 0) G4cout << "[" << z << "]" << dimage[(size_t)(dsize*0.55)] << ", ";
01945
01946 if(DEBUG || kVerbose > 0) {
01947 for(int j = 0; j < dsize; j++) {
01948 if(dimage[j] < 0)
01949 G4cout << "[" << j << "," << z << "]"
01950 << dimage[j] << ", ";
01951 }
01952 }
01953 }
01954 delete [] shimage;
01955 delete [] di;
01956 if(DEBUG || kVerbose > 0) G4cout << G4endl;
01957
01958 ifile.read((char *)ctmp, 3*4);
01959 convertEndian(ctmp, iCenter[0]);
01960 convertEndian(ctmp+4, iCenter[1]);
01961 convertEndian(ctmp+8, iCenter[2]);
01962 for(int i = 0; i < 3; i++) fCenter[i] = (float)iCenter[i];
01963 kDose[ndose].setCenterPosition(fCenter);
01964
01965 if(DEBUG || kVerbose > 0) {
01966 G4cout << "Dose dist. image relative location : ("
01967 << fCenter[0] << ", "
01968 << fCenter[1] << ", "
01969 << fCenter[2] << ")" << G4endl;
01970 }
01971
01972
01973
01974 char cname[81];
01975 ifile.read((char *)cname, 80);
01976 std::string dosename = cname;
01977 setDoseDistName(dosename, ndose);
01978 if(DEBUG || kVerbose > 0) {
01979 G4cout << "Dose dist. name : " << dosename << G4endl;
01980 }
01981
01982 }
01983
01984
01985 if(kPointerToROIData != 0) {
01986
01987 newROI();
01988
01989
01990 ifile.read((char *)ctmp, 3*sizeof(int));
01991 convertEndian(ctmp, size[0]);
01992 convertEndian(ctmp+sizeof(int), size[1]);
01993 convertEndian(ctmp+2*sizeof(int), size[2]);
01994 kRoi[0].setSize(size);
01995 if(DEBUG || kVerbose > 0) {
01996 G4cout << "ROI image size : ("
01997 << size[0] << ", "
01998 << size[1] << ", "
01999 << size[2] << ")"
02000 << G4endl;
02001 }
02002
02003
02004 ifile.read((char *)ctmp, sizeof(short)*2);
02005 convertEndian(ctmp, minmax[0]);
02006 convertEndian(ctmp+sizeof(short), minmax[1]);
02007 kRoi[0].setMinMax(minmax);
02008
02009
02010 ifile.read((char *)ctmp, sizeof(float));
02011 convertEndian(ctmp, scale);
02012 kRoi[0].setScale(dscale = scale);
02013 if(DEBUG || kVerbose > 0) {
02014 G4cout << "ROI image min., max., scale : "
02015 << minmax[0] << ", "
02016 << minmax[1] << ", "
02017 << scale << G4endl;
02018 }
02019
02020
02021 int rsize = size[0]*size[1];
02022 char * ri = new char[rsize*sizeof(short)];
02023 for(int i = 0; i < size[2]; i++) {
02024 ifile.read((char *)ri, rsize*sizeof(short));
02025 short * rimage = new short[rsize];
02026 for(int j = 0; j < rsize; j++) {
02027 convertEndian(ri+j*sizeof(short), rimage[j]);
02028 }
02029 kRoi[0].addImage(rimage);
02030
02031 }
02032 delete [] ri;
02033
02034
02035 ifile.read((char *)ctmp, 3*sizeof(int));
02036 convertEndian(ctmp, iCenter[0]);
02037 convertEndian(ctmp+sizeof(int), iCenter[1]);
02038 convertEndian(ctmp+2*sizeof(int), iCenter[2]);
02039 for(int i = 0; i < 3; i++) fCenter[i] = iCenter[i];
02040 kRoi[0].setCenterPosition(fCenter);
02041 if(DEBUG || kVerbose > 0) {
02042 G4cout << "ROI image relative location : ("
02043 << fCenter[0] << ", "
02044 << fCenter[1] << ", "
02045 << fCenter[2] << ")" << G4endl;
02046 }
02047
02048 }
02049
02050
02051 if(kPointerToTrackData != 0) {
02052
02053
02054 ifile.read((char *)ctmp, sizeof(int));
02055 int ntrk;
02056 convertEndian(ctmp, ntrk);
02057 if(DEBUG || kVerbose > 0) {
02058 G4cout << "# of tracks: " << ntrk << G4endl;
02059 }
02060
02061
02062 unsigned char rgb[3];
02063 for(int i = 0; i < ntrk; i++) {
02064
02065
02066
02067 ifile.read((char *)ctmp, sizeof(int));
02068 int nsteps;
02069 convertEndian(ctmp, nsteps);
02070
02071
02072 ifile.read((char *)rgb, 3);
02073
02074 std::vector<float *> steps;
02075
02076 for(int j = 0; j < nsteps; j++) {
02077
02078 float * steppoint = new float[6];
02079 ifile.read((char *)ctmp, sizeof(float)*6);
02080
02081 for(int k = 0; k < 6; k++) {
02082 convertEndian(ctmp+k*sizeof(float), steppoint[k]);
02083 }
02084
02085 steps.push_back(steppoint);
02086 }
02087
02088
02089 addTrack(steps, rgb);
02090
02091 if(DEBUG || kVerbose > 0) {
02092 if(i < 5) {
02093 G4cout << i << ": " ;
02094 for(int j = 0; j < 3; j++) G4cout << steps[0][j] << " ";
02095 int nstp = steps.size();
02096 G4cout << "<-> ";
02097 for(int j = 3; j < 6; j++) G4cout << steps[nstp-1][j] << " ";
02098 G4cout << " rgb( ";
02099 for(int j = 0; j < 3; j++) G4cout << (int)rgb[j] << " ";
02100 G4cout << ")" << G4endl;
02101 }
02102 }
02103 }
02104
02105
02106 }
02107
02108
02109
02110 if(kPointerToDetectorData != 0) {
02111
02112
02113 ifile.read((char *)ctmp, sizeof(int));
02114 int ndet;
02115 convertEndian(ctmp, ndet);
02116
02117 if(DEBUG || kVerbose > 0) {
02118 G4cout << "# of detectors : "
02119 << ndet << G4endl;
02120 }
02121
02122 for(int nd = 0; nd < ndet; nd++) {
02123
02124
02125 ifile.read((char *)ctmp, sizeof(int));
02126 int nedges;
02127 convertEndian(ctmp, nedges);
02128 if(DEBUG || kVerbose > 0) {
02129 G4cout << "# of edges in a detector : " << nedges << G4endl;
02130 }
02131
02132
02133 std::vector<float *> detector;
02134 char cftmp[24];
02135 for(int ne = 0; ne < nedges; ne++) {
02136
02137 ifile.read((char *)cftmp, sizeof(float)*6);
02138 float * edgePoints = new float[6];
02139 for(int j = 0; j < 6; j++) convertEndian(&cftmp[sizeof(float)*j], edgePoints[j]);
02140 detector.push_back(edgePoints);
02141
02142 }
02143
02144 if(DEBUG || kVerbose > 0) {
02145 G4cout << " first edge : (" << detector[0][0] << ", "
02146 << detector[0][1] << ", "
02147 << detector[0][2] << ") - ("
02148 << detector[0][3] << ", "
02149 << detector[0][4] << ", "
02150 << detector[0][5] << ")" << G4endl;
02151 }
02152
02153
02154 unsigned char dcolor[3];
02155 ifile.read((char *)dcolor, 3);
02156 if(DEBUG || kVerbose > 0) {
02157 G4cout << " detector color : rgb("
02158 << (int)dcolor[0] << ", "
02159 << (int)dcolor[1] << ", "
02160 << (int)dcolor[2] << G4endl;
02161 }
02162
02163
02164
02165 char cname[80];
02166 ifile.read((char *)cname, 80);
02167 std::string dname = cname;
02168 if(DEBUG || kVerbose > 0) {
02169 G4cout << " detector name : " << dname << G4endl;
02170 }
02171
02172
02173 addDetector(dname, detector, dcolor);
02174
02175 }
02176 }
02177
02178
02179 ifile.close();
02180
02181 return true;
02182 }
02183 bool G4GMocrenIO::retrieveData4(char * _filename) {
02184 kFileName = _filename;
02185 return retrieveData();
02186 }
02187
02188
02189 bool G4GMocrenIO::retrieveData3() {
02190
02191 bool DEBUG = false;
02192
02193
02194 std::ifstream ifile(kFileName.c_str(), std::ios_base::in|std::ios_base::binary);
02195 if(!ifile) {
02196 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
02197 G4cout << "Cannot open file: " << kFileName
02198 << " in G4GMocrenIO::retrieveData3()." << G4endl;
02199 return false;
02200 }
02201
02202
02203 char ctmp[12];
02204
02205
02206 char verid[9];
02207 ifile.read((char *)verid, 8);
02208
02209
02210 unsigned char ver;
02211 ifile.read((char *)&ver, 1);
02212 std::stringstream ss;
02213 ss << (int)ver;
02214 kVersion = ss.str();
02215 if(DEBUG || kVerbose > 0) G4cout << "File version : " << kVersion << G4endl;
02216
02217
02218 ifile.read((char *)&kLittleEndianInput, sizeof(char));
02219 if(DEBUG || kVerbose > 0) {
02220 G4cout << "Endian : ";
02221 if(kLittleEndianInput == 1)
02222 G4cout << " little" << G4endl;
02223 else {
02224 G4cout << " big" << G4endl;
02225 }
02226 }
02227
02228
02229 int clength;
02230 ifile.read((char *)ctmp, 4);
02231 convertEndian(ctmp, clength);
02232
02233 char cmt[1025];
02234 ifile.read((char *)cmt, clength);
02235 std::string scmt = cmt;
02236 setComment(scmt);
02237 if(DEBUG || kVerbose > 0) {
02238 G4cout << "Data comment : "
02239 << kComment << G4endl;
02240 }
02241
02242
02243 ifile.read((char *)ctmp, 12);
02244 convertEndian(ctmp, kVoxelSpacing[0]);
02245 convertEndian(ctmp+4, kVoxelSpacing[1]);
02246 convertEndian(ctmp+8, kVoxelSpacing[2]);
02247 if(DEBUG || kVerbose > 0) {
02248 G4cout << "Voxel spacing : ("
02249 << kVoxelSpacing[0] << ", "
02250 << kVoxelSpacing[1] << ", "
02251 << kVoxelSpacing[2]
02252 << ") mm " << G4endl;
02253 }
02254
02255
02256
02257 ifile.read((char *)ctmp, 4);
02258 convertEndian(ctmp, kPointerToModalityData);
02259
02260
02261 ifile.read((char *)ctmp, 4);
02262 int nDoseDist;
02263 convertEndian(ctmp, nDoseDist);
02264
02265
02266 for(int i = 0; i < nDoseDist; i++) {
02267 ifile.read((char *)ctmp, 4);
02268 unsigned int dptr;
02269 convertEndian(ctmp, dptr);
02270 addPointerToDoseDistData(dptr);
02271 }
02272
02273
02274 ifile.read((char *)ctmp, 4);
02275 convertEndian(ctmp, kPointerToROIData);
02276
02277
02278 ifile.read((char *)ctmp, 4);
02279 convertEndian(ctmp, kPointerToTrackData);
02280 if(DEBUG || kVerbose > 0) {
02281 G4cout << "Each pointer to data : "
02282 << kPointerToModalityData << ", ";
02283 for(int i = 0; i < nDoseDist; i++)
02284 G4cout << kPointerToDoseDistData[0] << ", ";
02285 G4cout << kPointerToROIData << ", "
02286 << kPointerToTrackData << G4endl;
02287 }
02288
02289 if(kPointerToModalityData == 0 && kPointerToDoseDistData.size() == 0 &&
02290 kPointerToROIData == 0 && kPointerToTrackData == 0) {
02291 if(DEBUG || kVerbose > 0) {
02292 G4cout << "No data." << G4endl;
02293 }
02294 return false;
02295 }
02296
02297
02298
02299
02300
02301
02302
02303 int size[3];
02304 float scale;
02305 double dscale;
02306 short minmax[2];
02307 float fCenter[3];
02308 int iCenter[3];
02309
02310
02311
02312 ifile.read(ctmp, 3*sizeof(int));
02313 convertEndian(ctmp, size[0]);
02314 convertEndian(ctmp+sizeof(int), size[1]);
02315 convertEndian(ctmp+2*sizeof(int), size[2]);
02316 if(DEBUG || kVerbose > 0) {
02317 G4cout << "Modality image size : ("
02318 << size[0] << ", "
02319 << size[1] << ", "
02320 << size[2] << ")"
02321 << G4endl;
02322 }
02323 kModality.setSize(size);
02324
02325
02326
02327
02328
02329
02330
02331
02332
02333 if(kPointerToModalityData != 0) {
02334
02335
02336 ifile.read((char *)ctmp, 4);
02337 convertEndian(ctmp, minmax[0]);
02338 convertEndian(ctmp+2, minmax[1]);
02339 kModality.setMinMax(minmax);
02340
02341
02342 char munit[13];
02343 ifile.read((char *)munit, 12);
02344 std::string smunit = munit;
02345 setModalityImageUnit(smunit);
02346
02347
02348 ifile.read((char *)ctmp, 4);
02349 convertEndian(ctmp, scale);
02350 kModality.setScale(dscale = scale);
02351 if(DEBUG || kVerbose > 0) {
02352 G4cout << "Modality image min., max., scale : "
02353 << minmax[0] << ", "
02354 << minmax[1] << ", "
02355 << scale << G4endl;
02356 }
02357
02358
02359 int psize = size[0]*size[1];
02360 if(DEBUG || kVerbose > 0) G4cout << "Modality image (" << psize << "): ";
02361 char * cimage = new char[psize*sizeof(short)];
02362 for(int i = 0; i < size[2]; i++) {
02363 ifile.read((char *)cimage, psize*sizeof(short));
02364 short * mimage = new short[psize];
02365 for(int j = 0; j < psize; j++) {
02366 convertEndian(cimage+j*sizeof(short), mimage[j]);
02367 }
02368 kModality.addImage(mimage);
02369
02370 if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << mimage[(size_t)(psize*0.55)] << ", ";
02371 }
02372 if(DEBUG || kVerbose > 0) G4cout << G4endl;
02373 delete [] cimage;
02374
02375
02376 size_t msize = minmax[1]-minmax[0]+1;
02377 if(DEBUG || kVerbose > 0) G4cout << "msize: " << msize << G4endl;
02378 char * pdmap = new char[msize*sizeof(float)];
02379 ifile.read((char *)pdmap, msize*sizeof(float));
02380 float ftmp;
02381 for(int i = 0; i < (int)msize; i++) {
02382 convertEndian(pdmap+i*sizeof(float), ftmp);
02383 kModalityImageDensityMap.push_back(ftmp);
02384 }
02385 delete [] pdmap;
02386 if(DEBUG || kVerbose > 0) {
02387 G4cout << "density map : " << std::ends;
02388 for(int i = 0; i < 10; i++)
02389 G4cout <<kModalityImageDensityMap[i] << ", ";
02390 G4cout << G4endl;
02391 for(int i = 0; i < 10; i++) G4cout << "..";
02392 G4cout << G4endl;
02393 for(size_t i =kModalityImageDensityMap.size() - 10; i <kModalityImageDensityMap.size(); i++)
02394 G4cout <<kModalityImageDensityMap[i] << ", ";
02395 G4cout << G4endl;
02396 }
02397
02398 }
02399
02400
02401
02402 for(int ndose = 0; ndose < nDoseDist; ndose++) {
02403
02404 newDoseDist();
02405
02406
02407 ifile.read((char *)ctmp, 3*sizeof(int));
02408 convertEndian(ctmp, size[0]);
02409 convertEndian(ctmp+sizeof(int), size[1]);
02410 convertEndian(ctmp+2*sizeof(int), size[2]);
02411 if(DEBUG || kVerbose > 0) {
02412 G4cout << "Dose dist. image size : ("
02413 << size[0] << ", "
02414 << size[1] << ", "
02415 << size[2] << ")"
02416 << G4endl;
02417 }
02418 kDose[ndose].setSize(size);
02419
02420
02421 ifile.read((char *)ctmp, sizeof(short)*2);
02422 convertEndian(ctmp, minmax[0]);
02423 convertEndian(ctmp+2, minmax[1]);
02424
02425
02426 char dunit[13];
02427 ifile.read((char *)dunit, 12);
02428 std::string sdunit = dunit;
02429 setDoseDistUnit(sdunit, ndose);
02430 if(DEBUG || kVerbose > 0) {
02431 G4cout << "Dose dist. unit : " << kDoseUnit << G4endl;
02432 }
02433
02434
02435 ifile.read((char *)ctmp, 4);
02436 convertEndian(ctmp, scale);
02437 kDose[ndose].setScale(dscale = scale);
02438
02439 double dminmax[2];
02440 for(int i = 0; i < 2; i++) dminmax[i] = minmax[i]*dscale;
02441 kDose[ndose].setMinMax(dminmax);
02442
02443 if(DEBUG || kVerbose > 0) {
02444 G4cout << "Dose dist. image min., max., scale : "
02445 << dminmax[0] << ", "
02446 << dminmax[1] << ", "
02447 << scale << G4endl;
02448 }
02449
02450
02451 int dsize = size[0]*size[1];
02452 if(DEBUG || kVerbose > 0) G4cout << "Dose dist. (" << dsize << "): ";
02453 char * di = new char[dsize*sizeof(short)];
02454 short * shimage = new short[dsize];
02455 for(int z = 0; z < size[2]; z++) {
02456 ifile.read((char *)di, dsize*sizeof(short));
02457 double * dimage = new double[dsize];
02458 for(int xy = 0; xy < dsize; xy++) {
02459 convertEndian(di+xy*sizeof(short), shimage[xy]);
02460 dimage[xy] = shimage[xy]*dscale;
02461 }
02462 kDose[ndose].addImage(dimage);
02463
02464 if(DEBUG || kVerbose > 0) G4cout << "[" << z << "]" << dimage[(size_t)(dsize*0.55)] << ", ";
02465
02466 if(DEBUG || kVerbose > 0) {
02467 for(int j = 0; j < dsize; j++) {
02468 if(dimage[j] < 0)
02469 G4cout << "[" << j << "," << z << "]"
02470 << dimage[j] << ", ";
02471 }
02472 }
02473 }
02474 delete [] shimage;
02475 delete [] di;
02476 if(DEBUG || kVerbose > 0) G4cout << G4endl;
02477
02478 ifile.read((char *)ctmp, 3*4);
02479 convertEndian(ctmp, iCenter[0]);
02480 convertEndian(ctmp+4, iCenter[1]);
02481 convertEndian(ctmp+8, iCenter[2]);
02482 for(int i = 0; i < 3; i++) fCenter[i] = (float)iCenter[i];
02483 kDose[ndose].setCenterPosition(fCenter);
02484
02485 if(DEBUG || kVerbose > 0) {
02486 G4cout << "Dose dist. image relative location : ("
02487 << fCenter[0] << ", "
02488 << fCenter[1] << ", "
02489 << fCenter[2] << ")" << G4endl;
02490 }
02491
02492
02493 }
02494
02495
02496 if(kPointerToROIData != 0) {
02497
02498 newROI();
02499
02500
02501 ifile.read((char *)ctmp, 3*sizeof(int));
02502 convertEndian(ctmp, size[0]);
02503 convertEndian(ctmp+sizeof(int), size[1]);
02504 convertEndian(ctmp+2*sizeof(int), size[2]);
02505 kRoi[0].setSize(size);
02506 if(DEBUG || kVerbose > 0) {
02507 G4cout << "ROI image size : ("
02508 << size[0] << ", "
02509 << size[1] << ", "
02510 << size[2] << ")"
02511 << G4endl;
02512 }
02513
02514
02515 ifile.read((char *)ctmp, sizeof(short)*2);
02516 convertEndian(ctmp, minmax[0]);
02517 convertEndian(ctmp+sizeof(short), minmax[1]);
02518 kRoi[0].setMinMax(minmax);
02519
02520
02521 ifile.read((char *)ctmp, sizeof(float));
02522 convertEndian(ctmp, scale);
02523 kRoi[0].setScale(dscale = scale);
02524 if(DEBUG || kVerbose > 0) {
02525 G4cout << "ROI image min., max., scale : "
02526 << minmax[0] << ", "
02527 << minmax[1] << ", "
02528 << scale << G4endl;
02529 }
02530
02531
02532 int rsize = size[0]*size[1];
02533 char * ri = new char[rsize*sizeof(short)];
02534 for(int i = 0; i < size[2]; i++) {
02535 ifile.read((char *)ri, rsize*sizeof(short));
02536 short * rimage = new short[rsize];
02537 for(int j = 0; j < rsize; j++) {
02538 convertEndian(ri+j*sizeof(short), rimage[j]);
02539 }
02540 kRoi[0].addImage(rimage);
02541
02542 }
02543 delete [] ri;
02544
02545
02546 ifile.read((char *)ctmp, 3*sizeof(int));
02547 convertEndian(ctmp, iCenter[0]);
02548 convertEndian(ctmp+sizeof(int), iCenter[1]);
02549 convertEndian(ctmp+2*sizeof(int), iCenter[2]);
02550 for(int i = 0; i < 3; i++) fCenter[i] = iCenter[i];
02551 kRoi[0].setCenterPosition(fCenter);
02552 if(DEBUG || kVerbose > 0) {
02553 G4cout << "ROI image relative location : ("
02554 << fCenter[0] << ", "
02555 << fCenter[1] << ", "
02556 << fCenter[2] << ")" << G4endl;
02557 }
02558
02559 }
02560
02561
02562 if(kPointerToTrackData != 0) {
02563
02564
02565 ifile.read((char *)ctmp, sizeof(int));
02566 int ntrk;
02567 convertEndian(ctmp, ntrk);
02568 if(DEBUG || kVerbose > 0) {
02569 G4cout << "# of tracks: " << ntrk << G4endl;
02570 }
02571
02572
02573 std::vector<float *> trkv4;
02574
02575
02576 for(int i = 0; i < ntrk; i++) {
02577 float * tp = new float[6];
02578
02579 ifile.read((char *)ctmp, sizeof(float)*3);
02580 if(DEBUG || kVerbose > 0) if(i < 10) G4cout << i << ": " ;
02581 for(int j = 0; j < 3; j++) {
02582 convertEndian(ctmp+j*sizeof(float), tp[j]);
02583 if(DEBUG || kVerbose > 0) if(i < 10) G4cout << tp[j] << ", ";
02584 }
02585
02586 ifile.read((char *)ctmp, sizeof(float)*3);
02587 for(int j = 0; j < 3; j++) {
02588 convertEndian(ctmp+j*sizeof(float), tp[j+3]);
02589 if(DEBUG || kVerbose > 0) if(i < 10) G4cout << tp[j+3] << ", ";
02590 }
02591 addTrack(tp);
02592 if(DEBUG || kVerbose > 0) if(i < 10) G4cout << G4endl;
02593
02594
02595 trkv4.push_back(tp);
02596 }
02597
02598
02599 unsigned char trkcolorv4[3];
02600
02601
02602 for(int i = 0; i < ntrk; i++) {
02603 unsigned char * rgb = new unsigned char[3];
02604 ifile.read((char *)rgb, 3);
02605 addTrackColor(rgb);
02606
02607
02608 for(int j = 0; j < 3; j++) trkcolorv4[j] = rgb[j];
02609 std::vector<float *> trk;
02610 trk.push_back(trkv4[i]);
02611 addTrack(trk, trkcolorv4);
02612
02613 }
02614
02615 }
02616
02617 ifile.close();
02618
02619 return true;
02620 }
02621 bool G4GMocrenIO::retrieveData3(char * _filename) {
02622 kFileName = _filename;
02623 return retrieveData();
02624 }
02625
02626
02627 bool G4GMocrenIO::retrieveData2() {
02628
02629 bool DEBUG = false;
02630
02631
02632 std::ifstream ifile(kFileName.c_str(), std::ios_base::in|std::ios_base::binary);
02633 if(!ifile) {
02634 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
02635 G4cout << "Cannot open file: " << kFileName
02636 << " in G4GMocrenIO::retrieveData2()." << G4endl;
02637 return false;
02638 }
02639
02640
02641 char ctmp[12];
02642
02643
02644 char verid[9];
02645 ifile.read((char *)verid, 8);
02646
02647
02648 unsigned char ver;
02649 ifile.read((char *)&ver, 1);
02650 std::stringstream ss;
02651 ss << (int)ver;
02652 kVersion = ss.str();
02653 if(DEBUG || kVerbose > 0) G4cout << "File version : " << kVersion << G4endl;
02654
02655
02656 char idtmp[IDLENGTH];
02657 ifile.read((char *)idtmp, IDLENGTH);
02658 kId = idtmp;
02659
02660 char vertmp[VERLENGTH];
02661 ifile.read((char *)vertmp, VERLENGTH);
02662
02663
02664 ifile.read((char *)&kLittleEndianInput, sizeof(char));
02665 if(DEBUG || kVerbose > 0) {
02666 G4cout << "Endian : ";
02667 if(kLittleEndianInput == 1)
02668 G4cout << " little" << G4endl;
02669 else {
02670 G4cout << " big" << G4endl;
02671 }
02672 }
02673
02674
02675 ifile.read((char *)ctmp, 12);
02676 convertEndian(ctmp, kVoxelSpacing[0]);
02677 convertEndian(ctmp+4, kVoxelSpacing[1]);
02678 convertEndian(ctmp+8, kVoxelSpacing[2]);
02679 if(DEBUG || kVerbose > 0) {
02680 G4cout << "Voxel spacing : ("
02681 << kVoxelSpacing[0] << ", "
02682 << kVoxelSpacing[1] << ", "
02683 << kVoxelSpacing[2]
02684 << ") mm " << G4endl;
02685 }
02686
02687
02688
02689 ifile.read((char *)ctmp, 4);
02690 convertEndian(ctmp, kPointerToModalityData);
02691
02692
02693 unsigned int ptddd;
02694 ifile.read((char *)ctmp, 4);
02695 convertEndian(ctmp, ptddd);
02696 kPointerToDoseDistData.push_back(ptddd);
02697
02698
02699 ifile.read((char *)ctmp, 4);
02700 convertEndian(ctmp, kPointerToROIData);
02701
02702
02703 ifile.read((char *)ctmp, 4);
02704 convertEndian(ctmp, kPointerToTrackData);
02705 if(DEBUG || kVerbose > 0) {
02706 G4cout << "Each pointer to data : "
02707 << kPointerToModalityData << ", "
02708 << kPointerToDoseDistData[0] << ", "
02709 << kPointerToROIData << ", "
02710 << kPointerToTrackData << G4endl;
02711 }
02712
02713 if(kPointerToModalityData == 0 && kPointerToDoseDistData.size() == 0 &&
02714 kPointerToROIData == 0 && kPointerToTrackData == 0) {
02715 if(DEBUG || kVerbose > 0) {
02716 G4cout << "No data." << G4endl;
02717 }
02718 return false;
02719 }
02720
02721
02722
02723
02724
02725
02726
02727 int size[3];
02728 float scale;
02729 double dscale;
02730 short minmax[2];
02731 float fCenter[3];
02732 int iCenter[3];
02733
02734
02735
02736 ifile.read(ctmp, 3*sizeof(int));
02737 convertEndian(ctmp, size[0]);
02738 convertEndian(ctmp+sizeof(int), size[1]);
02739 convertEndian(ctmp+2*sizeof(int), size[2]);
02740 if(DEBUG || kVerbose > 0) {
02741 G4cout << "Modality image size : ("
02742 << size[0] << ", "
02743 << size[1] << ", "
02744 << size[2] << ")"
02745 << G4endl;
02746 }
02747 kModality.setSize(size);
02748
02749
02750
02751
02752
02753
02754
02755
02756
02757 if(kPointerToModalityData != 0) {
02758
02759
02760 ifile.read((char *)ctmp, 4);
02761 convertEndian(ctmp, minmax[0]);
02762 convertEndian(ctmp+2, minmax[1]);
02763 kModality.setMinMax(minmax);
02764
02765
02766 ifile.read((char *)ctmp, 4);
02767 convertEndian(ctmp, scale);
02768 kModality.setScale(dscale = scale);
02769 if(DEBUG || kVerbose > 0) {
02770 G4cout << "Modality image min., max., scale : "
02771 << minmax[0] << ", "
02772 << minmax[1] << ", "
02773 << scale << G4endl;
02774 }
02775
02776
02777 int psize = size[0]*size[1];
02778 if(DEBUG || kVerbose > 0) G4cout << "Modality image (" << psize << "): ";
02779 char * cimage = new char[psize*sizeof(short)];
02780 for(int i = 0; i < size[2]; i++) {
02781 ifile.read((char *)cimage, psize*sizeof(short));
02782 short * mimage = new short[psize];
02783 for(int j = 0; j < psize; j++) {
02784 convertEndian(cimage+j*sizeof(short), mimage[j]);
02785 }
02786 kModality.addImage(mimage);
02787
02788 if(DEBUG || kVerbose > 0) G4cout << "[" << i << "]" << mimage[(size_t)(psize*0.55)] << ", ";
02789 }
02790 if(DEBUG || kVerbose > 0) G4cout << G4endl;
02791 delete [] cimage;
02792
02793
02794 size_t msize = minmax[1]-minmax[0]+1;
02795 if(DEBUG || kVerbose > 0) G4cout << "msize: " << msize << G4endl;
02796 char * pdmap = new char[msize*sizeof(float)];
02797 ifile.read((char *)pdmap, msize*sizeof(float));
02798 float ftmp;
02799 for(int i = 0; i < (int)msize; i++) {
02800 convertEndian(pdmap+i*sizeof(float), ftmp);
02801 kModalityImageDensityMap.push_back(ftmp);
02802 }
02803 delete [] pdmap;
02804 if(DEBUG || kVerbose > 0) {
02805 G4cout << "density map : " << std::ends;
02806 for(int i = 0; i < 10; i++)
02807 G4cout <<kModalityImageDensityMap[i] << ", ";
02808 G4cout << G4endl;
02809 for(int i = 0; i < 10; i++) G4cout << "..";
02810 G4cout << G4endl;
02811 for(size_t i =kModalityImageDensityMap.size() - 10; i <kModalityImageDensityMap.size(); i++)
02812 G4cout <<kModalityImageDensityMap[i] << ", ";
02813 G4cout << G4endl;
02814 }
02815
02816 }
02817
02818
02819
02820 if(kPointerToDoseDistData[0] != 0) {
02821
02822 newDoseDist();
02823
02824
02825 ifile.read((char *)ctmp, 3*sizeof(int));
02826 convertEndian(ctmp, size[0]);
02827 convertEndian(ctmp+sizeof(int), size[1]);
02828 convertEndian(ctmp+2*sizeof(int), size[2]);
02829 if(DEBUG || kVerbose > 0) {
02830 G4cout << "Dose dist. image size : ("
02831 << size[0] << ", "
02832 << size[1] << ", "
02833 << size[2] << ")"
02834 << G4endl;
02835 }
02836 kDose[0].setSize(size);
02837
02838
02839 ifile.read((char *)ctmp, sizeof(short)*2);
02840 convertEndian(ctmp, minmax[0]);
02841 convertEndian(ctmp+2, minmax[1]);
02842
02843 ifile.read((char *)ctmp, sizeof(float));
02844 convertEndian(ctmp, scale);
02845 kDose[0].setScale(dscale = scale);
02846
02847 double dminmax[2];
02848 for(int i = 0; i < 2; i++) dminmax[i] = minmax[i]*dscale;
02849 kDose[0].setMinMax(dminmax);
02850
02851 if(DEBUG || kVerbose > 0) {
02852 G4cout << "Dose dist. image min., max., scale : "
02853 << dminmax[0] << ", "
02854 << dminmax[1] << ", "
02855 << scale << G4endl;
02856 }
02857
02858
02859 int dsize = size[0]*size[1];
02860 if(DEBUG || kVerbose > 0) G4cout << "Dose dist. (" << dsize << "): ";
02861 char * di = new char[dsize*sizeof(short)];
02862 short * shimage = new short[dsize];
02863 for(int z = 0; z < size[2]; z++) {
02864 ifile.read((char *)di, dsize*sizeof(short));
02865 double * dimage = new double[dsize];
02866 for(int xy = 0; xy < dsize; xy++) {
02867 convertEndian(di+xy*sizeof(short), shimage[xy]);
02868 dimage[xy] = shimage[xy]*dscale;
02869 }
02870 kDose[0].addImage(dimage);
02871
02872 if(DEBUG || kVerbose > 0) G4cout << "[" << z << "]" << dimage[(size_t)(dsize*0.55)] << ", ";
02873
02874 if(DEBUG || kVerbose > 0) {
02875 for(int j = 0; j < dsize; j++) {
02876 if(dimage[j] < 0)
02877 G4cout << "[" << j << "," << z << "]"
02878 << dimage[j] << ", ";
02879 }
02880 }
02881 }
02882 delete [] shimage;
02883 delete [] di;
02884 if(DEBUG || kVerbose > 0) G4cout << G4endl;
02885
02886
02887
02888
02889
02890
02891
02892
02893
02894
02895
02896
02897
02898
02899
02900
02901
02902
02903
02904
02905
02906
02907
02908
02909
02910
02911
02912
02913
02914
02915
02916
02917
02918
02919
02920 ifile.read((char *)ctmp, 3*sizeof(int));
02921 convertEndian(ctmp, iCenter[0]);
02922 convertEndian(ctmp+sizeof(int), iCenter[1]);
02923 convertEndian(ctmp+2*sizeof(int), iCenter[2]);
02924 for(int i = 0; i < 3; i++) fCenter[i] = (float)iCenter[i];
02925 kDose[0].setCenterPosition(fCenter);
02926
02927 if(DEBUG || kVerbose > 0) {
02928 G4cout << "Dose dist. image relative location : ("
02929 << fCenter[0] << ", "
02930 << fCenter[1] << ", "
02931 << fCenter[2] << ")" << G4endl;
02932 }
02933
02934
02935 }
02936
02937
02938 if(kPointerToROIData != 0) {
02939
02940 newROI();
02941
02942
02943 ifile.read((char *)ctmp, 3*sizeof(int));
02944 convertEndian(ctmp, size[0]);
02945 convertEndian(ctmp+sizeof(int), size[1]);
02946 convertEndian(ctmp+2*sizeof(int), size[2]);
02947 kRoi[0].setSize(size);
02948 if(DEBUG || kVerbose > 0) {
02949 G4cout << "ROI image size : ("
02950 << size[0] << ", "
02951 << size[1] << ", "
02952 << size[2] << ")"
02953 << G4endl;
02954 }
02955
02956
02957 ifile.read((char *)ctmp, sizeof(short)*2);
02958 convertEndian(ctmp, minmax[0]);
02959 convertEndian(ctmp+sizeof(short), minmax[1]);
02960 kRoi[0].setMinMax(minmax);
02961
02962
02963 ifile.read((char *)ctmp, sizeof(float));
02964 convertEndian(ctmp, scale);
02965 kRoi[0].setScale(dscale = scale);
02966 if(DEBUG || kVerbose > 0) {
02967 G4cout << "ROI image min., max., scale : "
02968 << minmax[0] << ", "
02969 << minmax[1] << ", "
02970 << scale << G4endl;
02971 }
02972
02973
02974 int rsize = size[0]*size[1];
02975 char * ri = new char[rsize*sizeof(short)];
02976 for(int i = 0; i < size[2]; i++) {
02977 ifile.read((char *)ri, rsize*sizeof(short));
02978 short * rimage = new short[rsize];
02979 for(int j = 0; j < rsize; j++) {
02980 convertEndian(ri+j*sizeof(short), rimage[j]);
02981 }
02982 kRoi[0].addImage(rimage);
02983
02984 }
02985 delete [] ri;
02986
02987
02988 ifile.read((char *)ctmp, 3*sizeof(int));
02989 convertEndian(ctmp, iCenter[0]);
02990 convertEndian(ctmp+sizeof(int), iCenter[1]);
02991 convertEndian(ctmp+2*sizeof(int), iCenter[2]);
02992 for(int i = 0; i < 3; i++) fCenter[i] = iCenter[i];
02993 kRoi[0].setCenterPosition(fCenter);
02994 if(DEBUG || kVerbose > 0) {
02995 G4cout << "ROI image relative location : ("
02996 << fCenter[0] << ", "
02997 << fCenter[1] << ", "
02998 << fCenter[2] << ")" << G4endl;
02999 }
03000
03001 }
03002
03003
03004 if(kPointerToTrackData != 0) {
03005
03006
03007 ifile.read((char *)ctmp, sizeof(int));
03008 int ntrk;
03009 convertEndian(ctmp, ntrk);
03010 if(DEBUG || kVerbose > 0) {
03011 G4cout << "# of tracks: " << ntrk << G4endl;
03012 }
03013
03014
03015 unsigned char trkcolorv4[3] = {255, 0, 0};
03016
03017 for(int i = 0; i < ntrk; i++) {
03018 float * tp = new float[6];
03019
03020 std::vector<float *> trkv4;
03021
03022 ifile.read((char *)ctmp, sizeof(float)*3);
03023 if(DEBUG || kVerbose > 0) if(i < 10) G4cout << i << ": " ;
03024 for(int j = 0; j < 3; j++) {
03025 convertEndian(ctmp+j*sizeof(float), tp[j]);
03026 if(DEBUG || kVerbose > 0) if(i < 10) G4cout << tp[j] << ", ";
03027 }
03028
03029 ifile.read((char *)ctmp, sizeof(float)*3);
03030 for(int j = 0; j < 3; j++) {
03031 convertEndian(ctmp+j*sizeof(float), tp[j+3]);
03032 if(DEBUG || kVerbose > 0) if(i < 10) G4cout << tp[j+3] << ", ";
03033 }
03034
03035 kSteps.push_back(tp);
03036
03037 trkv4.push_back(tp);
03038 addTrack(trkv4, trkcolorv4);
03039
03040 if(DEBUG || kVerbose > 0) if(i < 10) G4cout << G4endl;
03041 }
03042
03043 }
03044
03045
03046
03047
03048
03049
03050
03051
03052
03053
03054
03055
03056
03057
03058
03059
03060
03061
03062
03063
03064
03065
03066
03067
03068
03069
03070
03071
03072
03073
03074
03075
03076
03077
03078
03079
03080
03081
03082
03083
03084
03085
03086
03087 ifile.close();
03088
03089 return true;
03090 }
03091
03092 bool G4GMocrenIO::retrieveData2(char * _filename) {
03093 kFileName = _filename;
03094 return retrieveData();
03095 }
03096
03097 void G4GMocrenIO::setID() {
03098 time_t t;
03099 time(&t);
03100
03101 tm * ti;
03102 ti = localtime(&t);
03103
03104 char cmonth[12][4] = {"Jan", "Feb", "Mar", "Apr",
03105 "May", "Jun", "Jul", "Aug",
03106 "Sep", "Oct", "Nov", "Dec"};
03107 std::stringstream ss;
03108 ss << std::setfill('0')
03109 << std::setw(2)
03110 << ti->tm_hour << ":"
03111 << std::setw(2)
03112 << ti->tm_min << ":"
03113 << std::setw(2)
03114 << ti->tm_sec << ","
03115 << cmonth[ti->tm_mon] << "."
03116 << std::setw(2)
03117 << ti->tm_mday << ","
03118 << ti->tm_year+1900;
03119
03120 kId = ss.str();
03121 }
03122
03123
03124 std::string & G4GMocrenIO::getVersion() {return kVersion;}
03125 void G4GMocrenIO::setVersion(std::string & _version) {kVersion = _version;}
03126
03127
03128 void G4GMocrenIO::setLittleEndianInput(bool _little) {kLittleEndianInput = _little;}
03129 void G4GMocrenIO::setLittleEndianOutput(bool _little) {kLittleEndianOutput = _little;}
03130
03131
03132 void G4GMocrenIO::setVoxelSpacing(float _spacing[3]) {
03133 for(int i = 0; i < 3; i++) kVoxelSpacing[i] = _spacing[i];
03134 }
03135 void G4GMocrenIO::getVoxelSpacing(float _spacing[3]) {
03136 for(int i = 0; i < 3; i++) _spacing[i] = kVoxelSpacing[i];
03137 }
03138
03139
03140 int & G4GMocrenIO::getNumberOfEvents() {
03141 return kNumberOfEvents;
03142 }
03143 void G4GMocrenIO::setNumberOfEvents(int & _numberOfEvents) {
03144 kNumberOfEvents = _numberOfEvents;
03145 }
03146 void G4GMocrenIO::addOneEvent() {
03147 kNumberOfEvents++;
03148 }
03149
03150
03151 void G4GMocrenIO::setPointerToModalityData(unsigned int & _pointer) {
03152 kPointerToModalityData = _pointer;
03153 }
03154 unsigned int G4GMocrenIO::getPointerToModalityData() {
03155 return kPointerToModalityData;
03156 }
03157
03158 void G4GMocrenIO::addPointerToDoseDistData(unsigned int & _pointer) {
03159 kPointerToDoseDistData.push_back(_pointer);
03160 }
03161 unsigned int G4GMocrenIO::getPointerToDoseDistData(int _elem) {
03162 if(kPointerToDoseDistData.size() == 0 ||
03163 kPointerToDoseDistData.size() < (size_t)_elem)
03164 return 0;
03165 else
03166 return kPointerToDoseDistData[_elem];
03167 }
03168
03169
03170 void G4GMocrenIO::setPointerToROIData(unsigned int & _pointer) {
03171 kPointerToROIData = _pointer;
03172 }
03173 unsigned int G4GMocrenIO::getPointerToROIData() {
03174 return kPointerToROIData;
03175 }
03176
03177 void G4GMocrenIO::setPointerToTrackData(unsigned int & _pointer) {
03178 kPointerToTrackData = _pointer;
03179 }
03180 unsigned int G4GMocrenIO::getPointerToTrackData() {
03181 return kPointerToTrackData;
03182 }
03183
03184
03185 void G4GMocrenIO::calcPointers4() {
03186
03187
03188 unsigned int pointer = 1070;
03189 int nDoseDist = getNumDoseDist();
03190 pointer += nDoseDist*4;
03191
03192 setPointerToModalityData(pointer);
03193
03194
03195
03196 int msize[3];
03197 getModalityImageSize(msize);
03198 short mminmax[2];
03199 getModalityImageMinMax(mminmax);
03200 int pmsize = 2*msize[0]*msize[1]*msize[2];
03201 int pmmap = 4*(mminmax[1] - mminmax[0] + 1);
03202 pointer += 32 + pmsize + pmmap;
03203
03204 kPointerToDoseDistData.clear();
03205 if(nDoseDist == 0) {
03206 unsigned int pointer0 = 0;
03207 addPointerToDoseDistData(pointer0);
03208 }
03209 for(int ndose = 0; ndose < nDoseDist; ndose++) {
03210 addPointerToDoseDistData(pointer);
03211 int dsize[3];
03212 getDoseDistSize(dsize);
03213 pointer += 44 + dsize[0]*dsize[1]*dsize[2]*2 + 80;
03214 }
03215
03216
03217 if(!isROIEmpty()) {
03218 setPointerToROIData(pointer);
03219
03220 int rsize[3];
03221 getROISize(rsize);
03222 int prsize = 2*rsize[0]*rsize[1]*rsize[2];
03223 pointer += 20 + prsize + 12;
03224 } else {
03225 unsigned int pointer0 = 0;
03226 setPointerToROIData(pointer0);
03227 }
03228
03229
03230 int ntrk = kTracks.size();
03231 if(ntrk != 0) {
03232 setPointerToTrackData(pointer);
03233
03234 pointer += 4;
03235 for(int nt = 0; nt < ntrk; nt++) {
03236 int nsteps = kTracks[nt].getNumberOfSteps();
03237 pointer += 4 + 3 + nsteps*(4*6);
03238 }
03239 } else {
03240 unsigned int pointer0 = 0;
03241 setPointerToTrackData(pointer0);
03242 }
03243 if(kVerbose > 0) G4cout << " pointer to the track data :"
03244 << kPointerToTrackData << G4endl;
03245
03246
03247 int ndet = kDetectors.size();
03248 if(ndet != 0) {
03249 kPointerToDetectorData = pointer;
03250 } else {
03251 kPointerToDetectorData = 0;
03252 }
03253 if(kVerbose > 0) G4cout << " pointer to the detector data :"
03254 << kPointerToDetectorData << G4endl;
03255
03256 }
03257
03258
03259 void G4GMocrenIO::calcPointers3() {
03260
03261
03262 unsigned int pointer = 1066;
03263 int nDoseDist = getNumDoseDist();
03264 pointer += nDoseDist*4;
03265
03266 setPointerToModalityData(pointer);
03267
03268
03269
03270 int msize[3];
03271 getModalityImageSize(msize);
03272 short mminmax[2];
03273 getModalityImageMinMax(mminmax);
03274 int pmsize = 2*msize[0]*msize[1]*msize[2];
03275 int pmmap = 4*(mminmax[1] - mminmax[0] + 1);
03276 pointer += 32 + pmsize + pmmap;
03277
03278 kPointerToDoseDistData.clear();
03279 if(nDoseDist == 0) {
03280 unsigned int pointer0 = 0;
03281 addPointerToDoseDistData(pointer0);
03282 }
03283 for(int ndose = 0; ndose < nDoseDist; ndose++) {
03284 addPointerToDoseDistData(pointer);
03285 int dsize[3];
03286 getDoseDistSize(dsize);
03287 pointer += 44 + dsize[0]*dsize[1]*dsize[2]*2;
03288 }
03289
03290
03291 if(!isROIEmpty()) {
03292 setPointerToROIData(pointer);
03293
03294 int rsize[3];
03295 getROISize(rsize);
03296 int prsize = 2*rsize[0]*rsize[1]*rsize[2];
03297 pointer += 20 + prsize + 12;
03298 } else {
03299 unsigned int pointer0 = 0;
03300 setPointerToROIData(pointer0);
03301 }
03302
03303
03304 if(getNumTracks() != 0)
03305 setPointerToTrackData(pointer);
03306 else {
03307 unsigned int pointer0 = 0;
03308 setPointerToTrackData(pointer0);
03309 }
03310
03311 }
03312
03313
03314 void G4GMocrenIO::calcPointers2() {
03315
03316
03317 unsigned int pointer = 65;
03318 setPointerToModalityData(pointer);
03319
03320
03321 int msize[3];
03322 getModalityImageSize(msize);
03323 short mminmax[2];
03324 getModalityImageMinMax(mminmax);
03325 int pmsize = 2*msize[0]*msize[1]*msize[2];
03326 int pmmap = 4*(mminmax[1] - mminmax[0] + 1);
03327 pointer += 20 + pmsize + pmmap;
03328 int dsize[3];
03329 getDoseDistSize(dsize);
03330 kPointerToDoseDistData.clear();
03331 if(dsize[0] != 0) {
03332 kPointerToDoseDistData.push_back(pointer);
03333
03334 int pdsize = 2*dsize[0]*dsize[1]*dsize[2];
03335 pointer += 20 + pdsize + 12;
03336 } else {
03337 unsigned int pointer0 = 0;
03338 kPointerToDoseDistData.push_back(pointer0);
03339 }
03340
03341
03342 if(!isROIEmpty()) {
03343 int rsize[3];
03344 getROISize(rsize);
03345 setPointerToROIData(pointer);
03346 int prsize = 2*rsize[0]*rsize[1]*rsize[2];
03347 pointer += 20 + prsize + 12;
03348
03349 } else {
03350 unsigned int pointer0 = 0;
03351 setPointerToROIData(pointer0);
03352 }
03353
03354
03355 if(getNumTracks() != 0)
03356 setPointerToTrackData(pointer);
03357 else {
03358 unsigned int pointer0 = 0;
03359 setPointerToTrackData(pointer0);
03360 }
03361
03362 }
03363
03364
03365
03366 void G4GMocrenIO::getModalityImageSize(int _size[3]) {
03367
03368 kModality.getSize(_size);
03369 }
03370 void G4GMocrenIO::setModalityImageSize(int _size[3]) {
03371
03372 kModality.setSize(_size);
03373 }
03374
03375
03376 void G4GMocrenIO::setModalityImageScale(double & _scale) {
03377
03378 kModality.setScale(_scale);
03379 }
03380 double G4GMocrenIO::getModalityImageScale() {
03381
03382 return kModality.getScale();
03383 }
03384
03385
03386 void G4GMocrenIO::setModalityImage(short * _image) {
03387
03388 kModality.addImage(_image);
03389 }
03390 short * G4GMocrenIO::getModalityImage(int _z) {
03391
03392 return kModality.getImage(_z);
03393 }
03394 void G4GMocrenIO::clearModalityImage() {
03395
03396 kModality.clearImage();
03397 }
03398
03399 void G4GMocrenIO::setModalityImageDensityMap(std::vector<float> & _map) {
03400 kModalityImageDensityMap = _map;
03401 }
03402 std::vector<float> & G4GMocrenIO::getModalityImageDensityMap() {
03403 return kModalityImageDensityMap;
03404 }
03405
03406 void G4GMocrenIO::setModalityImageMinMax(short _minmax[2]) {
03407
03408 kModality.setMinMax(_minmax);
03409 }
03410
03411 void G4GMocrenIO::getModalityImageMinMax(short _minmax[2]) {
03412
03413 short minmax[2];
03414 kModality.getMinMax(minmax);
03415 for(int i = 0; i < 2; i++) _minmax[i] = minmax[i];
03416 }
03417 short G4GMocrenIO::getModalityImageMax() {
03418
03419 short minmax[2];
03420 kModality.getMinMax(minmax);
03421 return minmax[1];
03422 }
03423 short G4GMocrenIO::getModalityImageMin() {
03424
03425 short minmax[2];
03426 kModality.getMinMax(minmax);
03427 return minmax[0];
03428 }
03429
03430 void G4GMocrenIO::setModalityCenterPosition(float _center[3]) {
03431
03432 kModality.setCenterPosition(_center);
03433 }
03434 void G4GMocrenIO::getModalityCenterPosition(float _center[3]) {
03435
03436 if(isROIEmpty())
03437 for(int i = 0; i < 3; i++) _center[i] = 0;
03438 else
03439 kModality.getCenterPosition(_center);
03440 }
03441
03442 std::string G4GMocrenIO::getModalityImageUnit() {
03443 return kModalityUnit;
03444 }
03445 void G4GMocrenIO::setModalityImageUnit(std::string & _unit) {
03446 kModalityUnit = _unit;
03447 }
03448
03449 short G4GMocrenIO::convertDensityToHU(float & _dens) {
03450 short rval = -1024;
03451 int nmap = (int)kModalityImageDensityMap.size();
03452 if(nmap != 0) {
03453 short minmax[2];
03454 kModality.getMinMax(minmax);
03455 rval = minmax[1];
03456 for(int i = 0; i < nmap; i++) {
03457
03458 if(_dens <= kModalityImageDensityMap[i]) {
03459 rval = i + minmax[0];
03460 break;
03461 }
03462 }
03463 }
03464 return rval;
03465 }
03466
03467
03468
03469
03470 void G4GMocrenIO::newDoseDist() {
03471 GMocrenDataPrimitive<double> doseData;
03472 kDose.push_back(doseData);
03473 }
03474 int G4GMocrenIO::getNumDoseDist() {
03475 return (int)kDose.size();
03476 }
03477
03478
03479 std::string G4GMocrenIO::getDoseDistUnit(int _num) {
03480
03481 if(kDoseUnit.size() > static_cast<size_t>(_num)) return kDoseUnit;
03482
03483 return kDoseUnit;
03484 }
03485 void G4GMocrenIO::setDoseDistUnit(std::string & _unit, int _num) {
03486
03487 if(_unit.size() > static_cast<size_t>(_num)) kDoseUnit = _unit;
03488
03489
03490
03491
03492 kDoseUnit = _unit;
03493 }
03494
03495 void G4GMocrenIO::getDoseDistSize(int _size[3], int _num) {
03496 if(isDoseEmpty())
03497 for(int i = 0; i < 3; i++) _size[i] = 0;
03498 else
03499 kDose[_num].getSize(_size);
03500 }
03501 void G4GMocrenIO::setDoseDistSize(int _size[3], int _num) {
03502
03503 kDose[_num].setSize(_size);
03504
03505
03506 }
03507
03508 void G4GMocrenIO::setDoseDistMinMax(short _minmax[2], int _num) {
03509
03510 double minmax[2];
03511 double scale = kDose[_num].getScale();
03512 for(int i = 0; i < 2; i++) minmax[i] = (double)_minmax[i]*scale;
03513 kDose[_num].setMinMax(minmax);
03514 }
03515 void G4GMocrenIO::getDoseDistMinMax(short _minmax[2], int _num) {
03516
03517 if(isDoseEmpty())
03518 for(int i = 0; i < 2; i++) _minmax[i] = 0;
03519 else {
03520 double minmax[2];
03521 kDose[_num].getMinMax(minmax);
03522 double scale = kDose[_num].getScale();
03523 for(int i = 0; i < 2; i++) _minmax[i] = (short)(minmax[i]/scale+0.5);
03524 }
03525 }
03526 void G4GMocrenIO::setDoseDistMinMax(double _minmax[2], int _num) {
03527
03528 kDose[_num].setMinMax(_minmax);
03529 }
03530 void G4GMocrenIO::getDoseDistMinMax(double _minmax[2], int _num) {
03531
03532 if(isDoseEmpty())
03533 for(int i = 0; i < 2; i++) _minmax[i] = 0.;
03534 else
03535 kDose[_num].getMinMax(_minmax);
03536 }
03537
03538
03539 void G4GMocrenIO::setDoseDistScale(double & _scale, int _num) {
03540
03541 kDose[_num].setScale(_scale);
03542 }
03543 double G4GMocrenIO::getDoseDistScale(int _num) {
03544
03545 if(isDoseEmpty())
03546 return 0.;
03547 else
03548 return kDose[_num].getScale();
03549 }
03550
03551
03552
03553
03554
03555
03556
03557
03558
03559
03560 void G4GMocrenIO::setShortDoseDist(short * _image, int _num) {
03561
03562 int size[3];
03563 kDose[_num].getSize(size);
03564 int dsize = size[0]*size[1];
03565 double * ddata = new double[dsize];
03566 double scale = kDose[_num].getScale();
03567 double minmax[2];
03568 kDose[_num].getMinMax(minmax);
03569 for(int xy = 0; xy < dsize; xy++) {
03570 ddata[xy] = _image[xy]*scale;
03571 if(ddata[xy] < minmax[0]) minmax[0] = ddata[xy];
03572 if(ddata[xy] > minmax[1]) minmax[1] = ddata[xy];
03573 }
03574 kDose[_num].addImage(ddata);
03575
03576
03577 kDose[_num].setMinMax(minmax);
03578 }
03579 void G4GMocrenIO::getShortDoseDist(short * _data, int _z, int _num) {
03580
03581 if(_data == NULL) {
03582 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
03583 G4cout << "In G4GMocrenIO::getShortDoseDist(), "
03584 << "first argument is NULL pointer. "
03585 << "The argument must be allocated array."
03586 << G4endl;
03587 std::exit(-1);
03588 }
03589
03590 int size[3];
03591 kDose[_num].getSize(size);
03592
03593 double * ddata = kDose[_num].getImage(_z);
03594 double scale = kDose[_num].getScale();
03595 for(int xy = 0; xy < size[0]*size[1]; xy++) {
03596 _data[xy] = (short)(ddata[xy]/scale+0.5);
03597 }
03598 }
03599 void G4GMocrenIO::getShortDoseDistMinMax(short _minmax[2], int _num) {
03600 double scale = kDose[_num].getScale();
03601 double minmax[2];
03602 kDose[_num].getMinMax(minmax);
03603 for(int i = 0; i < 2; i++)
03604 _minmax[i] = (short)(minmax[i]/scale+0.5);
03605 }
03606
03607 void G4GMocrenIO::setDoseDist(double * _image, int _num) {
03608
03609 kDose[_num].addImage(_image);
03610 }
03611 double * G4GMocrenIO::getDoseDist(int _z, int _num) {
03612
03613 double * image;
03614 if(isDoseEmpty()) {
03615 image = 0;
03616 } else {
03617 image = kDose[_num].getImage(_z);
03618 }
03619 return image;
03620 }
03621
03622
03623
03624
03625
03626
03627
03628
03629
03630
03631
03632
03633
03634 bool G4GMocrenIO::addDoseDist(std::vector<double *> & _image, int _num) {
03635
03636 int size[3];
03637 getDoseDistSize(size, _num);
03638 std::vector<double *> dosedist = kDose[_num].getImage();
03639
03640 int nimg = size[0]*size[1];
03641 for(int z = 0; z < size[2]; z++) {
03642 for(int xy = 0; xy < nimg; xy++) {
03643 dosedist[z][xy] += _image[z][xy];
03644 }
03645 }
03646
03647 return true;
03648 }
03649
03650
03651 void G4GMocrenIO::setDoseDistCenterPosition(float _center[3], int _num) {
03652
03653 kDose[_num].setCenterPosition(_center);
03654 }
03655 void G4GMocrenIO::getDoseDistCenterPosition(float _center[3], int _num) {
03656
03657 if(isDoseEmpty())
03658 for(int i = 0; i < 3; i++) _center[i] = 0;
03659 else
03660 kDose[_num].getCenterPosition(_center);
03661 }
03662
03663 void G4GMocrenIO::setDoseDistName(std::string _name, int _num) {
03664
03665 kDose[_num].setName(_name);
03666 }
03667 std::string G4GMocrenIO::getDoseDistName(int _num) {
03668
03669 std::string name;
03670 if(isDoseEmpty())
03671 return name;
03672 else
03673 return kDose[_num].getName();
03674 }
03675
03676 void G4GMocrenIO::copyDoseDist(std::vector<class GMocrenDataPrimitive<double> > & _dose) {
03677 std::vector<class GMocrenDataPrimitive<double> >::iterator itr;
03678 for(itr = kDose.begin(); itr != kDose.end(); itr++) {
03679 _dose.push_back(*itr);
03680 }
03681 }
03682
03683 bool G4GMocrenIO::mergeDoseDist(std::vector<class GMocrenDataPrimitive<double> > & _dose) {
03684 if(kDose.size() != _dose.size()) {
03685 if (G4VisManager::GetVerbosity() >= G4VisManager::errors) {
03686 G4cout << "G4GMocrenIO::mergeDoseDist() : Error" << G4endl;
03687 G4cout << " Unable to merge the dose distributions,"<< G4endl;
03688 G4cout << " because of different size of dose maps."<< G4endl;
03689 }
03690 return false;
03691 }
03692
03693 int num = kDose.size();
03694 std::vector<class GMocrenDataPrimitive<double> >::iterator itr1 = kDose.begin();
03695 std::vector<class GMocrenDataPrimitive<double> >::iterator itr2 = _dose.begin();
03696 for(int i = 0; i < num; i++, itr1++, itr2++) {
03697 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
03698 if(kVerbose > 0)
03699 G4cout << "merged dose distribution [" << i << "]" << G4endl;
03700 *itr1 += *itr2;
03701 }
03702
03703 return true;
03704 }
03705
03706 void G4GMocrenIO::clearDoseDistAll() {
03707
03708 if(!isDoseEmpty()) {
03709 for(int i = 0; i < getNumDoseDist(); i++) {
03710 kDose[i].clear();
03711 }
03712 kDose.clear();
03713 }
03714 }
03715
03716 bool G4GMocrenIO::isDoseEmpty() {
03717 if(kDose.empty()) {
03718
03719
03720 return true;
03721 } else {
03722 return false;
03723 }
03724 }
03725
03726
03727 void G4GMocrenIO::calcDoseDistScale() {
03728
03729 double scale;
03730 double minmax[2];
03731
03732 for(int i = 0; i < (int)kDose.size(); i++) {
03733 kDose[i].getMinMax(minmax);
03734 scale = minmax[1]/DOSERANGE;
03735 kDose[i].setScale(scale);
03736 }
03737 }
03738
03739
03740
03741
03742
03743 void G4GMocrenIO::newROI() {
03744 GMocrenDataPrimitive<short> roiData;
03745 kRoi.push_back(roiData);
03746 }
03747 int G4GMocrenIO::getNumROI() {
03748 return (int)kRoi.size();
03749 }
03750
03751
03752 void G4GMocrenIO::setROIScale(double & _scale, int _num) {
03753
03754 kRoi[_num].setScale(_scale);
03755 }
03756 double G4GMocrenIO::getROIScale(int _num) {
03757
03758 if(isROIEmpty())
03759 return 0.;
03760 else
03761 return kRoi[_num].getScale();
03762 }
03763
03764 void G4GMocrenIO::setROI(short * _image, int _num) {
03765
03766 kRoi[_num].addImage(_image);
03767 }
03768 short * G4GMocrenIO::getROI(int _z, int _num) {
03769
03770 if(isROIEmpty())
03771 return 0;
03772 else
03773 return kRoi[_num].getImage(_z);
03774 }
03775
03776 void G4GMocrenIO::setROISize(int _size[3], int _num) {
03777
03778 return kRoi[_num].setSize(_size);
03779 }
03780 void G4GMocrenIO::getROISize(int _size[3], int _num) {
03781
03782 if(isROIEmpty())
03783 for(int i = 0; i < 3; i++) _size[i] = 0;
03784 else
03785 return kRoi[_num].getSize(_size);
03786 }
03787
03788 void G4GMocrenIO::setROIMinMax(short _minmax[2], int _num) {
03789
03790 kRoi[_num].setMinMax(_minmax);
03791 }
03792 void G4GMocrenIO::getROIMinMax(short _minmax[2], int _num) {
03793
03794 if(isROIEmpty())
03795 for(int i = 0; i < 2; i++) _minmax[i] = 0;
03796 else
03797 kRoi[_num].getMinMax(_minmax);
03798 }
03799
03800 void G4GMocrenIO::setROICenterPosition(float _center[3], int _num) {
03801
03802 kRoi[_num].setCenterPosition(_center);
03803 }
03804 void G4GMocrenIO::getROICenterPosition(float _center[3], int _num) {
03805
03806 if(isROIEmpty())
03807 for(int i = 0; i < 3; i++) _center[i] = 0;
03808 else
03809 kRoi[_num].getCenterPosition(_center);
03810 }
03811
03812 void G4GMocrenIO::clearROIAll() {
03813
03814 if(!isROIEmpty()) {
03815 for(int i = 0; i < getNumROI(); i++) {
03816 kRoi[i].clear();
03817 }
03818 kRoi.clear();
03819 }
03820 }
03821
03822 bool G4GMocrenIO::isROIEmpty() {
03823 if(kRoi.empty()) {
03824
03825
03826 return true;
03827 } else {
03828 return false;
03829 }
03830 }
03831
03832
03833
03834
03835
03836 int G4GMocrenIO::getNumTracks() {
03837 return (int)kSteps.size();
03838 }
03839 int G4GMocrenIO::getNumTracks4() {
03840 return (int)kTracks.size();
03841 }
03842 void G4GMocrenIO::addTrack(float * _tracks) {
03843 kSteps.push_back(_tracks);
03844 }
03845 void G4GMocrenIO::setTracks(std::vector<float *> & _tracks) {
03846 kSteps = _tracks;
03847 }
03848 std::vector<float *> & G4GMocrenIO::getTracks() {
03849 return kSteps;
03850 }
03851 void G4GMocrenIO::addTrackColor(unsigned char * _colors) {
03852 kStepColors.push_back(_colors);
03853 }
03854 void G4GMocrenIO::setTrackColors(std::vector<unsigned char *> & _trackColors) {
03855 kStepColors = _trackColors;
03856 }
03857 std::vector<unsigned char *> & G4GMocrenIO::getTrackColors() {
03858 return kStepColors;
03859 }
03860 void G4GMocrenIO::copyTracks(std::vector<float *> & _tracks,
03861 std::vector<unsigned char *> & _colors) {
03862 std::vector<float *>::iterator titr;
03863 for(titr = kSteps.begin(); titr != kSteps.end(); titr++) {
03864 float * pts = new float[6];
03865 for(int i = 0; i < 6; i++) {
03866 pts[i] = (*titr)[i];
03867 }
03868 _tracks.push_back(pts);
03869 }
03870
03871 std::vector<unsigned char *>::iterator citr;
03872 for(citr = kStepColors.begin(); citr != kStepColors.end(); citr++) {
03873 unsigned char * pts = new unsigned char[3];
03874 for(int i = 0; i < 3; i++) {
03875 pts[i] = (*citr)[i];
03876 }
03877 _colors.push_back(pts);
03878 }
03879 }
03880 void G4GMocrenIO::mergeTracks(std::vector<float *> & _tracks,
03881 std::vector<unsigned char *> & _colors) {
03882 std::vector<float *>::iterator titr;
03883 for(titr = _tracks.begin(); titr != _tracks.end(); titr++) {
03884 addTrack(*titr);
03885 }
03886
03887 std::vector<unsigned char *>::iterator citr;
03888 for(citr = _colors.begin(); citr != _colors.end(); citr++) {
03889 addTrackColor(*citr);
03890 }
03891 }
03892 void G4GMocrenIO::addTrack(std::vector<float *> & _steps, unsigned char _color[3]) {
03893
03894 std::vector<float *>::iterator itr = _steps.begin();
03895 std::vector<struct GMocrenTrack::Step> steps;
03896 for(; itr != _steps.end(); itr++) {
03897 struct GMocrenTrack::Step step;
03898 for(int i = 0; i < 3; i++) {
03899 step.startPoint[i] = (*itr)[i];
03900 step.endPoint[i] = (*itr)[i+3];
03901 }
03902 steps.push_back(step);
03903 }
03904 GMocrenTrack track;
03905 track.setTrack(steps);
03906 track.setColor(_color);
03907 kTracks.push_back(track);
03908
03909 }
03910 void G4GMocrenIO::getTrack(int _num, std::vector<float *> & _steps,
03911 std::vector<unsigned char *> & _color) {
03912
03913 if(_num > (int)kTracks.size()) {
03914 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
03915 G4cout << "ERROR in getTrack() : " << G4endl;
03916 std::exit(-1);
03917 }
03918 unsigned char * color = new unsigned char[3];
03919 kTracks[_num].getColor(color);
03920 _color.push_back(color);
03921
03922
03923 int nsteps = kTracks[_num].getNumberOfSteps();
03924 for(int isteps = 0; isteps < nsteps; isteps++) {
03925 float * stepPoints = new float[6];
03926 kTracks[_num].getStep(stepPoints[0], stepPoints[1], stepPoints[2],
03927 stepPoints[3], stepPoints[4], stepPoints[5],
03928 isteps);
03929 _steps.push_back(stepPoints);
03930 }
03931 }
03932
03933 void G4GMocrenIO::translateTracks(std::vector<float> & _translate) {
03934 std::vector<class GMocrenTrack>::iterator itr = kTracks.begin();
03935 for(; itr != kTracks.end(); itr++) {
03936 itr->translate(_translate);
03937 }
03938 }
03939
03940
03941
03942
03943
03944 int G4GMocrenIO::getNumberOfDetectors() {
03945 return (int)kDetectors.size();
03946 }
03947 void G4GMocrenIO::addDetector(std::string & _name,
03948 std::vector<float *> & _det,
03949 unsigned char _color[3]) {
03950
03951 std::vector<float *>::iterator itr = _det.begin();
03952 std::vector<struct GMocrenDetector::Edge> edges;
03953 for(; itr != _det.end(); itr++) {
03954 struct GMocrenDetector::Edge edge;
03955 for(int i = 0; i < 3; i++) {
03956 edge.startPoint[i] = (*itr)[i];
03957 edge.endPoint[i] = (*itr)[i+3];
03958 }
03959 edges.push_back(edge);
03960 }
03961 GMocrenDetector detector;
03962 detector.setDetector(edges);
03963 detector.setColor(_color);
03964 detector.setName(_name);
03965 kDetectors.push_back(detector);
03966
03967 }
03968
03969 void G4GMocrenIO::getDetector(int _num, std::vector<float *> & _edges,
03970 std::vector<unsigned char *> & _color,
03971 std::string & _detName) {
03972
03973 if(_num > (int)kDetectors.size()) {
03974 if (G4VisManager::GetVerbosity() >= G4VisManager::errors)
03975 G4cout << "ERROR in getDetector() : " << G4endl;
03976 std::exit(-1);
03977 }
03978
03979 _detName = kDetectors[_num].getName();
03980
03981 unsigned char * color = new unsigned char[3];
03982 kDetectors[_num].getColor(color);
03983 _color.push_back(color);
03984
03985
03986 int nedges = kDetectors[_num].getNumberOfEdges();
03987 for(int ne = 0; ne < nedges; ne++) {
03988 float * edgePoints = new float[6];
03989 kDetectors[_num].getEdge(edgePoints[0], edgePoints[1], edgePoints[2],
03990 edgePoints[3], edgePoints[4], edgePoints[5],
03991 ne);
03992 _edges.push_back(edgePoints);
03993 }
03994 }
03995
03996 void G4GMocrenIO::translateDetector(std::vector<float> & _translate) {
03997 std::vector<class GMocrenDetector>::iterator itr = kDetectors.begin();
03998 for(; itr != kDetectors.end(); itr++) {
03999 itr->translate(_translate);
04000 }
04001 }
04002
04003
04004 template <typename T>
04005 void G4GMocrenIO::convertEndian(char * _val, T & _rval) {
04006
04007 if((kLittleEndianOutput && !kLittleEndianInput) ||
04008 (!kLittleEndianOutput && kLittleEndianInput)) {
04009
04010 const int SIZE = sizeof(_rval);
04011 char ctemp;
04012 for(int i = 0; i < SIZE/2; i++) {
04013 ctemp = _val[i];
04014 _val[i] = _val[SIZE - 1 - i];
04015 _val[SIZE - 1 - i] = ctemp;
04016 }
04017 }
04018 _rval = *(T *)_val;
04019 }
04020
04021
04022 template <typename T>
04023 void G4GMocrenIO::invertByteOrder(char * _val, T & _rval) {
04024
04025 const int SIZE = sizeof(_rval);
04026
04027 union {
04028 char cu[16];
04029 T tu;
04030 } uni;
04031 for(int i = 0; i < SIZE; i++) {
04032 uni.cu[i] = _val[SIZE-1-i];
04033
04034 }
04035
04036 _rval = uni.tu;
04037
04038 }
04039
04040
04041 void G4GMocrenIO::setVerboseLevel(int _level) {
04042 kVerbose = _level;
04043 }
04044