Example/exUnMosaicStack.cxx

   1 /*=========================================================================
   2
   3   Program:   gdcm
   4   Module:    $RCSfile: exUnMosaicStack.cxx,v $
   5   Language:  C++
   6   Date:      $Date: 2011/03/29 07:45:00 $
   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 //  UNFINISHED !  JPRx
  21
  22
  23 #include <stdlib.h> // for exit
  24 #include <math.h> // for sqrt, ceil
  25 #include "gdcmFile.h"
  26 #include "gdcmFileHelper.h"
  27
  28 bool reorganize_mosaic(const unsigned short *input, const unsigned int *inputdims, unsigned int square,
  29   const unsigned int *outputdims, unsigned short *output );
  30
  31
  32
  33 bool reorganize_mosaic(const unsigned short *input, const unsigned int *inputdims, unsigned int nbImagesPerRow,
  34   const unsigned int *outputdims, unsigned short *output )
  35 {
  36   for(unsigned x = 0; x < outputdims[0]; ++x)
  37     {
  38     for(unsigned y = 0; y < outputdims[1]; ++y)
  39       {
  40       for(unsigned z = 0; z < outputdims[2]; ++z)
  41         {
  42         output[ x + y*outputdims[0] + z*outputdims[0]*outputdims[1] ] =
  43           input[ (x + z * outputdims[0]) + (y + (z/nbImagesPerRow)*outputdims[0])*inputdims[0] ];
  44         }
  45       }
  46     }
  47   return true;
  48 }
  49
  50
  51 bool reorganize_mosaic_stack(const unsigned short *input, const unsigned int *inputdims, unsigned int nbImagesPerRow,
  52   const unsigned int *outputdims, unsigned short *output )
  53 {
  54   for(unsigned x = 0; x < outputdims[0]; ++x)
  55     {
  56     for(unsigned y = 0; y < outputdims[1]; ++y)
  57       {
  58       for(unsigned z = 0; z < outputdims[2]; ++z)
  59         {
  60         output[ x + y*outputdims[0] + z*outputdims[0]*outputdims[1] ] =
  61           input[ (x + z * outputdims[0]) + (y + (z/nbImagesPerRow)*outputdims[0])*inputdims[0] ];
  62         }
  63       }
  64     }
  65   return true;
  66 }
  67
  68
  69 int main(int argc, char *argv[])
  70 {
  71
  72    std::string fileName       = argv[1];
  73    std::string outputFileName = argv[2];
  74    unsigned int nbImagesPerRow = 8;
  75    unsigned int numberOfImagesInMosaic = 52;
  76    unsigned int numberMosaics = 52;
  77
  78    unsigned int inputdims[2];
  79    inputdims[0] = 1024;
  80    inputdims[1] = 1024;
  81
  82    unsigned int outputdims[3];
  83
  84 // ============================================================
  85 //   Read the input image.
  86 // ============================================================
  87
  88    std::cout << argv[1] << std::endl;
  89
  90       std::cout << "dirName [" << dirName << "]" << std::endl;
  91
  92       GDCM_NAME_SPACE::DirList dirList(dirName,rec); // gets recursively (or not) the file list
  93       GDCM_NAME_SPACE::DirListType fileList = dirList.GetFilenames();
  94       GDCM_NAME_SPACE::File *f;
  95       bool res;
  96
  97       if (fileList.size() == 0)
  98       {
  99          std::cout << "No file found in : [" << dirName << "]" << std::endl;
 100       }
 101
 102
 103
 104
 105    GDCM_NAME_SPACE::File *f = GDCM_NAME_SPACE::File::New();
 106    f->SetLoadMode( GDCM_NAME_SPACE::LD_ALL);
 107    f->SetFileName( fileName );
 108    bool res = f->Load();
 109
 110    if (!res) {
 111        std::cerr << "Sorry, " << fileName <<"  not a gdcm-readable "
 112                  << "DICOM / ACR File"
 113                  <<std::endl;
 114        f->Delete();
 115        return 0;
 116    }
 117    std::cout << " ... is readable " << std::endl;
 118
 119
 120 // ============================================================
 121 //   Read the Pixels
 122 //
 123 // ============================================================
 124
 125    // Pixel Reading must be done here, to be sure
 126    // to load the Palettes Color (if any)
 127
 128    // First, create a GDCM_NAME_SPACE::FileHelper
 129    GDCM_NAME_SPACE::FileHelper *fh = GDCM_NAME_SPACE::FileHelper::New(f);
 130
 131    // Load the pixels, DO NOT transform LUT (if any) into RGB Pixels
 132    uint8_t *imageDataRaw = fh->GetImageDataRaw();
 133    // Get the image data size
 134    size_t dataRawSize    = fh->GetImageDataRawSize();
 135
 136  // ============================================================
 137 //   Create a new GDCM_NAME_SPACE::Filehelper, to hold new image.
 138 // ============================================================
 139
 140    GDCM_NAME_SPACE::FileHelper *copy = GDCM_NAME_SPACE::FileHelper::New( );
 141    copy->SetFileName( outputFileName );
 142   // copy->Load();
 143
 144   unsigned int div = (unsigned int)ceil(sqrt( (double)numberOfImagesInMosaic ) );
 145   outputdims[0] = inputdims[0]/ div;
 146   outputdims[1] = inputdims[1] / div;
 147   outputdims[2] = numberOfImagesInMosaic;
 148
 149 /*
 150   if (outputdims[0] *  outputdims[1] * outputdims[2] * 2 != dataRawSize)
 151   {
 152      std::cout << "outputdims[0] : " << outputdims[0] << " outputdims[1] : " <<  outputdims[1] << " numberOfImagesInMosaic : " << numberOfImagesInMosaic
 153                << " ( = " << outputdims[0] *  outputdims[1] * outputdims[2] * 2 << " ) doesn't match with dataRawSize : " << dataRawSize << std::endl;
 154      //exit(0);
 155   }
 156 */
 157
 158
 159   unsigned short *input = (unsigned short *)imageDataRaw;
 160
 161   unsigned short *output = (unsigned short *) malloc(outputdims[0] * outputdims[1] * numberOfImagesInMosaic * 2);
 162
 163   reorganize_mosaic(input, inputdims,
 164                     nbImagesPerRow,
 165                     outputdims, output);
 166
 167   FILE *fPixels;
 168   fPixels = fopen(outputFileName.c_str(), "w");
 169   fwrite(output,2, outputdims[0] * outputdims[1] * numberOfImagesInMosaic, fPixels);
 170 }
 171