Example/Dense2007ToDicom.cxx

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