Example/PcpdenseToDicom.cxx

   1 /*=========================================================================
   2
   3   Program:   gdcm
   4   Module:    $RCSfile: PcpdenseToDicom.cxx,v $
   5   Language:  C++
   6   Date:      $Date: 2011/09/06 16:08:07 $
   7   Version:   $Revision: 1.9 $
   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 #include "gdcmUtil.h"
  32 #include "gdcmArgMgr.h"
  33
  34 /**
  35   * \brief
  36   *  Converts the "pcpdense" ".txt" (2008 version)  files into 16 bits Dicom Files,
  37   * Hope they don't change soon!
  38   */
  39
  40 void MakeDicomImage(unsigned short int *tabVal, int X, int Y, std::string dcmImageName, const char * patientName, int nbFrames,
  41                     std::string studyUID, std::string serieUID, std::string SerieDescr, int imgNum, bool m );
  42
  43 void LoadImage(std::ifstream &from,  unsigned short int * );
  44
  45 void LoadImageX2(std::ifstream &from, unsigned short int * );
  46 void RotateImage(unsigned short int *image, unsigned short int *image2, int NX, int NY);
  47 void FlipImage  (unsigned short int *image, unsigned short int *image2, int NX, int NY);
  48
  49 bool verbose;
  50
  51 int main(int argc, char *argv[])
  52 {
  53    START_USAGE(usage)
  54    " \n pcpdenseToDicom :\n                                                   ",
  55    "        Converts the '.txt' files into 16 bits Dicom-like Files,          ",
  56    " usage:                                                                   ",
  57    " pcpdenseToDicom rootfilename=...                                         ",
  58    "                 (e.g.. :   meas_MID380_DENSE_stacked_slices_aif_FID81637)",
  59    "                 numberOfSlices =  (default : 3)                          ",
  60    "                 X2 : multiply x 2 image size                             ",
  61    "                 m :create multiframe files instead of image stacks       ",
  62    "                 [patientname = Patient's name]                           ",
  63    "                 [verbose] [debug]                                        ",
  64    "                                                                          ",
  65    " verbose  : user wants to run the program in 'verbose mode'               ",
  66    " debug    : *developer*  wants to run the program in 'debug mode'         ",
  67    FINISH_USAGE
  68
  69    // ----- Initialize Arguments Manager ------
  70
  71    GDCM_NAME_SPACE::ArgMgr *am = new GDCM_NAME_SPACE::ArgMgr(argc, argv);
  72
  73    if (argc == 1 || am->ArgMgrDefined("usage"))
  74    {
  75       am->ArgMgrUsage(usage); // Display 'usage'
  76       delete am;
  77       return 0;
  78    }
  79    // Seems that ArgMgrWantString doesn't work on MacOS
  80    if(!am->ArgMgrDefined("rootfilename"))
  81    {
  82       std::cout << "'rootfilename' is mandatory" << std::endl;
  83       exit(0);
  84    }
  85
  86    const char *rootfilename = am->ArgMgrWantString("rootfilename",usage);
  87    int numberOfSlices       = am->ArgMgrGetInt("numberOfSlices",3);
  88    const char *patientName  = am->ArgMgrGetString("patientname", "Patient^Name");
  89
  90    if (am->ArgMgrDefined("debug"))
  91       GDCM_NAME_SPACE::Debug::DebugOn();
  92
  93    verbose         = ( 0 != am->ArgMgrDefined("verbose") );
  94    bool X2         = ( 0 != am->ArgMgrDefined("X2") );
  95    bool multiframe = ( 0 != am->ArgMgrDefined("m") );
  96
  97    // if unused Param we give up
  98    if ( am->ArgMgrPrintUnusedLabels() )
  99    {
 100       am->ArgMgrUsage(usage);
 101       delete am;
 102       return 0;
 103    }
 104    delete am;  // we don't need Argument Manager any longer
 105
 106    // ----- Begin Processing -----
 107
 108    std::ostringstream Ecc;
 109    std::ostringstream perf;
 110    std::ostringstream WashoutTc;
 111
 112    std::string strSerieUID;
 113    std::string strStudyUID =  GDCM_NAME_SPACE::Util::CreateUniqueUID();
 114    std::string dcmImageName, textFileName, patientname, serieDescr;
 115    std::string deb(rootfilename);
 116
 117    unsigned short int *image;
 118    unsigned short int *image2;
 119
 120    int NX, NY;
 121
 122   // Get some info
 123   // -------------
 124
 125    {
 126         Ecc.str(rootfilename);
 127         Ecc   << Ecc.str() << "_s0" << "_Ecc.txt";
 128
 129         std::ifstream fromEcc( Ecc.str().c_str() );
 130         if ( !fromEcc )
 131         {
 132            std::cout << "Can't open file [" << Ecc.str() << "]" << std::endl;
 133            exit(0);
 134         }
 135         std::string str1;
 136
 137          fromEcc >> str1;
 138          fromEcc >> str1;
 139
 140          fromEcc >> NX;
 141          fromEcc >> NY;
 142          std::cout << "NX, NY : " << NX << ", " << NY << std::endl;
 143    }
 144
 145    int mult;
 146    if (X2)
 147       mult=4;
 148    else
 149       mult=1;
 150
 151    if (multiframe) {
 152       image  = new unsigned short int[NX*NY*mult*numberOfSlices];
 153       image2 = new unsigned short int[NX*NY*mult*numberOfSlices];
 154    } else {
 155       image  = new unsigned short int[NX*NY*mult];
 156       image2 = new unsigned short int[NX*NY*mult];
 157    }
 158
 159    // === Ecc ===
 160
 161    strSerieUID =  GDCM_NAME_SPACE::Util::CreateUniqueUID();
 162
 163    serieDescr = "Ecc";
 164
 165    if (!multiframe) {
 166      for (int i=0; i<numberOfSlices; i++)
 167      {
 168         Ecc.str(rootfilename);
 169         Ecc << Ecc.str() << "_s" << i << "_Ecc.txt";
 170
 171         std::ifstream fromEcc( Ecc.str().c_str() );
 172         if ( !fromEcc )
 173         {
 174            std::cout << "Can't open file [" << Ecc.str() << "]" << std::endl;
 175            exit(0);
 176         }
 177         std::cout << "Open file [" << Ecc.str() << "] : OK" << std::endl;
 178         dcmImageName = Ecc.str() + ".dcm";
 179
 180         if (X2)
 181         {
 182            LoadImageX2(fromEcc, image);
 183            RotateImage(image, image2, NX*2, NY*2);
 184            FlipImage(image2, image2, NY*2, NX*2);
 185            MakeDicomImage(image, NY*2, NX*2, dcmImageName, patientName, 1, strStudyUID, strSerieUID, serieDescr, i, multiframe );
 186         }
 187         else
 188         {
 189            LoadImage(fromEcc, image);
 190            RotateImage(image, image2, NX, NY);
 191            FlipImage(image2, image2, NY, NX);
 192            MakeDicomImage(image2, NY, NX,     dcmImageName, patientName, 1, strStudyUID, strSerieUID, serieDescr, i, multiframe );
 193         }
 194         fromEcc.close();
 195
 196      } // end : for (int i=0; i<numberOfSlices
 197    }
 198
 199    if (multiframe) {
 200      for (int i=0; i<numberOfSlices; i++)
 201      {
 202         Ecc.str(rootfilename);
 203         Ecc   << Ecc.str()    << "_s" << i << "_Ecc.txt";
 204
 205         std::ifstream fromEcc( Ecc.str().c_str() );
 206         if ( !fromEcc )
 207         {
 208            std::cout << "Can't open file [" << Ecc.str() << "]" << std::endl;
 209            exit(0);
 210         }
 211
 212         std::cout << "Open file [" << Ecc.str() << "] : OK" << std::endl;
 213         if (X2)
 214         {
 215           LoadImageX2(fromEcc ,&image[NX*NY*4*i] );
 216           RotateImage(&image[NX*NY*4*i], &image2[NX*NY*4*i] , NX*2, NY*2);
 217           FlipImage(&image2[NX*NY*4*i],  &image2[NX*NY*4*i] , NY*2, NX*2);
 218         }
 219         else
 220         {
 221           LoadImage(fromEcc, &image[NX*NY*i] );
 222           RotateImage(&image[NX*NY*i],&image2[NX*NY*i], NX, NY);
 223           FlipImage(&image2[NX*NY*i], &image2[NX*NY*i], NY, NX);
 224         }
 225
 226         fromEcc.close();
 227      } // end : for (int i=0; i<numberOfSlices
 228
 229      dcmImageName = deb + "_Ecc.dcm";
 230      if (X2)
 231         MakeDicomImage(image2, NY*2, NX*2, dcmImageName, patientName, numberOfSlices, strStudyUID, strSerieUID, serieDescr, 0, multiframe );
 232      else
 233         MakeDicomImage(image2, NY, NX,     dcmImageName, patientName, numberOfSlices, strStudyUID, strSerieUID, serieDescr, 0, multiframe );
 234    }  // end : if (multiframe)
 235
 236
 237
 238 // === perf ===
 239
 240    strSerieUID =  GDCM_NAME_SPACE::Util::CreateUniqueUID();
 241
 242    serieDescr = "perf";
 243
 244    if (!multiframe) {
 245      for (int i=0; i<numberOfSlices; i++)
 246      {
 247         perf.str(rootfilename);
 248         perf  << perf.str()    << "_s" << i << "_perf.txt";
 249
 250         std::ifstream fromperf( perf.str().c_str() );
 251         if ( !fromperf )
 252         {
 253            std::cout << "Can't open file [" << perf.str() << "]" << std::endl;
 254            exit(0);
 255         }
 256         std::cout << "Open file [" << perf.str() << "] : OK" << std::endl;
 257
 258         if (X2)
 259           LoadImageX2(fromperf ,image );
 260         else
 261           LoadImage(fromperf, image );
 262
 263         fromperf.close();
 264
 265         dcmImageName = perf.str() + ".dcm";
 266
 267         if (X2)
 268            MakeDicomImage(image, NX*2, NY*2, dcmImageName, patientName, 1, strStudyUID, strSerieUID, serieDescr, i, multiframe );
 269         else
 270            MakeDicomImage(image, NX, NY,     dcmImageName, patientName, 1, strStudyUID, strSerieUID, serieDescr, i, multiframe );
 271
 272      } // end : for (int i=0; i<numberOfSlices
 273    }
 274
 275
 276    if (multiframe) {
 277      for (int i=0; i<numberOfSlices; i++)
 278      {
 279         perf.str(rootfilename);
 280         perf  << perf.str()    << "_s" << i << "_perf.txt";
 281
 282         std::ifstream fromperf( perf.str().c_str() );
 283         if ( !fromperf )
 284         {
 285            std::cout << "Can't open file [" << perf.str() << "]" << std::endl;
 286            exit(0);
 287         }
 288         std::cout << "Open file [" << perf.str() << "] : OK" << std::endl;
 289
 290       if (X2)
 291          LoadImageX2(fromperf, &image[NX*NY*i*4]);
 292       else
 293          LoadImage(fromperf, &image[NX*NY*i] );
 294
 295       fromperf.close();
 296    } // end : for (int i=0; i<numberOfSlices
 297
 298      dcmImageName = deb + "_perf.dcm";
 299      if (X2)
 300         MakeDicomImage(image, NX*2, NY*2, dcmImageName, patientName, numberOfSlices, strStudyUID, strSerieUID, serieDescr, 0, multiframe );
 301      else
 302         MakeDicomImage(image, NX, NY,     dcmImageName, patientName, numberOfSlices, strStudyUID, strSerieUID, serieDescr, 0, multiframe );
 303  }
 304
 305
 306
 307 // === WashoutTc ===
 308
 309
 310    strSerieUID =  GDCM_NAME_SPACE::Util::CreateUniqueUID();
 311
 312    serieDescr = "WashoutTc";
 313
 314    if (!multiframe) {
 315
 316    for (int i1=0; i1<numberOfSlices; i1++)
 317    {
 318       WashoutTc.str(rootfilename);
 319       WashoutTc <<  WashoutTc.str() << "_s" << i1 << "_WashoutTc.txt";
 320
 321       std::ifstream fromWashoutTc( WashoutTc.str().c_str() );
 322       if ( !fromWashoutTc )
 323       {
 324          std::cout << "Can't open file [" << WashoutTc.str() << "]" << std::endl;
 325          exit(0);
 326       }
 327
 328       std::cout << "Open file [" << WashoutTc.str() << "] : OK" << std::endl;
 329       if (X2)
 330          LoadImageX2(fromWashoutTc,image );
 331        else
 332          LoadImage(fromWashoutTc,image );
 333       fromWashoutTc.close();
 334
 335       dcmImageName = WashoutTc.str() + ".dcm";
 336       if (X2)
 337          MakeDicomImage(image, NX*2, NY*2, dcmImageName, patientName, 1, strStudyUID, strSerieUID, serieDescr, i1, multiframe );
 338       else
 339          MakeDicomImage(image, NX, NY,     dcmImageName, patientName, 1, strStudyUID, strSerieUID, serieDescr, i1, multiframe );
 340     } // end : for (int i=0; i<numberOfSlices
 341   }
 342
 343    if (multiframe) {
 344      for (int i=0; i<numberOfSlices; i++)
 345      {
 346         WashoutTc.str(rootfilename);
 347         WashoutTc <<  WashoutTc.str() << "_s" << i << "_WashoutTc.txt";
 348
 349         std::ifstream fromWashoutTc( WashoutTc.str().c_str() );
 350         if ( !fromWashoutTc )
 351         {
 352            std::cout << "Can't open file [" << WashoutTc.str() << "]" << std::endl;
 353            exit(0);
 354         }
 355         std::cout << "Open file [" << WashoutTc.str() << "] : OK" << std::endl;
 356
 357         if (X2)
 358           LoadImageX2(fromWashoutTc,  &image[NX*NY*4*i]);
 359         else
 360           LoadImage(fromWashoutTc, &image[NX*NY*i] );
 361         fromWashoutTc.close();
 362
 363      } // end : for (int i=0; i<numberOfSlices
 364
 365      dcmImageName = deb + "_WashoutTc.dcm";
 366      if (X2)
 367         MakeDicomImage(image, NX*2, NY*2, dcmImageName, patientName, numberOfSlices, strStudyUID, strSerieUID, serieDescr, 0, multiframe );
 368      else
 369         MakeDicomImage(image, NX, NY,     dcmImageName, patientName, numberOfSlices, strStudyUID, strSerieUID, serieDescr, 0, multiframe );
 370    }
 371
 372    delete []image;
 373    return 1;
 374 }
 375
 376 // =====================================================================================================================
 377
 378 void LoadImage(std::ifstream &from,  unsigned short int *image)
 379 {
 380 // in any file ".txt" :
 381
 382 /*
 383 XY Dimensions           47          50
 384      0.000000     0.000000     0.000000     0.000000     0.000000     0.000000
 385      ...
 386 */
 387
 388    if (!from)
 389       return;
 390
 391    std::string str1;
 392
 393    from >> str1;
 394    from >> str1;
 395
 396    int NX, NY;
 397    from >> NX;
 398    from >> NY;
 399    std::cout << "NX, NY : " << NX << ", " << NY << std::endl;
 400
 401    float pixelValue;
 402
 403     int i, j;
 404     for( i=0;i<NY;i++) {
 405         for(j=0;j<NX;j++) {
 406            from >> pixelValue;
 407            image[i*NX+j] = (unsigned short int)(pixelValue * 1000.); // Why do we multiply by 1000? // JPR
 408         }
 409      }
 410 }
 411
 412
 413 // =====================================================================================================================
 414
 415 void LoadImageX2(std::ifstream &from,  unsigned short int *image )
 416 {
 417 // in any file ".txt" :
 418
 419 /*
 420 XY Dimensions           47          50
 421      0.000000     0.000000     0.000000     0.000000     0.000000     0.000000
 422      ...
 423 */
 424
 425    if (!from)
 426       return;
 427
 428    std::string str1;
 429
 430     from >> str1;
 431     from >> str1;
 432
 433     int NX, NY;
 434     from >> NX;
 435     from >> NY;
 436     std::cout << "NX, NY : " << NX << ", " << NY << std::endl;
 437
 438      int k;
 439      for( k=0;k<NX*NY*4;k++)
 440           image[k] = 0;
 441
 442     float pixelValue;
 443
 444      int i, j;
 445  /*
 446      for( i=0;i<NY;i++) {
 447            for(j=0;j<NX;j++) {
 448               from >> pixelValue;
 449               pixelValue*=1000.;  // Why do we multiply by 1000? // JPR
 450               image[i*4*NX + j*2] = image[i*4*NX + j*2+1] =  image[(i*4+2)*NX + j*2] =  image[(i*4+2)*NX + j*2+1] = (unsigned short int)(pixelValue);
 451            }
 452      }
 453  */
 454
 455      int lgrLigneNvlleImage = NX+NX;
 456      int debLigneNvlleImage = 0;
 457      for( i=0;i<NY;i++) {
 458            for(j=0;j<NX;j++) {
 459               from >> pixelValue;
 460               pixelValue*=1000.;  // Why do we multiply by 1000? // JPR
 461               image[debLigneNvlleImage + j+j] =
 462               image[debLigneNvlleImage + j+j +1] =
 463               image[debLigneNvlleImage + lgrLigneNvlleImage +j+j] =
 464               image[debLigneNvlleImage + lgrLigneNvlleImage + j+j +1] =
 465                  (unsigned short int)(pixelValue);
 466            }
 467            debLigneNvlleImage += 2*lgrLigneNvlleImage;
 468      }
 469 }
 470 // =====================================================================================================================================
 471 void RotateImage(unsigned short int *image, unsigned short int *image2, int NX, int NY)
 472 {
 473      int k = 0;
 474      for( int i=0;i<NY;i++) {
 475            for(int j=0;j<NX;j++) {
 476               image2[NY*j + i] = image[k];
 477               k++;
 478            }
 479      }
 480 }
 481
 482
 483 // =====================================================================================================================================
 484 void FlipImage(unsigned short int *image, unsigned short int *image2, int NX, int NY)
 485 {
 486      unsigned short int temp;
 487      for(int i=0;i<NY/2;i++) {
 488            for(int j=0;j<NX;j++) {
 489               temp = image[NX*i + j];
 490               image2[NX*i + j] = image[NX*(NY-i-1) + j];
 491               image2[NX*(NY-i-1) + j] = temp;
 492            }
 493      }
 494 }
 495
 496
 497 // =====================================================================================================================================
 498 void MakeDicomImage(unsigned short int *tabVal, int X, int Y, std::string dcmImageName, const char * patientName, int nbFrames, std::string studyUID, std::string serieUID, std::string SerieDescr, int imgNum, bool m)
 499 {
 500
 501 std::cout << "========================> in MakeDicomImage : dcmImageName = [" << dcmImageName << "] NX= " << X << " NY= " << Y << std::endl;
 502  // GDCM_NAME_SPACE::Debug::DebugOn();
 503
 504    std::ostringstream str;
 505
 506    GDCM_NAME_SPACE::File *file;
 507    file = GDCM_NAME_SPACE::File::New();
 508
 509   // Set the image size
 510    str.str("");
 511    str << X;
 512    file->InsertEntryString(str.str(),0x0028,0x0011,"US"); // Columns
 513    str.str("");
 514    str << Y;
 515    file->InsertEntryString(str.str(),0x0028,0x0010,"US"); // Rows
 516
 517   // Set the pixel type
 518   //      16; //8, 16, 32
 519    file->InsertEntryString("16",0x0028,0x0100,"US"); // Bits Allocated
 520
 521    str.str("");
 522    str << 16; // may be 12 or 16 if componentSize =16
 523    file->InsertEntryString("16",0x0028,0x0101,"US"); // Bits Stored
 524    file->InsertEntryString("15",0x0028,0x0102,"US"); // High Bit
 525
 526   // Set the pixel representation // 0/1 , 0=unsigned
 527    file->InsertEntryString("1",0x0028,0x0103, "US"); // Pixel Representation
 528
 529   // Set the samples per pixel // 1:Grey level, 3:RGB
 530    file->InsertEntryString("1",0x0028,0x0002, "US"); // Samples per Pixel
 531
 532
 533    if (nbFrames != 1)
 534    {
 535       str.str("");
 536       str << nbFrames;
 537       file->InsertEntryString(str.str(),0x0028,0x0008,"IS"); // Number of Frames
 538    }
 539
 540    if (strlen(patientName) != 0)
 541       file->InsertEntryString(patientName,0x0010,0x0010, "PN"); // Patient's Name
 542
 543    file->InsertEntryString(studyUID, 0x0020, 0x000d, "UI");
 544    file->InsertEntryString(serieUID, 0x0020, 0x000e, "UI");
 545
 546    file->InsertEntryString(SerieDescr,0x0008,0x103e, "LO");  // Series Description
 547
 548 // 0020 0037 DS 6 Image Orientation (Patient)
 549    file->InsertEntryString("1.0\\0.0\\0.0\\0.0\\1.0\\0.0",0x0020,0x0037, "DS"); //[1\0\0\0\1\0] : Axial   (Tant pis!)
 550
 551 // 0020 0032 DS 3 Image Position (Patient)
 552    char charImagePosition[256];
 553    sprintf(charImagePosition,"0.0\\0.0\\%f",(float)imgNum);
 554    file->InsertEntryString(charImagePosition,0x0020,0x0032, "DS");  //0020 0032 DS 3 Image Position (Patient)
 555
 556 // 0020 0x1041 DS 1 Slice Location
 557         sprintf(charImagePosition,"%f",float(imgNum));
 558         file->InsertEntryString(charImagePosition,0x0020,0x1041, "DS");
 559 /*
 560   // Set Rescale Intercept
 561         str.str("");
 562         str << div;
 563         file->InsertEntryString(str.str(),0x0028,0x1052,"DS");
 564
 565   // Set Rescale Slope
 566         str.str("");
 567         str << mini;
 568         file->InsertEntryString(str.str(),0x0028,0x1053,"DS");
 569 */
 570
 571    GDCM_NAME_SPACE::FileHelper *fileH;
 572    fileH = GDCM_NAME_SPACE::FileHelper::New(file);
 573    // cast is just to avoid warnings (*no* conversion is performed)
 574    //fileH->SetImageData((uint8_t *)img,int(maxX*maxY)*sizeof(uint16_t)); // troubles when maxX, mayY are *actually* float!
 575
 576    fileH->SetImageData((uint8_t *)tabVal,X*Y*nbFrames*sizeof(uint16_t));
 577    fileH->SetWriteModeToRaw();
 578    fileH->SetWriteTypeToDcmExplVR();
 579
 580    if( !fileH->Write(dcmImageName))
 581       std::cout << "Failed for [" << dcmImageName << "]\n"
 582                 << "           File is unwrittable" << std::endl;
 583
 584   // file->Print();
 585
 586   // delete img;
 587    file->Delete();
 588    fileH->Delete();
 589 }
 590
 591
 592
 593 // =====================================================================================================================
 594