Example/RawToDicomStack.cxx

   1 /*=========================================================================
   2
   3   Program:   gdcm
   4   Module:    $RCSfile: RawToDicomStack.cxx,v $
   5   Language:  C++
   6   Date:      $Date: 2008/06/12 13:26:39 $
   7   Version:   $Revision: 1.1 $
   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 /**
  20  * Writes a Dicom file from a Raw File
  21  * User has to supply parameters.
  22  *
  23  */
  24 #include "gdcmFile.h"
  25 #include "gdcmFileHelper.h"
  26 #include "gdcmDebug.h"
  27 #include "gdcmUtil.h"
  28 #include "gdcmArgMgr.h"
  29
  30 #include <iostream>
  31 #include <sstream>
  32
  33 typedef char * PS8;
  34 typedef unsigned char * PU8;
  35 typedef short int * PS16;
  36 typedef unsigned short int * PU16;
  37 typedef int * PS32;
  38 typedef unsigned int * PU32;
  39
  40 #define CRR(t1,t2)     { for(int i=0;i<nY;i++)                   \
  41                            for(int j=0;j<nX;j++)                 \
  42                              for(int k=0;k<samplesPerPixel;k++)  \
  43                                  *((t2)planePixels + i*nX +j+k) = *(((t1)pixels)+ nbPlanes*nX*nY + i*nX +j+k);\
  44                        }
  45
  46 #define CFR(PPt)  switch ( pixelTypeOutCode ) {     \
  47           case -8   : CRR(PPt,PS8); break;          \
  48           case 8    : CRR(PPt,PU8); break;          \
  49           case -16  : CRR(PPt,PS16); break;         \
  50           case 16   : CRR(PPt,PU16); break;         \
  51           case -32  : CRR(PPt,PS32); break;         \
  52           case 32   : CRR(PPt,PU32); break;         \
  53           }
  54
  55
  56 void ConvertSwapZone(int pixelSize, void *Raw, size_t RawSize);
  57
  58 void ConvertSwapZone(int pixelSize, void *Raw, size_t RawSize)
  59 {
  60    unsigned int i;
  61    if ( pixelSize == 2 )
  62    {
  63       uint16_t *im16 = (uint16_t*)Raw;
  64       for( i = 0; i < RawSize / 2; i++ )
  65       {
  66          im16[i]= (im16[i] >> 8) | (im16[i] << 8 );
  67       }
  68    }
  69    else if ( pixelSize == 4 )
  70    {
  71       uint32_t s32;
  72       uint16_t high;
  73       uint16_t low;
  74       uint32_t *im32 = (uint32_t*)Raw;
  75
  76       for( i = 0; i < RawSize / 4; i++ )
  77       {
  78          low     = im32[i] & 0x0000ffff; // 3412
  79          high    = im32[i] >> 16;
  80          s32     = low;
  81          im32[i] = ( s32 << 16 ) | high;
  82       }
  83
  84    }
  85 }
  86
  87 int main(int argc, char *argv[])
  88 {
  89    START_USAGE(usage)
  90    " \n RawToDicom : \n                                                       ",
  91    " Writes a Dicom file from a Raw File                                      ",
  92    " usage: RawToDicom filein=inputFileName                                   ",
  93    "                   fileout=outputSkeletonName                             ",
  94    "                   rows=nb of Rows                                        ",
  95    "                   lines=nb of Lines,                                     ",
  96    "                   [frames = nb of Frames] //defaulted to 1               ",
  97    "                   pixeltype={8U|8S|16U|16S|32U|32S}                      ",
  98    "                   pixeltypeout={8U|8S|16U|16S|32U|32S}                   ",
  99    "                   [{b|l}] b:BigEndian,l:LittleEndian default : l         ",
 100    "                   [samples = {1|3}}       //(1:Gray,3:RGB) defaulted to 1",
 101    "                   [monochrome1]                                          ",
 102    "                   [studyid = ] [patientname = Patient's name]            ",
 103    "                   [debug]                                                ",
 104    "                                                                          ",
 105    "  monochrome1 = user wants MONOCHROME1 photom. interp. (0=white)          ",
 106    "  studyUID   : *aware* user wants to add the serie                        ",
 107    "                                             to an already existing study ",
 108    "  debug      : developper wants to run the program in 'debug mode'        ",
 109    FINISH_USAGE
 110
 111
 112    // ------------ Initialize Arguments Manager ----------------
 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 1;
 120    }
 121
 122    const char *inputFileName          = am->ArgMgrGetString("filein");
 123    const char *outputSkeletonFileName = am->ArgMgrGetString("fileout");
 124
 125   const char *patientName = am->ArgMgrGetString("patientname", "g^Fantomas");
 126
 127    int nX = am->ArgMgrWantInt("rows", usage);
 128    int nY = am->ArgMgrWantInt("lines", usage);
 129    int nZ = am->ArgMgrGetInt("frames", 1);
 130    int samplesPerPixel = am->ArgMgrGetInt("samples", 1);
 131
 132    int b = am->ArgMgrDefined("b");
 133    int l = am->ArgMgrDefined("l");
 134
 135    char *pixelType = am->ArgMgrWantString("pixeltype", usage);
 136
 137    const char *pixelTypeOut = am->ArgMgrGetString("pixeltypeout", pixelType);
 138
 139    bool monochrome1 = ( 0 != am->ArgMgrDefined("monochrome1") );
 140
 141    if (am->ArgMgrDefined("debug"))
 142       GDCM_NAME_SPACE::Debug::DebugOn();
 143
 144    bool userDefinedStudy = ( 0 != am->ArgMgrDefined("studyUID") );
 145    const char *studyUID;
 146    if (userDefinedStudy)
 147       studyUID  = am->ArgMgrGetString("studyUID");
 148
 149    // not described *on purpose* in the Usage !
 150    bool userDefinedSerie = ( 0 != am->ArgMgrDefined("serieUID") );
 151    const char *serieUID;
 152    if(userDefinedSerie)
 153       serieUID = am->ArgMgrGetString("serieUID");
 154
 155    /* if unused Param we give up */
 156    if ( am->ArgMgrPrintUnusedLabels() )
 157    {
 158       am->ArgMgrUsage(usage);
 159       delete am;
 160       return 1;
 161    }
 162
 163    delete am;  // we don't need Argument Manager any longer
 164
 165    // ----------- End Arguments Manager ---------
 166
 167  /// \TODO Deal with all the images of a directory
 168
 169  // Read the Raw file
 170    std::ifstream *Fp = new std::ifstream(inputFileName, std::ios::in | std::ios::binary);
 171    if ( ! *Fp )
 172    {
 173       std::cout << "Cannot open file: " << inputFileName;
 174       delete Fp;
 175       Fp = 0;
 176       return 0;
 177    }
 178
 179    bool bigEndian = GDCM_NAME_SPACE::Util::IsCurrentProcessorBigEndian();
 180
 181    std::string strPixelType(pixelType);
 182    int pixelSign;
 183    int pixelSize;
 184    int pixelTypeCode; // for the switch case
 185
 186    if (strPixelType == "8S")
 187    {
 188       pixelSize = 1;
 189       pixelSign = 1;
 190       pixelTypeCode = -8;
 191    }
 192    else if (strPixelType == "8U")
 193    {
 194       pixelSize = 1;
 195       pixelSign = 0;
 196       pixelTypeCode = 8;
 197    }
 198    else if (strPixelType == "16S")
 199    {
 200       pixelSize = 2;
 201       pixelSign = 1;
 202       pixelTypeCode = -16;
 203    }
 204    else if (strPixelType == "16U")
 205    {
 206       pixelSize = 2;
 207       pixelSign = 0;
 208       pixelTypeCode = 16;
 209    }
 210    else if (strPixelType == "32S")
 211    {
 212       pixelSize = 4;
 213       pixelSign = 1;
 214       pixelTypeCode = -32;
 215    }
 216    else if (strPixelType == "32U")
 217    {
 218       pixelSize = 4;
 219       pixelSign = 0;
 220       pixelTypeCode = 32;
 221    }
 222    else
 223    {
 224       std::cout << "Wrong 'pixeltype' (" << strPixelType << ")" << std::endl;
 225       return 1;
 226    }
 227
 228    std::string strPixelTypeOut(pixelTypeOut);
 229    int pixelSignOut;
 230    int pixelSizeOut;
 231    int pixelTypeOutCode; // for the swirch case
 232
 233    if (strPixelTypeOut == "8S")
 234    {
 235       pixelSizeOut = 1;
 236       pixelSignOut = 1;
 237       pixelTypeOutCode = -8;
 238    }
 239    else if (strPixelTypeOut == "8U")
 240    {
 241       pixelSizeOut = 1;
 242       pixelSignOut = 0;
 243       pixelTypeOutCode = 8;
 244    }
 245    else if (strPixelTypeOut == "16S")
 246    {
 247       pixelSizeOut = 2;
 248       pixelSignOut = 1;
 249       pixelTypeOutCode = -16;
 250    }
 251    else if (strPixelTypeOut == "16U")
 252    {
 253       pixelSizeOut = 2;
 254       pixelSignOut = 0;
 255       pixelTypeOutCode = 16;
 256    }
 257    else if (strPixelTypeOut == "32S")
 258    {
 259       pixelSizeOut = 4;
 260       pixelSignOut = 1;
 261       pixelTypeOutCode = -32;
 262    }
 263    else if (strPixelTypeOut == "32U")
 264    {
 265       pixelSizeOut = 4;
 266       pixelSignOut = 0;
 267       pixelTypeOutCode = 32;
 268    }
 269    else
 270    {
 271       std::cout << "Wrong 'pixeltypeout' (" << strPixelTypeOut << ")" << std::endl;
 272       return 1;
 273    }
 274
 275    std::string strStudyUID;
 276    std::string strSerieUID;
 277
 278    if (userDefinedStudy)
 279       strSerieUID = studyUID;
 280    else
 281       strStudyUID = GDCM_NAME_SPACE::Util::CreateUniqueUID();
 282
 283    if (userDefinedSerie)
 284      strSerieUID = serieUID;
 285    else
 286       strSerieUID = GDCM_NAME_SPACE::Util::CreateUniqueUID();
 287
 288    // Read the pixels
 289
 290    int singlePlaneDataSize =  nX*nY*samplesPerPixel*pixelSizeOut;
 291    int dataSize            =  singlePlaneDataSize*nZ;
 292
 293    uint8_t *pixels = new uint8_t[dataSize];
 294    uint8_t *planePixels = new uint8_t[singlePlaneDataSize];
 295
 296    Fp->read((char*)pixels, (size_t)dataSize);
 297
 298    if ( pixelSize !=1 && ( (l && bigEndian) || (b && ! bigEndian) ) )
 299    {
 300       ConvertSwapZone(pixelSize, pixels, dataSize);
 301    }
 302
 303 // iterate on the planes.
 304
 305    char outputFileName[200];
 306
 307 for(int nbPlanes=0; nbPlanes<nZ; nbPlanes++)
 308 {
 309
 310 sprintf(outputFileName, "%s.%04d.dcm", outputSkeletonFileName,nbPlanes);
 311
 312 // Copy (and convert) pixels of a single plane
 313
 314   switch ( pixelTypeCode ) {
 315     case 8    : CFR(PU8);  break;
 316     case -8   : CFR(PS8);  break;
 317     case -16  : CFR(PU16); break;
 318     case 16   : CFR(PS16); break;
 319     case -32  : CFR(PS32); break;
 320     case 32   : CFR(PU32); break;
 321   }
 322
 323 // Create an empty FileHelper
 324
 325    GDCM_NAME_SPACE::FileHelper *fileH = GDCM_NAME_SPACE::FileHelper::New();
 326
 327  // Get the (empty) image header.
 328    GDCM_NAME_SPACE::File *fileToBuild = fileH->GetFile();
 329
 330    // 'Study Instance UID'
 331    // The user is allowed to create his own Study,
 332    //          keeping the same 'Study Instance UID' for various images
 333    // The user may add images to a 'Manufacturer Study',
 334    //          adding new Series to an already existing Study
 335
 336    fileToBuild->InsertEntryString(strStudyUID,0x0020,0x000d,"UI");  //  Study UID
 337
 338    // 'Serie Instance UID'
 339    // The user is allowed to create his own Series,
 340    // keeping the same 'Serie Instance UID' for various images
 341    // The user shouldn't add any image to a 'Manufacturer Serie'
 342    // but there is no way no to prevent him for doing that
 343
 344    fileToBuild->InsertEntryString(strSerieUID,0x0020,0x000e,"UI");  //  Serie UID
 345
 346    std::ostringstream str;
 347
 348    // Set the image size
 349    str.str("");
 350    str << nX;
 351    fileToBuild->InsertEntryString(str.str(),0x0028,0x0011, "US"); // Columns
 352    str.str("");
 353    str << nY;
 354    fileToBuild->InsertEntryString(str.str(),0x0028,0x0010, "US"); // Rows
 355
 356    // Set the pixel type
 357
 358    str.str("");
 359    str << pixelSize*8;
 360    fileToBuild->InsertEntryString(str.str(),0x0028,0x0100, "US"); // Bits Allocated
 361
 362    str.str("");
 363    str << pixelSize*8;
 364    fileToBuild->InsertEntryString(str.str(),0x0028,0x0101, "US"); // Bits Stored
 365
 366    str.str("");
 367    str << ( pixelSize*8 - 1 );
 368    fileToBuild->InsertEntryString(str.str(),0x0028,0x0102, "US"); // High Bit
 369
 370    str.str("");
 371    str << pixelSign;
 372    fileToBuild->InsertEntryString(str.str(),0x0028,0x0103, "US"); // Pixel Representation
 373
 374 // If you deal with a Serie of images, as slices of a volume,
 375 // it up to you to tell gdcm, for each image, what are the values of :
 376 //
 377 // 0020 0032 DS 3 Image Position (Patient)
 378 // 0020 0037 DS 6 Image Orientation (Patient)
 379
 380    str.str("");
 381    str << "0.0 \\ 0.0 \\" <<  nbPlanes;
 382    fileToBuild->InsertEntryString(str.str(),0x0020,0x0032, "DS");
 383
 384    str.str("");
 385    str << samplesPerPixel;
 386    fileToBuild->InsertEntryString(str.str(),0x0028,0x0002, "US"); // Samples per Pixel
 387
 388    if (strlen(patientName) != 0)
 389       fileToBuild->InsertEntryString(patientName,0x0010,0x0010, "PN"); // Patient's Name
 390
 391    //  0=white
 392    if(monochrome1)
 393       fileH->SetPhotometricInterpretationToMonochrome1();
 394
 395 // Set the image Pixel Data (plane)
 396    fileH->SetImageData(planePixels,singlePlaneDataSize);
 397
 398 // Set the writting mode and write the image
 399    fileH->SetWriteModeToRaw();
 400
 401  // Write a DICOM Explicit VR file
 402    fileH->SetWriteTypeToDcmExplVR();
 403
 404    if( !fileH->Write(outputFileName ) )
 405    {
 406       std::cout << "Failed for [" << outputFileName << "]\n"
 407                 << "           File is unwrittable\n";
 408    }
 409
 410    fileH->Delete();
 411
 412 } // end iterate on the planes
 413
 414    delete[] pixels;
 415    delete[] planePixels;
 416    return 1;
 417 }