Example/DenseMultiFramesToDicom.cxx

   1 /*=========================================================================
   2
   3   Program:   gdcm
   4   Module:    $RCSfile: DenseMultiFramesToDicom.cxx,v $
   5   Language:  C++
   6   Date:      $Date: 2007/06/27 08:38:44 $
   7   Version:   $Revision: 1.7 $
   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 #include "gdcmUtil.h"
  28
  29 #include "gdcmArgMgr.h"
  30
  31 /**
  32   * \brief
  33   *          - explores recursively the (single Patient, single Study) directory
  34   *          - examines the ".txt" Dense multiframe files
  35   *          - Converts the files into 16 bits Dicom Files
  36   *           WARNING : directory must contain ONLY .txt files
  37   */
  38
  39 void Load(std::ifstream &from, std::string imageName, std::string patName,
  40           std::string strStudyUID, std::string strSerieUID, int serieNumber, int verbose);
  41
  42 //std::ifstream& eatwhite(std::ifstream& is);
  43
  44 /*
  45 window center (level) and width, are defined in the DICOM
  46 standard very precisely as follows (see PS 3.3 C.11.2.1.2):
  47 >
  48 >"These Attributes are applied according to the following pseudo-code,
  49 >where :
  50 x is the input value,
  51 y is an output value with a range from ymin to ymax,
  52 c is Window Center (0028,1050)
  53 w is Window Width (0028,1051):
  54 >
  55 >           if      (x  <= c - 0.5 - (w-1)/2), then y = ymin
  56 >           else if (x > c - 0.5 + (w-1)/2), then y = ymax,
  57 >           else    y = ((x - (c - 0.5)) / (w-1) + 0.5) * (ymax - ymin)+ ymin
  58
  59 */
  60 /*
  61 From:   David Clunie - view profile
  62 Date:   Thurs, Jun 1 2006 3:03 pm
  63 Email:  David Clunie <dclu...@dclunie.com>
  64 Groups: comp.protocols.dicom
  65
  66 The value of x is the output of the preceding step, the so-called
  67 "modality LUT", which may either be:
  68
  69 - identity (no rescale values or Modality LUT, or the SOP Class is
  70   PET and rescale values are ignored), in which case x is the stored
  71   pixel value in (7FE0,0010)
  72
  73 - Rescale Slope and Intercept (as typically used in CT), in which
  74   case x is the value obtained from applying the rescale values to
  75   the stored pixel value
  76
  77 - an actual LUT, in which case x is the value stored in the LUT
  78   corresponding to the LUT index value that is the stored pixel
  79   value
  80
  81 The ymin and ymax are intended to represent the output range; for
  82 example, if the hypothetical Presentation LUT step that follows
  83 the VOI LUT (window) stage is an identity operation, then the
  84 ymin and ymax represent P-Values, which might be the range of
  85 digital driving levels for your display (calibrated to the GSDF),
  86 in the 8-bit wide output case ranging from 0 to 255, for example.
  87
  88 The terms brightness and contrast are not used in radiology imaging
  89 -the window center and width are used instead.
  90 */
  91
  92 int main(int argc, char *argv[])
  93 {
  94    START_USAGE(usage)
  95    " \n DenseMultiFramesToDicom :                                           \n",
  96    " - explores recursively the given (single Patient, single Study) directory",
  97    "         - examines the '.txt' files                                      ",
  98    "         - Converts the files into 16 bits Dicom files,                   ",
  99    " WARNING : directory must contain ONLY .txt files                         ",
 100    " usage:                                                                   ",
 101    " DenseMultiFramesToDicom dirin=rootDirectoryName                          ",
 102    "              [studyUID = ] [patName = ]                                  ",
 103    "              [listonly] [verbose] [debug]                                ",
 104    "                                                                          ",
 105    " studyUID   : *aware* user wants to add the serie                         ",
 106    "                                             to an already existing study ",
 107    " verbose  : user wants to run the program in 'verbose mode'               ",
 108    " debug    : *developer*  wants to run the program in 'debug mode'         ",
 109    FINISH_USAGE
 110
 111    // ----- Initialize Arguments Manager ------
 112
 113    GDCM_NAME_SPACE::ArgMgr *am = new GDCM_NAME_SPACE::ArgMgr(argc, argv);
 114
 115    if (argc == 1 || am->ArgMgrDefined("usage"))
 116    {
 117       am->ArgMgrUsage(usage); // Display 'usage'
 118       delete am;
 119       return 0;
 120    }
 121
 122    const char *dirNamein;
 123    dirNamein  = am->ArgMgrGetString("dirin",".");
 124
 125    if (am->ArgMgrDefined("debug"))
 126       GDCM_NAME_SPACE::Debug::DebugOn();
 127
 128    int verbose  = am->ArgMgrDefined("verbose");
 129    int listonly = am->ArgMgrDefined("listonly");
 130    std::string patName = am->ArgMgrGetString("patname", dirNamein);
 131
 132    bool userDefinedStudy = am->ArgMgrDefined("studyUID");
 133    const char *studyUID;
 134    if (userDefinedStudy)
 135       studyUID  = am->ArgMgrGetString("studyUID");
 136
 137    // not described *on purpose* in the Usage !
 138
 139    bool userDefinedSerie = am->ArgMgrDefined("serieUID");
 140    const char *serieUID;
 141    if(userDefinedSerie)
 142       serieUID = am->ArgMgrGetString("serieUID");
 143
 144    // if unused Param we give up
 145    if ( am->ArgMgrPrintUnusedLabels() )
 146    {
 147       am->ArgMgrUsage(usage);
 148       delete am;
 149       return 0;
 150    }
 151    delete am;  // we don't need Argument Manager any longer
 152
 153    // ----- Begin Processing -----
 154
 155    if ( ! GDCM_NAME_SPACE::DirList::IsDirectory(dirNamein) )
 156    {
 157       std::cout << "KO : [" << dirNamein << "] is not a Directory."
 158                 << std::endl;
 159       return 0;
 160    }
 161    else
 162    {
 163       if (verbose)
 164          std::cout << "OK : [" << dirNamein << "] is a Directory." << std::endl;
 165    }
 166
 167    std::string strDirNamein(dirNamein);
 168    GDCM_NAME_SPACE::DirList dirList(strDirNamein, true); // (recursively) the list of files
 169
 170    if (listonly)
 171    {
 172       std::cout << "------------List of found files ------------" << std::endl;
 173       dirList.Print();
 174       std::cout << std::endl;
 175       return 0;
 176     }
 177
 178    std::string filenameout;
 179
 180
 181
 182
 183    std::string strStudyUID;
 184    std::string strSerieUID;
 185
 186    if (userDefinedStudy)
 187       strSerieUID =  studyUID;
 188    else
 189       strStudyUID =  GDCM_NAME_SPACE::Util::CreateUniqueUID();
 190
 191    if (userDefinedStudy)
 192      strSerieUID =  serieUID;
 193    else
 194       strStudyUID =  GDCM_NAME_SPACE::Util::CreateUniqueUID();
 195
 196
 197    int serieNumber =0;
 198    GDCM_NAME_SPACE::DirListType fileNames;
 199    fileNames = dirList.GetFilenames();
 200    for (GDCM_NAME_SPACE::DirListType::iterator it = fileNames.begin();
 201                                     it != fileNames.end();
 202                                   ++it)
 203    {
 204       if ( GDCM_NAME_SPACE::Util::GetName((*it)).c_str()[0] == '.' )
 205       {
 206       // skip hidden files
 207          continue;
 208       }
 209       int sz = (*it).size();
 210       if ( (*it).c_str()[sz-1] != 't')
 211       {
 212          // skip non .txt files
 213          continue;
 214       }
 215       std::ifstream from( (*it).c_str() );
 216       if ( !from )
 217       {
 218          std::cout << "Can't open file" << *it << std::endl;
 219          //return 0;
 220       }
 221       else
 222       {
 223          if (verbose)
 224            std::cout << "Success in open file" << *it << std::endl;
 225          Load(from, *it, patName, strStudyUID, strSerieUID, serieNumber, verbose);
 226          serieNumber+=2;
 227          //return 0;
 228       }
 229    }
 230 }
 231
 232
 233 void Load(std::ifstream &from, std::string imageName, std::string patName,
 234           std::string strStudyUID, std::string strSerieUID, int serieNumber, int verbose)
 235 {
 236    if (verbose)
 237       std::cout << " ========= Deal with file [" << imageName << "]" << std::endl;
 238    if (!from)
 239       return;
 240    if (verbose)
 241      std::cout << " ========= Create Parametric images" << std::endl;
 242 /* was OK for single frames
 243 eg :
 244 ---------------------------
 245  Array dimensions = 58 x 56
 246 The following is the 2D array of peak circumferential strain map, all zero value
 247 pixels are outside the mask
 248      0.000000     0.000000     0.000000     0.000000     0.000000     0.000000
 249 -----------------------------
 250    std::string str1;
 251    int nx, ny;
 252    // Get nx, ny
 253    from >> str1;
 254    from >> str1;
 255    from >> str1;
 256    from >> nx;
 257    from >> str1;
 258    from >> ny;
 259
 260    std::cout << "nx " << nx << " ny " << ny << std::endl;
 261
 262    // Skip 2 lines.
 263    std::getline(from, str1);
 264    std::cout << "[" << str1 << "]" << std::endl;
 265    std::getline(from, str1);
 266    std::cout << "[" << str1 << "]" << std::endl;
 267    std::getline(from, str1);
 268    std::cout << "[" << str1 << "]" << std::endl;
 269  */
 270
 271  /* now, we deal with multiframes
 272  eg :
 273  ------------------------------------------
 274 X dim, Y dim, N of frames = 52x59x14
 275 Temporal resolution = 31.9200 ms
 276 First frame starts at 47.9600 ms
 277 The following are the 3D arrays of peak circumferential strain map and the magnitude images,
 278 All pixels with zero strain values are outside the masks.
 279      0.000000     0.000000     0.000000     0.000000     0.000000     0.000000
 280 --------------------------------------------
 281 */
 282    std::string str1;
 283    int nx, ny, nf;
 284     from >> str1; // X dim,
 285     from >> str1;
 286     from >> str1; // Y dim,
 287     from >> str1;
 288     from >> str1; // N of frames =
 289     from >> str1;
 290     from >> str1;
 291     from >> str1;
 292     from >> str1; // 52x59x14
 293
 294     if(verbose)
 295        std::cout << "[" << str1 << "]" << std::endl;
 296
 297     sscanf( str1.c_str(),"%dx%dx%d", &nx,&ny,&nf);
 298     std::cout << nx << " " << ny << " " << nf << std::endl;
 299
 300     std::getline(from, str1);
 301
 302     from >> str1; // Temporal
 303     from >> str1; // Resolution
 304     from >> str1; // =
 305
 306     from >> str1;
 307
 308     float temporalResolution;
 309     sscanf( str1.c_str(),"%f",&temporalResolution);
 310     if(verbose)
 311       std::cout << "temporal Resolution = " << temporalResolution << std::endl;
 312     std::getline(from, str1);
 313
 314     from >> str1; // First
 315     from >> str1; // frame
 316     from >> str1; // starts
 317     from >> str1; // at
 318
 319     from >> str1;
 320     float timeStart;
 321     sscanf( str1.c_str(),"%f",&timeStart);
 322     std::cout << "time Start = " << timeStart << std::endl;
 323     std::getline(from, str1);
 324
 325        // Skip 3 lines.
 326     for (int k=0; k<2; k++) //
 327     {
 328        std::getline(from, str1);
 329        std::cout << str1 << std::endl;
 330     }
 331
 332   //float *f = new float(nx*ny);
 333   // -->float *f = new float[nx*ny]; // Would be better !
 334    float *f = (float *) malloc(nx*ny*nf*sizeof(float));
 335   // float mini = FLT_MAX, maxi = FLT_MIN;
 336    float val;
 337
 338    int imageNumber = 0;
 339    float currentTime;
 340    currentTime = timeStart;
 341    int l1;
 342    for (l1=0; l1<nf; l1++)  // Loop on the frames
 343    {
 344      //std::cout << "Frame nb " << l1 << std::endl;
 345      for( int j=0; j<ny; j++)
 346      {
 347       for (int i=0; i<nx; i++)
 348       {
 349          str1="";
 350          //eatwhite(from);
 351          char c;
 352          for (;;) //eatwhite(from);
 353          {
 354             if (!from.get(c))
 355             {
 356                std::cout << " !from.get(c) ";
 357                break;
 358             }
 359             if (!isspace(c) )
 360             {
 361                //std::cout << " !isspace(c) ";
 362                from.putback(c);
 363                break;
 364             }
 365          } //end eatwhite(from);
 366
 367         // trouble : when space is missing "-0.0990263-8.8778"
 368         // is not interpreted as TWO values  :-(
 369         // from >> str1;
 370
 371          int first = 1;
 372          char previous = 'z';
 373          for(;;)
 374          {
 375             from.get(c);
 376             if ( c == ' ')
 377                break;
 378             if ( first != 1 && c == '-' && previous != 'e')
 379             {
 380                from.putback(c);
 381                //std::cout << " One more gauffre in frame:" << std::dec << l
 382                //         << ", line : " << j << " element " << i << std::endl;
 383                break;
 384              }
 385
 386              first=0;
 387              previous = c;
 388              str1=str1+c;
 389          }
 390
 391          val = (float)atof(str1.c_str());
 392          //std::cout << "  " << val;
 393          *(f+ /*l*nx*ny + */j*nx+i) = val;
 394
 395         if(from.eof())
 396         {
 397             std::cout << "Missing values at [" << std::dec << j <<"," << i << "]"
 398                       << std::endl;
 399            break;
 400          }
 401       }
 402
 403       //std::cout << std::endl;
 404       //std::cout << std::endl << " line nb : " << j
 405       //          << " line length : " << l << std::endl;
 406     }
 407
 408    // std::cout << "mini : "<< mini  << " maxi : " << maxi << std::endl;
 409 /*
 410 // values are expressed as % as a fraction, actually!)
 411 // It's useless to rescale them as uint16_t : just *100
 412     uint16_t *img = new uint16_t[ny*nx];
 413     uint16_t *ptr = img;
 414     float *tmp = f;
 415     float div = maxi-mini;
 416     std::cout << "div : " << div << " mini : " << mini << std::endl;
 417     for( int k=0; k<ny*nx; k++)
 418     {
 419        *ptr = (uint16_t)((*tmp * 65535.0) / div);
 420        tmp ++;
 421        ptr++;
 422     }
 423 */
 424     int16_t *img = new int16_t[ny*nx];
 425     int16_t *ptr = img;
 426     float *tmp = f /* + l*ny*nx */ ; // start on image nb l.
 427     for( int k=0; k<ny*nx; k++)
 428     {
 429        if(*tmp > 1.0) // artefacted pixel
 430           *ptr = 0;
 431        else        /// \todo FIXME : what about max threshold ?
 432         *ptr = (int16_t)(*tmp *100);
 433
 434       // std::cout << std::dec << "[" << *tmp <<" : " << *ptr << "] ";
 435        tmp ++;
 436        ptr++;
 437     }
 438
 439  // GDCM_NAME_SPACE::Debug::DebugOn();
 440
 441         std::ostringstream str;
 442         GDCM_NAME_SPACE::File *file;
 443         file = GDCM_NAME_SPACE::File::New();
 444
 445   // Set the image size
 446         str.str("");
 447         str << nx;
 448         file->InsertEntryString(str.str(),0x0028,0x0011,"US"); // Columns
 449         str.str("");
 450         str << ny;
 451         file->InsertEntryString(str.str(),0x0028,0x0010,"US"); // Rows
 452
 453   // Set the pixel type
 454   //      16; //8, 16, 32
 455         file->InsertEntryString("16",0x0028,0x0100,"US"); // Bits Allocated
 456         str.str("");
 457         str << 16; // may be 12 or 16 if componentSize =16
 458         file->InsertEntryString("16",0x0028,0x0101,"US"); // Bits Stored
 459
 460         file->InsertEntryString("15",0x0028,0x0102,"US"); // High Bit
 461
 462   // Set the pixel representation // 0/1 1 : signed
 463         file->InsertEntryString("1",0x0028,0x0103, "US"); // Pixel Representation
 464
 465   // Set the samples per pixel // 1:Grey level, 3:RGB
 466         file->InsertEntryString("1",0x0028,0x0002, "US"); // Samples per Pixel
 467
 468 /*
 469   // Set Rescale Intercept
 470         str.str("");
 471         str << div;
 472         file->InsertEntryString(str.str(),0x0028,0x1052,"DS");
 473
 474   // Set Rescale Slope
 475         str.str("");
 476         str << mini;
 477         file->InsertEntryString(str.str(),0x0028,0x1053,"DS");
 478 */
 479
 480 // 0020 0037 DS 6 Image Orientation (Patient)
 481          file->InsertEntryString("1.0\\0.0\\0.0\\0.0\\1.0\\0.0",0x0020,0x0037, "DS"); //[1\0\0\0\1\0] : Axial
 482
 483 // 0020 0032 DS 3 Image Position (Patient)
 484         char charImagePosition[256];
 485         sprintf(charImagePosition,"0.0\\0.0\\%f",float(l1%nf));
 486
 487 // 0020 0x1041 DS 1 Slice Location
 488         sprintf(charImagePosition,"%f",float(l1%nf));
 489         file->InsertEntryString(charImagePosition,0x0020,0x1041, "DS");
 490
 491 //0008 103e LO 1 Series Description
 492         file->InsertEntryString(imageName,0x0008,0x103e, "LO");
 493
 494         file->InsertEntryString(strStudyUID,0x0020,0x000d,"UI");
 495         file->InsertEntryString(strSerieUID,0x0020,0x000e,"UI");
 496         file->InsertEntryString(patName,0x0010,0x0010, "PN");   // Patient's Name
 497
 498 //0020 0011 "IS" Series Number
 499          sprintf(charImagePosition,"%04d",serieNumber);
 500          file->InsertEntryString(charImagePosition,0x0020,0x0011, "IS");
 501
 502 //0020 0011 "IS" Instance Number
 503          sprintf(charImagePosition,"%04d",imageNumber);
 504          file->InsertEntryString(charImagePosition,0x0020,0x0013, "IS");
 505
 506 //0018 1060 "DS" Time Trigger
 507          sprintf(charImagePosition,"%f",currentTime);
 508          file->InsertEntryString(charImagePosition,0x0018,0x1060, "DS");
 509
 510    // file->Print();
 511
 512     GDCM_NAME_SPACE::FileHelper *fh;
 513     fh = GDCM_NAME_SPACE::FileHelper::New(file);
 514     // cast is just to avoid warnings (*no* conversion)
 515     fh->SetImageData((uint8_t *)img,nx*ny*sizeof(uint16_t));
 516     fh->SetWriteModeToRaw();
 517     fh->SetWriteTypeToDcmExplVR();
 518
 519     fh->SetWriteTypeToDcmExplVR();
 520
 521     char numero[10];
 522     sprintf(numero, "%02d", l1);
 523     std::string fileName = imageName + "." + numero + ".dcm";
 524     if(verbose)
 525       std::cout << "fileName " << fileName << std::endl;
 526
 527     if( !fh->Write(fileName))
 528        std::cout << "Failed for [" << fileName << "]\n"
 529                  << "           File is unwrittable" << std::endl;
 530
 531     delete img;
 532     currentTime += temporalResolution;
 533     imageNumber ++;
 534
 535   } // end loop on frames
 536
 537
 538    // Anatomical Images.
 539   std::cout << " ========= Create Anatomical images" << std::endl;
 540
 541   strSerieUID =  GDCM_NAME_SPACE::Util::CreateUniqueUID();
 542   imageNumber = 0;
 543   currentTime = timeStart;
 544
 545   for (int fr=0; fr<nf; fr++)  // Loop on the frames
 546   {
 547    //std::cout << "Frame nb " << l << std::endl;
 548      for( int j=0; j<ny; j++)
 549      {
 550       int l =0;
 551       for (int i=0; i<nx; i++)
 552       {
 553          //eatwhite(from);
 554          char c;
 555          for (;;)
 556          {
 557             if (!from.get(c))
 558                break;
 559             if (!isspace(c))
 560             {
 561                from.putback(c);
 562                break;
 563             }
 564          }
 565          from >> str1;
 566          val = (float)atof(str1.c_str());
 567         // std::cout << "  " << val;
 568          *(f+ /*l*nx*ny + */j*nx+i) = val;
 569
 570         if(from.eof())
 571         {
 572             std::cout << "Missing values at [" << std::dec << j <<"," << i << "]"
 573                       << std::endl;
 574            break;
 575          }
 576          l++;
 577       }
 578      }
 579
 580     uint16_t *img = new uint16_t[ny*nx];
 581     uint16_t *ptr = img;
 582     float *tmp = f /* + l*ny*nx */ ; // start on image nb l.
 583     for( int k=0; k<ny*nx; k++)
 584     {
 585        *ptr = (int16_t)(*tmp);
 586        tmp ++;
 587        ptr++;
 588     }
 589         std::ostringstream str;
 590         GDCM_NAME_SPACE::File *file;
 591         file = GDCM_NAME_SPACE::File::New();
 592
 593   // Set the image size
 594         str.str("");
 595         str << nx;
 596         file->InsertEntryString(str.str(),0x0028,0x0011,"US"); // Columns
 597         str.str("");
 598         str << ny;
 599         file->InsertEntryString(str.str(),0x0028,0x0010,"US"); // Rows
 600
 601   // Set the pixel type
 602   //      16; //8, 16, 32
 603         file->InsertEntryString("16",0x0028,0x0100,"US"); // Bits Allocated
 604         str.str("");
 605         str << 16; // may be 12 or 16 if componentSize =16
 606         file->InsertEntryString("16",0x0028,0x0101,"US"); // Bits Stored
 607
 608         file->InsertEntryString("15",0x0028,0x0102,"US"); // High Bit
 609
 610   // Set the pixel representation // 0/1 1 : signed
 611         file->InsertEntryString("0",0x0028,0x0103, "US"); // Pixel Representation
 612
 613   // Set the samples per pixel // 1:Grey level, 3:RGB
 614         file->InsertEntryString("1",0x0028,0x0002, "US"); // Samples per Pixel
 615
 616 /*
 617   // Set Rescale Intercept
 618         str.str("");
 619         str << div;
 620         file->InsertEntryString(str.str(),0x0028,0x1052,"DS");
 621
 622   // Set Rescale Slope
 623         str.str("");
 624         str << mini;
 625         file->InsertEntryString(str.str(),0x0028,0x1053,"DS");
 626 */
 627
 628 // 0020 0037 DS 6 Image Orientation (Patient)
 629          file->InsertEntryString("1.0\\0.0\\0.0\\0.0\\1.0\\0.0",0x0020,0x0037, "DS"); //[1\0\0\0\1\0] : Axial
 630
 631 // 0020 0032 DS 3 Image Position (Patient)
 632         char charImagePosition[256];
 633         sprintf(charImagePosition,"0.0\\0.0\\%f",float(l1%nf));
 634
 635 // 0020 0x1041 DS 1 Slice Location
 636         sprintf(charImagePosition,"%f",float(l1%nf));
 637         file->InsertEntryString(charImagePosition,0x0020,0x1041, "DS");
 638
 639 //0008 103e LO 1 Series Description
 640         file->InsertEntryString(imageName,0x0008,0x103e, "LO");
 641
 642         file->InsertEntryString(strStudyUID,0x0020,0x000d,"UI");
 643         file->InsertEntryString(strSerieUID,0x0020,0x000e,"UI");
 644         file->InsertEntryString(patName,0x0010,0x0010, "PN");   // Patient's Name
 645
 646 //0020 0011 "IS" Series Number
 647          sprintf(charImagePosition,"%04d",serieNumber+1);
 648          file->InsertEntryString(charImagePosition,0x0020,0x0011, "IS");
 649
 650 //0020 0011 "IS" Instance Number
 651          sprintf(charImagePosition,"%04d",imageNumber);
 652          file->InsertEntryString(charImagePosition,0x0020,0x0013, "IS");
 653
 654 //0018 1060 "DS" Time Trigger
 655          sprintf(charImagePosition,"%f",currentTime);
 656          file->InsertEntryString(charImagePosition,0x0018,0x1060, "DS");
 657    // file->Print();
 658
 659     GDCM_NAME_SPACE::FileHelper *fh;
 660     fh = GDCM_NAME_SPACE::FileHelper::New(file);
 661     // cast is just to avoid warnings (*no* conversion)
 662     fh->SetImageData((uint8_t *)img,nx*ny*sizeof(uint16_t));
 663     fh->SetWriteModeToRaw();
 664     fh->SetWriteTypeToDcmExplVR();
 665
 666     fh->SetWriteTypeToDcmExplVR();
 667
 668     char numero[10];
 669     sprintf(numero, "%02d", l1);
 670     std::string fileName = imageName + ".Anatomical." + numero + ".dcm";
 671     std::cout << "fileName " << fileName << std::endl;
 672
 673     if( !fh->Write(fileName))
 674        std::cout << "Failed for [" << fileName << "]\n"
 675                  << "           File is unwrittable" << std::endl;
 676
 677     delete img;
 678     currentTime += temporalResolution;
 679     imageNumber ++;
 680
 681   } // end loop on frames
 682
 683    from.close();
 684 } // end void Load(