Example/Dense2007ToDicom.cxx

   1 /*=========================================================================
   2
   3   Program:   gdcm
   4   Module:    $RCSfile: Dense2007ToDicom.cxx,v $
   5   Language:  C++
   6   Date:      $Date: 2007/10/24 08:03:10 $
   7   Version:   $Revision: 1.4 $
   8
   9   Copyright (c) CREATIS (Centre de Recherche et d'Applications en Traitement de
  10   l'Image). All rights reserved. See Doc/License.txt or
  11   http://www.creatis.insa-lyon.fr/Public/Gdcm/License.html for details.
  12
  13      This software is distributed WITHOUT ANY WARRANTY; without even
  14      the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  15      PURPOSE.  See the above copyright notices for more information.
  16
  17 =========================================================================*/
  18
  19 #include <fstream>
  20 #include <iostream>
  21 //#include <values.h>
  22
  23 #include "gdcmFile.h"
  24 #include "gdcmFileHelper.h"
  25 #include "gdcmDebug.h"
  26 #include "gdcmDirList.h"
  27
  28 #include "gdcmArgMgr.h"
  29
  30 /**
  31   * \brief
  32   *  Converts the "Dense" ".txt" (2007 version)  files into 16 bits Dicom Files,
  33   * Hope they don't change soon!
  34   */
  35
  36
  37 void LoadPeakStrain(std::ifstream &from, std::string imageName);
  38 void LoadStrain(std::ifstream &from, std::string imageName);
  39 void MakeDicomImage(float *tabVal, float *X, float *Y, float *Z, int NP, std::string dcmImageName);
  40
  41 int verbose;
  42
  43 int main(int argc, char *argv[])
  44 {
  45    START_USAGE(usage)
  46    " \n Dense2007ToDicom :\n                                                  ",
  47    "        Converts the '.txt' files into 16 bits Dicom Files,               ",
  48    " usage:                                                                   ",
  49    " Dense2007ToDicom strain=...strain.txt  peak_strain=...peak_strain.txt    ",
  50    "                 [verbose] [debug]                                        ",
  51    "                                                                          ",
  52    " verbose  : user wants to run the program in 'verbose mode'               ",
  53    " debug    : *developer*  wants to run the program in 'debug mode'         ",
  54    FINISH_USAGE
  55
  56    // ----- Initialize Arguments Manager ------
  57
  58    GDCM_NAME_SPACE::ArgMgr *am = new GDCM_NAME_SPACE::ArgMgr(argc, argv);
  59
  60    if (argc == 1 || am->ArgMgrDefined("usage"))
  61    {
  62       am->ArgMgrUsage(usage); // Display 'usage'
  63       delete am;
  64       return 0;
  65    }
  66
  67    const char *strain      = am->ArgMgrWantString("strain",usage);
  68    const char *peak_strain = am->ArgMgrWantString("peak_strain",usage);
  69
  70    if (am->ArgMgrDefined("debug"))
  71       GDCM_NAME_SPACE::Debug::DebugOn();
  72
  73    verbose  = am->ArgMgrDefined("verbose");
  74
  75    // if unused Param we give up
  76    if ( am->ArgMgrPrintUnusedLabels() )
  77    {
  78       am->ArgMgrUsage(usage);
  79       delete am;
  80       return 0;
  81    }
  82    delete am;  // we don't need Argument Manager any longer
  83
  84    // ----- Begin Processing -----
  85
  86
  87    std::ifstream fromPeakStrain( peak_strain );
  88    if ( !fromPeakStrain )
  89    {
  90       std::cout << "Can't open file" << peak_strain << std::endl;
  91       exit(0);
  92    }
  93
  94    std::ifstream fromStrain( strain );
  95    if ( !fromStrain )
  96    {
  97       std::cout << "Can't open file" << strain << std::endl;
  98       exit(0);
  99    }
 100
 101    std::cout << "Success in open file" << peak_strain << std::endl;
 102    LoadPeakStrain(fromPeakStrain, peak_strain);
 103    fromPeakStrain.close();
 104
 105    std::cout << "Success in open file" << strain << std::endl;
 106    LoadStrain(fromStrain, strain);
 107    fromStrain.close();
 108    return 1;
 109 }
 110
 111 // =====================================================================================================================
 112
 113 void LoadPeakStrain(std::ifstream &from, std::string textFileName)
 114 {
 115 // in sax_base_slice0_peak_strain.txt :
 116
 117 /*
 118 Number of material points (NP) = 181
 119 Origin of (readout, phase enc, slice sel) coordinates in 3D =  87.3243 3.19392
 120 88.2381
 121 Readout direction in 3D =  -0.162314 -0.0771294 -0.983720
 122 Phase Enc. direction in 3D =  -0.540606 -0.827052 0.154046
 123 Slice select direction in 3D =  0.825469 -0.556809 -0.0925458
 124 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
 125 followed by their peak Ecc strain, an array of NP elements,
 126 followed by their peak Err strain, an array of NP elements,
 127 followed by their peak E11 strain, an array of NP elements,
 128 followed by their Peak E22 strain, an array of NP elements,
 129       42.0000      10.0000     0.000000
 130       ...
 131     -0.154905   -0.0840482    -0.157350    -0.221403    -0.168118    -0.131331
 132     -0.153781    -0.148481    -0.166602    -0.232858    -0.222650    -0.213712
 133     ...
 134 */
 135
 136    if (!from)
 137       return;
 138
 139    std::string str1;
 140    int NP;
 141
 142    //Number of material points (NP) = 181
 143     from >> str1;
 144     from >> str1;
 145     from >> str1;
 146     from >> str1;
 147     from >> str1;
 148     from >> str1;
 149     from >> NP;
 150
 151     std::cout << "NP : " << NP << std::endl;
 152
 153    //Origin of (readout, phase enc, slice sel) coordinates in 3D =  87.3243 3.19392 88.2381
 154      from >> str1;
 155      from >> str1;
 156      from >> str1;
 157      from >> str1;
 158      from >> str1;
 159      from >> str1;
 160      from >> str1;
 161      from >> str1;
 162      from >> str1;
 163      from >> str1;
 164      from >> str1;
 165
 166      float readout,  phase_enc, slice_sel;
 167      from >> readout;
 168      from >> phase_enc;
 169      from >> slice_sel;
 170      std::cout << " readout " << readout << " phase_enc " << phase_enc << " slice_sel " << slice_sel << std::endl;
 171
 172     // Readout direction in 3D =  -0.162314 -0.0771294 -0.983720
 173
 174     from >> str1;
 175     from >> str1;
 176     from >> str1;
 177     from >> str1;
 178     from >> str1;
 179
 180     float readoutX, readoutY, readoutZ;
 181     from >> readoutX;
 182     from >> readoutY;
 183     from >> readoutZ;
 184     std::cout << " readoutX " << readoutX <<  " readoutY " << readoutY <<  " readoutZ " << readoutZ << std::endl;
 185
 186 // Phase Enc. direction in 3D =  -0.540606 -0.827052 0.154046
 187
 188      from >> str1;
 189      from >> str1;
 190      from >> str1;
 191      from >> str1;
 192      from >> str1;
 193      from >> str1;
 194
 195     float phase_encX, phase_encY, phase_encZ;
 196     from >> phase_encX;
 197     from >> phase_encY;
 198     from >> phase_encZ;
 199     std::cout << " phase_encX " << phase_encX <<  " phase_encY " << phase_encY <<  " phase_encZ " << phase_encZ << std::endl;
 200
 201 // Slice select direction in 3D =  0.825469 -0.556809 -0.0925458
 202      from >> str1;
 203      from >> str1;
 204      from >> str1;
 205      from >> str1;
 206      from >> str1;
 207      from >> str1;
 208
 209     float slice_selX, slice_selY, slice_selZ;
 210     from >> slice_selX;
 211     from >> slice_selY;
 212     from >> slice_selZ;
 213     std::cout << " slice_selX " << slice_selX <<  " slice_selY " << slice_selY <<  " slice_selZ " << slice_selZ << std::endl;
 214
 215
 216
 217 // Skip 5 lines :
 218 /*
 219 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
 220 followed by their peak Ecc strain, an array of NP elements,
 221 followed by their peak Err strain, an array of NP elements,
 222 followed by their peak E11 strain, an array of NP elements,
 223 followed by their Peak E22 strain, an array of NP elements,
 224 */
 225
 226 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 227    std::getline(from, str1);
 228    std::cout << "[" << str1 << "]" << std::endl;
 229 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 230    std::getline(from, str1);
 231    std::cout << "[" << str1 << "]" << std::endl;
 232 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 233    std::getline(from, str1);
 234    std::cout << "[" << str1 << "]" << std::endl;
 235 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 236    std::getline(from, str1);
 237    std::cout << "[" << str1 << "]" << std::endl;
 238 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 239    std::getline(from, str1);
 240    std::cout << "[" << str1 << "]" << std::endl;
 241 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 242    std::getline(from, str1);
 243    std::cout << "[" << str1 << "]" << std::endl;
 244 std::cout << "------------stop skipping ---------------- " << std::endl;
 245
 246    float *X = new float[NP];
 247    float *Y = new float[NP];
 248    float *Z = new float[NP];
 249
 250    char c;
 251    int i;
 252    for (i=0; i<NP; i++) {
 253
 254       from >> X[i];
 255       for (;;) {
 256         if (!from.get(c))
 257           break;
 258         if (!isspace(c)) {
 259           from.putback(c);
 260           break;
 261         }
 262      }
 263
 264       from >> Y[i];
 265       for (;;) {
 266         if (!from.get(c))
 267           break;
 268         if (!isspace(c)) {
 269           from.putback(c);
 270           break;
 271         }
 272      }
 273       from >> Z[i];
 274
 275    } // end for i<NP
 276
 277    std::cout << "--------------- Ecc_strain ------------------" << std::endl;
 278    float *ecc_strain = new float[NP];
 279    for (i=0; i<NP; i++) {
 280        from >> ecc_strain[i];
 281        if (verbose)
 282        std::cout <<  ecc_strain[i] <<  std::endl;
 283    }
 284
 285    std::cout << "--------------- Err_strain ------------------" << std::endl;
 286    float *err_strain = new float[NP];
 287    for (i=0; i<NP; i++) {
 288        from >> err_strain[i];
 289        if (verbose)
 290        std::cout <<  err_strain[i] <<  std::endl;
 291    }
 292
 293    std::cout << "--------------- E11_strain ------------------" << std::endl;
 294    float *e11_strain = new float[NP];
 295    for (i=0; i<NP; i++) {
 296        from >> e11_strain[i];
 297        if (verbose)
 298        std::cout <<  e11_strain[i] <<  std::endl;
 299    }
 300
 301    std::cout << "--------------- E22_strain ------------------" << std::endl;
 302    float *e22_strain = new float[NP];
 303    for (i=0; i<NP; i++) {
 304        from >> e22_strain[i];
 305        if (verbose)
 306        std::cout <<  e22_strain[i] <<  std::endl;
 307    }
 308
 309    std::string dcmImageName;
 310
 311 //followed by their peak Ecc strain, an array of NP elements,
 312    dcmImageName = textFileName + "_peak_Ecc_strain.dcm";
 313    MakeDicomImage(ecc_strain, X, Y, Z, NP, dcmImageName);
 314
 315 //followed by their peak Err strain, an array of NP elements,
 316    dcmImageName = textFileName + "_peak_Err_strain.dcm";
 317    MakeDicomImage(err_strain, X, Y, Z, NP, dcmImageName);
 318
 319 //followed by their peak E11 strain, an array of NP elements,
 320    dcmImageName = textFileName + "_peak_E11_strain.dcm";
 321    MakeDicomImage(e11_strain, X, Y, Z, NP, dcmImageName);
 322
 323 //followed by their Peak E22 strain, an array of NP elements,
 324    dcmImageName = textFileName + "_peak_E22_strain.dcm";
 325    MakeDicomImage(e22_strain, X, Y, Z, NP, dcmImageName);
 326 }
 327
 328 // =====================================================================================================================
 329
 330 void LoadStrain(std::ifstream &from, std::string textFileName)
 331 {
 332
 333 // in sax_base_slice0_strain.txt :
 334 /*
 335 Number of cine frames = 18
 336 Temporal resolution = 32.0000 ms
 337 First frame starts at 48.0000 ms
 338 Number of material points (NP) = 181
 339 Origin of (readout, phase enc, slice sel) coordinates in 3D =  87.324341 3.193918 88.238113
 340 Readout direction in 3D =  -0.162314 -0.0771294 -0.983720
 341 Phase Enc. direction in 3D =  -0.540606 -0.827052 0.154046
 342 Slice select direction in 3D =  0.825469 -0.556809 -0.0925458
 343 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
 344 followed by their Ecc strain, an array of dimensions(NP, number of cine frames),
 345 followed by their Err strain, an array of dimensions(NP, number of cine frames),
 346 followed by their E11 strain, an array of dimensions(NP, number of cine frames),
 347 followed by their E22 strain, an array of dimensions(NP, number of cine frames),
 348 Note that RV Err, E11 and E22 strains are not calculated due to the small thickness.
 349       42.0000      10.0000     0.000000
 350       44.0000      10.0000     0.000000
 351       ...
 352    -0.0622793   -0.0840482    -0.157350    -0.196722    -0.105844    -0.131331
 353     -0.153781  -0.00940573   -0.0542236    -0.100403   -0.0369671   -0.0696840
 354 */
 355
 356    if (!from)
 357       return;
 358
 359    std::string str1;
 360    int NP;    // Number of Pints
 361    int NCF;   // Number of cine frames
 362    float TR;  // Temporal resolution
 363    float FFS; // First frame starts
 364
 365    // Number of cine frames = 18
 366     from >> str1;
 367     from >> str1;
 368     from >> str1;
 369     from >> str1;
 370     from >> str1;
 371     from >> NCF;
 372
 373    // Temporal resolution = 32.0000 ms
 374     from >> str1;
 375     from >> str1;
 376     from >> str1;
 377     from >> TR;
 378     from >> str1;
 379
 380    // First frame starts at 48.0000 ms
 381     from >> str1;
 382     from >> str1;
 383     from >> str1;
 384     from >> str1;
 385     from >> FFS;
 386     from >> str1;
 387
 388    //Number of material points (NP) = 181
 389     from >> str1;
 390     from >> str1;
 391     from >> str1;
 392     from >> str1;
 393     from >> str1;
 394     from >> str1;
 395     from >> NP;
 396
 397     std::cout << "NP : " << NP << std::endl;
 398
 399    //Origin of (readout, phase enc, slice sel) coordinates in 3D = 87.324341 3.193918 88.238113
 400      from >> str1;
 401      from >> str1;
 402      from >> str1;
 403      from >> str1;
 404      from >> str1;
 405      from >> str1;
 406      from >> str1;
 407      from >> str1;
 408      from >> str1;
 409      from >> str1;
 410      from >> str1;
 411
 412      float readout,  phase_enc, slice_sel;
 413      from >> readout;
 414      from >> phase_enc;
 415      from >> slice_sel;
 416      std::cout << " readout " << readout << " phase_enc " << phase_enc << " slice_sel " << slice_sel << std::endl;
 417
 418     // Readout direction in 3D =  -0.162314 -0.0771294 -0.983720
 419
 420     from >> str1;
 421     from >> str1;
 422     from >> str1;
 423     from >> str1;
 424     from >> str1;
 425
 426     float readoutX, readoutY, readoutZ;
 427     from >> readoutX;
 428     from >> readoutY;
 429     from >> readoutZ;
 430     std::cout << " readoutX " << readoutX <<  " readoutY " << readoutY <<  " readoutZ " << readoutZ << std::endl;
 431
 432 // Phase Enc. direction in 3D =  -0.540606 -0.827052 0.154046
 433
 434      from >> str1;
 435      from >> str1;
 436      from >> str1;
 437      from >> str1;
 438      from >> str1;
 439      from >> str1;
 440
 441     float phase_encX, phase_encY, phase_encZ;
 442     from >> phase_encX;
 443     from >> phase_encY;
 444     from >> phase_encZ;
 445     std::cout << " phase_encX " << phase_encX <<  " phase_encY " << phase_encY <<  " phase_encZ " << phase_encZ << std::endl;
 446
 447 // Slice select direction in 3D =  0.825469 -0.556809 -0.0925458
 448      from >> str1;
 449      from >> str1;
 450      from >> str1;
 451      from >> str1;
 452      from >> str1;
 453      from >> str1;
 454
 455     float slice_selX, slice_selY, slice_selZ;
 456     from >> slice_selX;
 457     from >> slice_selY;
 458     from >> slice_selZ;
 459     std::cout << " slice_selX " << slice_selX <<  " slice_selY " << slice_selY <<  " slice_selZ " << slice_selZ << std::endl;
 460
 461
 462
 463 // Skip 6 lines :
 464 /*
 465 The following are the (readout, phase enc, slice sel) coordinates (mm) of the grid points for which strains are calculated,
 466 followed by their Ecc strain, an array of dimensions(NP, number of cine frames),
 467 followed by their Err strain, an array of dimensions(NP, number of cine frames),
 468 followed by their E11 strain, an array of dimensions(NP, number of cine frames),
 469 followed by their E22 strain, an array of dimensions(NP, number of cine frames),
 470 Note that RV Err, E11 and E22 strains are not calculated due to the small thickness.
 471 */
 472 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 473    std::getline(from, str1);
 474    std::cout << "[" << str1 << "]" << std::endl;
 475 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 476    std::getline(from, str1);
 477    std::cout << "[" << str1 << "]" << std::endl;
 478 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 479    std::getline(from, str1);
 480    std::cout << "[" << str1 << "]" << std::endl;
 481 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 482    std::getline(from, str1);
 483    std::cout << "[" << str1 << "]" << std::endl;
 484 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 485    std::getline(from, str1);
 486    std::cout << "[" << str1 << "]" << std::endl;
 487 std::cout << "------------start skipping 1 line---------------- " << std::endl;
 488    std::getline(from, str1);
 489    std::cout << "[" << str1 << "]" << std::endl;
 490 std::cout << "------------stop skipping ---------------- " << std::endl;
 491    std::getline(from, str1);
 492    std::cout << "[" << str1 << "]" << std::endl;
 493 std::cout << "------------stop skipping ---------------- " << std::endl;
 494
 495    float *X = new float[NP];
 496    float *Y = new float[NP];
 497    float *Z = new float[NP];
 498
 499    char c;
 500    int i;
 501    for (i=0; i<NP; i++) {
 502
 503       from >> X[i];
 504       for (;;) {
 505         if (!from.get(c))
 506           break;
 507         if (!isspace(c)) {
 508           from.putback(c);
 509           break;
 510         }
 511      }
 512       from >> Y[i];
 513       for (;;) {
 514         if (!from.get(c))
 515           break;
 516         if (!isspace(c)) {
 517           from.putback(c);
 518           break;
 519         }
 520      }
 521       from >> Z[i];
 522    }
 523
 524    std::cout << "--------------- Ecc_strain ------------------" << std::endl;
 525    float *ecc_strain = new float[NP];
 526    for (i=0; i<NP; i++) {
 527        from >> ecc_strain[i];
 528        if (verbose)
 529        std::cout <<  ecc_strain[i] <<  std::endl;
 530    }
 531
 532    std::cout << "--------------- Err_strain ------------------" << std::endl;
 533    float *err_strain = new float[NP];
 534    for (i=0; i<NP; i++) {
 535        from >> err_strain[i];
 536        if (verbose)
 537        std::cout <<  err_strain[i] <<  std::endl;
 538    }
 539
 540    std::cout << "--------------- E11_strain ------------------" << std::endl;
 541    float *e11_strain = new float[NP];
 542    for (i=0; i<NP; i++) {
 543        from >> e11_strain[i];
 544        if (verbose)
 545        std::cout <<  e11_strain[i] <<  std::endl;
 546    }
 547
 548    std::cout << "--------------- E22_strain ------------------" << std::endl;
 549    float *e22_strain = new float[NP];
 550    for (i=0; i<NP; i++) {
 551        from >> e22_strain[i];
 552        if (verbose)
 553        std::cout <<  e22_strain[i] <<  std::endl;
 554    }
 555
 556    std::string dcmImageName;
 557
 558 //followed by their Ecc strain, an array of NP elements,
 559    dcmImageName = textFileName + "_Ecc_strain.dcm";
 560    MakeDicomImage(ecc_strain, X, Y, Z, NP, dcmImageName);
 561
 562 //followed by their  Err strain, an array of NP elements,
 563    dcmImageName = textFileName + "_Err_strain.dcm";
 564    MakeDicomImage(err_strain, X, Y, Z, NP, dcmImageName);
 565
 566 //followed by their E11 strain, an array of NP elements,
 567    dcmImageName = textFileName + "_E11_strain.dcm";
 568    MakeDicomImage(e11_strain, X, Y, Z, NP, dcmImageName);
 569
 570 //followed by their E22 strain, an array of NP elements,
 571    dcmImageName = textFileName + "_E22_strain.dcm";
 572    MakeDicomImage(e22_strain, X, Y, Z, NP, dcmImageName);
 573
 574 }
 575
 576
 577 // =====================================================================================================================
 578
 579
 580 void MakeDicomImage(float *tabVal, float *X, float *Y, float *Z, int NP, std::string dcmImageName)
 581 {
 582    float minX = 99999., minY = 99999., minZ = 99999.;
 583    float maxX = 0., maxY = 0., maxZ = 0.;
 584    int i;
 585
 586    for (i=0; i<NP; i++) {
 587       // std::cout << X[i] << " " << Y[i] << " " << Z[i] <<  std::endl;
 588       if(maxX < X[i])
 589          maxX = X[i];
 590       if(maxY < Y[i])
 591          maxY = Y[i];
 592       if(maxZ < Z[i])
 593          maxZ = Z[i];
 594
 595       if(minX > X[i])
 596          minX = X[i];
 597       if(minY > Y[i])
 598          minY = Y[i];
 599       if(minZ > Z[i])
 600          minZ = Z[i];
 601    }
 602    std::cout << "Min X,Y,Z " << minX << " " << minY << " " << minZ <<  std::endl;
 603    std::cout << "Max X,Y,Z " << maxX << " " << maxY << " " << maxZ <<  std::endl;
 604    std::cout << "Size X,Y,Z " << maxX-minX << " " << maxY-minY << " " << maxZ-minZ <<  std::endl;
 605
 606 //   uint16_t *img = new uint16_t[int(maxX+0.5)*int(maxY+0.5)];
 607    uint16_t *img = new uint16_t[int(maxX*4.)*int(maxY*4.)];
 608
 609    // Set whole image to 0
 610    for(int i3=0;i3<int(maxX*4.)*int(maxY*4.);i3++)
 611 //   for(int i3=0;i3<int(maxX+0.5)*int(maxY+0.5);i3++)
 612       img[i3] = 0;
 613
 614    for(int i2=0; i2<NP; i2++) {
 615
 616       int ordX = int(X[i2]*4.-30);
 617       int ordY = int(maxY*4.) - int(Y[i2]*4.)+30;
 618 //   img[ /*int(maxX) -*/ int(X[i2]+0.5-1)   +  (int(maxY+0.5) - int(Y[i2]+0.5))   * int(maxX+0.5) ] = int(tabVal[i2]*100);
 619
 620 //   img[ /*int(maxX) -*/ int(X[i2]*4.-30)   +  (int(maxY*4.) - int(Y[i2]*4.)+30)   * int(maxX*4.) ] = int(tabVal[i2]*100);
 621    img[ /*int(maxX) -*/ ordX   +  ordY   * int(maxX*4.) ] = int(tabVal[i2]*100);
 622
 623       // Try to round up, just to see.
 624        for(int iii=ordY-3; iii<ordY+4; iii++) for(int jjj=ordX-3; jjj<ordX+4; jjj++)
 625           img[  jjj  +  iii   * int(maxX*4.) ] = int(tabVal[i2]*100);
 626       std::cout << int(X[i2]*4.) << " " << int(Y[i2]*4.) << " = " << int(tabVal[i2]*100) << std::endl;
 627    }
 628
 629  // GDCM_NAME_SPACE::Debug::DebugOn();
 630
 631    std::ostringstream str;
 632
 633    GDCM_NAME_SPACE::File *file;
 634    file = GDCM_NAME_SPACE::File::New();
 635
 636   // Set the image size
 637    str.str("");
 638    str << (int)(maxX*4.);
 639    file->InsertEntryString(str.str(),0x0028,0x0011,"US"); // Columns
 640    str.str("");
 641    str << (int)(maxY*4.);
 642    file->InsertEntryString(str.str(),0x0028,0x0010,"US"); // Rows
 643
 644   // Set the pixel type
 645   //      16; //8, 16, 32
 646    file->InsertEntryString("16",0x0028,0x0100,"US"); // Bits Allocated
 647    str.str("");
 648    str << 16; // may be 12 or 16 if componentSize =16
 649    file->InsertEntryString("16",0x0028,0x0101,"US"); // Bits Stored
 650    file->InsertEntryString("15",0x0028,0x0102,"US"); // High Bit
 651
 652   // Set the pixel representation // 0/1 , 0=unsigned
 653    file->InsertEntryString("1",0x0028,0x0103, "US"); // Pixel Representation
 654   // Set the samples per pixel // 1:Grey level, 3:RGB
 655    file->InsertEntryString("1",0x0028,0x0002, "US"); // Samples per Pixel
 656
 657 /*
 658   // Set Rescale Intercept
 659         str.str("");
 660         str << div;
 661         file->InsertEntryString(str.str(),0x0028,0x1052,"DS");
 662
 663   // Set Rescale Slope
 664         str.str("");
 665         str << mini;
 666         file->InsertEntryString(str.str(),0x0028,0x1053,"DS");
 667 */
 668
 669    GDCM_NAME_SPACE::FileHelper *fileH;
 670    fileH = GDCM_NAME_SPACE::FileHelper::New(file);
 671    // cast is just to avoid warnings (*no* conversion)
 672    //fileH->SetImageData((uint8_t *)img,int(maxX*maxY)*sizeof(uint16_t)); // troubles when maxX, mayY are *actually* float!
 673    fileH->SetImageData((uint8_t *)img,int(maxX*4.)*int(maxY*4.)*sizeof(uint16_t));
 674    fileH->SetWriteModeToRaw();
 675    fileH->SetWriteTypeToDcmExplVR();
 676
 677    if( !fileH->Write(dcmImageName))
 678       std::cout << "Failed for [" << dcmImageName << "]\n"
 679                 << "           File is unwrittable" << std::endl;
 680
 681    //file->Print();
 682
 683    delete img;
 684    file->Delete();
 685    fileH->Delete();
 686 }