src/gdcmPixelReadConvert.cxx

   1 /*=========================================================================
   2
   3   Program:   gdcm
   4   Module:    $RCSfile: gdcmPixelReadConvert.cxx,v $
   5   Language:  C++
   6   Date:      $Date: 2005/06/13 15:43:48 $
   7   Version:   $Revision: 1.64 $
   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 "gdcmPixelReadConvert.h"
  20 #include "gdcmDebug.h"
  21 #include "gdcmFile.h"
  22 #include "gdcmGlobal.h"
  23 #include "gdcmTS.h"
  24 #include "gdcmDocEntry.h"
  25 #include "gdcmRLEFramesInfo.h"
  26 #include "gdcmJPEGFragmentsInfo.h"
  27
  28 #include <fstream>
  29 #include <stdio.h> //for sscanf
  30
  31 namespace gdcm
  32 {
  33
  34 //bool ReadMPEGFile (std::ifstream *fp, void *image_buffer, size_t lenght);
  35 bool gdcm_read_JPEG2000_file (void* raw,
  36                               char *inputdata, size_t inputlength);
  37 //-----------------------------------------------------------------------------
  38 #define str2num(str, typeNum) *((typeNum *)(str))
  39
  40 //-----------------------------------------------------------------------------
  41 // Constructor / Destructor
  42 /// Constructor
  43 PixelReadConvert::PixelReadConvert()
  44 {
  45    RGB          = 0;
  46    RGBSize      = 0;
  47    Raw          = 0;
  48    RawSize      = 0;
  49    LutRGBA      = 0;
  50    LutRedData   = 0;
  51    LutGreenData = 0;
  52    LutBlueData  = 0;
  53 }
  54
  55 /// Canonical Destructor
  56 PixelReadConvert::~PixelReadConvert()
  57 {
  58    Squeeze();
  59 }
  60
  61 //-----------------------------------------------------------------------------
  62 // Public
  63 /**
  64  * \brief Predicate to know whether the image[s] (once Raw) is RGB.
  65  * \note See comments of \ref ConvertHandleColor
  66  */
  67 bool PixelReadConvert::IsRawRGB()
  68 {
  69    if (   IsMonochrome
  70        || PlanarConfiguration == 2
  71        || IsPaletteColor )
  72    {
  73       return false;
  74    }
  75    return true;
  76 }
  77 /**
  78  * \brief Gets various usefull informations from the file header
  79  * @param file gdcm::File pointer
  80  */
  81 void PixelReadConvert::GrabInformationsFromFile( File *file )
  82 {
  83    // Number of Bits Allocated for storing a Pixel is defaulted to 16
  84    // when absent from the file.
  85    BitsAllocated = file->GetBitsAllocated();
  86    if ( BitsAllocated == 0 )
  87    {
  88       BitsAllocated = 16;
  89    }
  90
  91    // Number of "Bits Stored", defaulted to number of "Bits Allocated"
  92    // when absent from the file.
  93    BitsStored = file->GetBitsStored();
  94    if ( BitsStored == 0 )
  95    {
  96       BitsStored = BitsAllocated;
  97    }
  98
  99    // High Bit Position, defaulted to "Bits Allocated" - 1
 100    HighBitPosition = file->GetHighBitPosition();
 101    if ( HighBitPosition == 0 )
 102    {
 103       HighBitPosition = BitsAllocated - 1;
 104    }
 105
 106    XSize           = file->GetXSize();
 107    YSize           = file->GetYSize();
 108    ZSize           = file->GetZSize();
 109    SamplesPerPixel = file->GetSamplesPerPixel();
 110    PixelSize       = file->GetPixelSize();
 111    PixelSign       = file->IsSignedPixelData();
 112    SwapCode        = file->GetSwapCode();
 113    std::string ts  = file->GetTransferSyntax();
 114    IsRaw =
 115         ( ! file->IsDicomV3() )
 116      || Global::GetTS()->GetSpecialTransferSyntax(ts) == TS::ImplicitVRLittleEndian
 117      || Global::GetTS()->GetSpecialTransferSyntax(ts) == TS::ImplicitVRLittleEndianDLXGE
 118      || Global::GetTS()->GetSpecialTransferSyntax(ts) == TS::ExplicitVRLittleEndian
 119      || Global::GetTS()->GetSpecialTransferSyntax(ts) == TS::ExplicitVRBigEndian
 120      || Global::GetTS()->GetSpecialTransferSyntax(ts) == TS::DeflatedExplicitVRLittleEndian;
 121
 122    IsMPEG          = Global::GetTS()->IsMPEG(ts);
 123    IsJPEG2000      = Global::GetTS()->IsJPEG2000(ts);
 124    IsJPEGLS        = Global::GetTS()->IsJPEGLS(ts);
 125    IsJPEGLossy     = Global::GetTS()->IsJPEGLossy(ts);
 126    IsJPEGLossless  = Global::GetTS()->IsJPEGLossless(ts);
 127    IsRLELossless   = Global::GetTS()->IsRLELossless(ts);
 128
 129    PixelOffset     = file->GetPixelOffset();
 130    PixelDataLength = file->GetPixelAreaLength();
 131    RLEInfo         = file->GetRLEInfo();
 132    JPEGInfo        = file->GetJPEGInfo();
 133
 134    IsMonochrome    = file->IsMonochrome();
 135    IsMonochrome1   = file->IsMonochrome1();
 136    IsPaletteColor  = file->IsPaletteColor();
 137    IsYBRFull       = file->IsYBRFull();
 138
 139    PlanarConfiguration = file->GetPlanarConfiguration();
 140
 141    /////////////////////////////////////////////////////////////////
 142    // LUT section:
 143    HasLUT = file->HasLUT();
 144    if ( HasLUT )
 145    {
 146       // Just in case some access to a File element requires disk access.
 147       LutRedDescriptor   = file->GetEntryValue( 0x0028, 0x1101 );
 148       LutGreenDescriptor = file->GetEntryValue( 0x0028, 0x1102 );
 149       LutBlueDescriptor  = file->GetEntryValue( 0x0028, 0x1103 );
 150
 151       // Depending on the value of Document::MAX_SIZE_LOAD_ELEMENT_VALUE
 152       // [ refer to invocation of Document::SetMaxSizeLoadEntry() in
 153       // Document::Document() ], the loading of the value (content) of a
 154       // [Bin|Val]Entry occurence migth have been hindered (read simply NOT
 155       // loaded). Hence, we first try to obtain the LUTs data from the file
 156       // and when this fails we read the LUTs data directly from disk.
 157       // \TODO Reading a [Bin|Val]Entry directly from disk is a kludge.
 158       //       We should NOT bypass the [Bin|Val]Entry class. Instead
 159       //       an access to an UNLOADED content of a [Bin|Val]Entry occurence
 160       //       (e.g. BinEntry::GetBinArea()) should force disk access from
 161       //       within the [Bin|Val]Entry class itself. The only problem
 162       //       is that the [Bin|Val]Entry is unaware of the FILE* is was
 163       //       parsed from. Fix that. FIXME.
 164
 165       // //// Red round
 166       file->LoadEntryBinArea(0x0028, 0x1201);
 167       LutRedData = (uint8_t*)file->GetEntryBinArea( 0x0028, 0x1201 );
 168       if ( ! LutRedData )
 169       {
 170          gdcmWarningMacro( "Unable to read Red LUT data" );
 171       }
 172
 173       // //// Green round:
 174       file->LoadEntryBinArea(0x0028, 0x1202);
 175       LutGreenData = (uint8_t*)file->GetEntryBinArea(0x0028, 0x1202 );
 176       if ( ! LutGreenData)
 177       {
 178          gdcmWarningMacro( "Unable to read Green LUT data" );
 179       }
 180
 181       // //// Blue round:
 182       file->LoadEntryBinArea(0x0028, 0x1203);
 183       LutBlueData = (uint8_t*)file->GetEntryBinArea( 0x0028, 0x1203 );
 184       if ( ! LutBlueData )
 185       {
 186          gdcmWarningMacro( "Unable to read Blue LUT data" );
 187       }
 188    }
 189
 190    ComputeRawAndRGBSizes();
 191 }
 192
 193 /// \brief Reads from disk and decompresses Pixels
 194 bool PixelReadConvert::ReadAndDecompressPixelData( std::ifstream *fp )
 195 {
 196    // ComputeRawAndRGBSizes is already made by
 197    // ::GrabInformationsFromfile. So, the structure sizes are
 198    // correct
 199    Squeeze();
 200
 201    //////////////////////////////////////////////////
 202    //// First stage: get our hands on the Pixel Data.
 203    if ( !fp )
 204    {
 205       gdcmWarningMacro( "Unavailable file pointer." );
 206       return false;
 207    }
 208
 209    fp->seekg( PixelOffset, std::ios::beg );
 210    if( fp->fail() || fp->eof())
 211    {
 212       gdcmWarningMacro( "Unable to find PixelOffset in file." );
 213       return false;
 214    }
 215
 216    AllocateRaw();
 217
 218    //////////////////////////////////////////////////
 219    //// Second stage: read from disk dans decompress.
 220    if ( BitsAllocated == 12 )
 221    {
 222       ReadAndDecompress12BitsTo16Bits( fp);
 223    }
 224    else if ( IsRaw )
 225    {
 226       // This problem can be found when some obvious informations are found
 227       // after the field containing the image data. In this case, these
 228       // bad data are added to the size of the image (in the PixelDataLength
 229       // variable). But RawSize is the right size of the image !
 230       if( PixelDataLength != RawSize)
 231       {
 232          gdcmWarningMacro( "Mismatch between PixelReadConvert : "
 233                             << PixelDataLength << " and RawSize : " << RawSize );
 234       }
 235       if( PixelDataLength > RawSize)
 236       {
 237          fp->read( (char*)Raw, RawSize);
 238       }
 239       else
 240       {
 241          fp->read( (char*)Raw, PixelDataLength);
 242       }
 243
 244       if ( fp->fail() || fp->eof())
 245       {
 246          gdcmWarningMacro( "Reading of Raw pixel data failed." );
 247          return false;
 248       }
 249    }
 250    else if ( IsRLELossless )
 251    {
 252       if ( ! RLEInfo->DecompressRLEFile( fp, Raw, XSize, YSize, ZSize, BitsAllocated ) )
 253       {
 254          gdcmWarningMacro( "RLE decompressor failed." );
 255          return false;
 256       }
 257    }
 258    else if ( IsMPEG )
 259    {
 260       //gdcmWarningMacro( "Sorry, MPEG not yet taken into account" );
 261       //return false;
 262 //      ReadMPEGFile(fp, Raw, PixelDataLength); // fp has already been seek to start of mpeg
 263       return true;
 264    }
 265    else
 266    {
 267       // Default case concerns JPEG family
 268       if ( ! ReadAndDecompressJPEGFile( fp ) )
 269       {
 270          gdcmWarningMacro( "JPEG decompressor failed." );
 271          return false;
 272       }
 273    }
 274
 275    ////////////////////////////////////////////
 276    //// Third stage: twigle the bytes and bits.
 277    ConvertReorderEndianity();
 278    ConvertReArrangeBits();
 279    ConvertFixGreyLevels();
 280    ConvertHandleColor();
 281
 282    return true;
 283 }
 284
 285 /// Deletes Pixels Area
 286 void PixelReadConvert::Squeeze()
 287 {
 288    if ( RGB )
 289       delete [] RGB;
 290    RGB = 0;
 291
 292    if ( Raw )
 293       delete [] Raw;
 294    Raw = 0;
 295
 296    if ( LutRGBA )
 297       delete [] LutRGBA;
 298    LutRGBA = 0;
 299 }
 300
 301 /**
 302  * \brief Build the RGB image from the Raw imagage and the LUTs.
 303  */
 304 bool PixelReadConvert::BuildRGBImage()
 305 {
 306    if ( RGB )
 307    {
 308       // The job is already done.
 309       return true;
 310    }
 311
 312    if ( ! Raw )
 313    {
 314       // The job can't be done
 315       return false;
 316    }
 317
 318    BuildLUTRGBA();
 319    if ( ! LutRGBA )
 320    {
 321       // The job can't be done
 322       return false;
 323    }
 324
 325    // Build RGB Pixels
 326    AllocateRGB();
 327    uint8_t *localRGB = RGB;
 328    for (size_t i = 0; i < RawSize; ++i )
 329    {
 330       int j  = Raw[i] * 4;
 331       *localRGB++ = LutRGBA[j];
 332       *localRGB++ = LutRGBA[j+1];
 333       *localRGB++ = LutRGBA[j+2];
 334    }
 335    return true;
 336 }
 337
 338 //-----------------------------------------------------------------------------
 339 // Protected
 340
 341 //-----------------------------------------------------------------------------
 342 // Private
 343 /**
 344  * \brief Read from file a 12 bits per pixel image and decompress it
 345  *        into a 16 bits per pixel image.
 346  */
 347 void PixelReadConvert::ReadAndDecompress12BitsTo16Bits( std::ifstream *fp )
 348                throw ( FormatError )
 349 {
 350    int nbPixels = XSize * YSize;
 351    uint16_t *localDecompres = (uint16_t*)Raw;
 352
 353    for( int p = 0; p < nbPixels; p += 2 )
 354    {
 355       uint8_t b0, b1, b2;
 356
 357       fp->read( (char*)&b0, 1);
 358       if ( fp->fail() || fp->eof() )
 359       {
 360          throw FormatError( "PixelReadConvert::ReadAndDecompress12BitsTo16Bits()",
 361                                 "Unfound first block" );
 362       }
 363
 364       fp->read( (char*)&b1, 1 );
 365       if ( fp->fail() || fp->eof())
 366       {
 367          throw FormatError( "PixelReadConvert::ReadAndDecompress12BitsTo16Bits()",
 368                                 "Unfound second block" );
 369       }
 370
 371       fp->read( (char*)&b2, 1 );
 372       if ( fp->fail() || fp->eof())
 373       {
 374          throw FormatError( "PixelReadConvert::ReadAndDecompress12BitsTo16Bits()",
 375                                 "Unfound second block" );
 376       }
 377
 378       // Two steps are necessary to please VC++
 379       //
 380       // 2 pixels 12bit =     [0xABCDEF]
 381       // 2 pixels 16bit = [0x0ABD] + [0x0FCE]
 382       //                        A                     B                 D
 383       *localDecompres++ =  ((b0 >> 4) << 8) + ((b0 & 0x0f) << 4) + (b1 & 0x0f);
 384       //                        F                     C                 E
 385       *localDecompres++ =  ((b2 & 0x0f) << 8) + ((b1 >> 4) << 4) + (b2 >> 4);
 386
 387       /// \todo JPR Troubles expected on Big-Endian processors ?
 388    }
 389 }
 390
 391 /**
 392  * \brief     Reads from disk the Pixel Data of JPEG Dicom encapsulated
 393  *            file and decompress it.
 394  * @param     fp File Pointer
 395  * @return    Boolean
 396  */
 397 bool PixelReadConvert::ReadAndDecompressJPEGFile( std::ifstream *fp )
 398 {
 399    if ( IsJPEG2000 )
 400    {
 401      // FIXME this is really ugly but it seems I have to load the complete
 402      // jpeg2000 stream to use jasper:
 403       // I don't think we'll ever be able to deal with multiple fragments properly
 404
 405       unsigned long inputlength = 0;
 406       JPEGFragment *jpegfrag = JPEGInfo->GetFirstFragment();
 407       while( jpegfrag )
 408       {
 409          inputlength += jpegfrag->GetLength();
 410          jpegfrag = JPEGInfo->GetNextFragment();
 411       }
 412       gdcmAssertMacro( inputlength != 0);
 413       uint8_t *inputdata = new uint8_t[inputlength];
 414       char *pinputdata = (char*)inputdata;
 415       jpegfrag = JPEGInfo->GetFirstFragment();
 416       while( jpegfrag )
 417       {
 418          fp->seekg( jpegfrag->GetOffset(), std::ios::beg);
 419          fp->read(pinputdata, jpegfrag->GetLength());
 420          pinputdata += jpegfrag->GetLength();
 421          jpegfrag = JPEGInfo->GetNextFragment();
 422       }
 423       // Warning the inputdata buffer is delete in the function
 424       if ( ! gdcm_read_JPEG2000_file( Raw,
 425           (char*)inputdata, inputlength ) )
 426       {
 427          return true;
 428       }
 429    }
 430
 431    if ( IsJPEGLS )
 432    {
 433    // WARNING : JPEG-LS is NOT the 'classical' Jpeg Lossless :
 434    // [JPEG-LS is the basis for new lossless/near-lossless compression
 435    // standard for continuous-tone images intended for JPEG2000. The standard
 436    // is based on the LOCO-I algorithm (LOw COmplexity LOssless COmpression
 437    // for Images) developed at Hewlett-Packard Laboratories]
 438    //
 439    // see http://datacompression.info/JPEGLS.shtml
 440    //
 441
 442       gdcmWarningMacro( "Sorry, JPEG-LS not yet taken into account" );
 443       fp->seekg( JPEGInfo->GetFirstFragment()->GetOffset(), std::ios::beg);
 444 //    if ( ! gdcm_read_JPEGLS_file( fp,Raw ) )
 445          return false;
 446    }
 447    else
 448      {
 449      // Precompute the offset localRaw will be shifted with
 450      int length = XSize * YSize * SamplesPerPixel;
 451      int numberBytes = BitsAllocated / 8;
 452
 453      JPEGInfo->DecompressFromFile(fp, Raw, BitsStored, numberBytes, length );
 454      return true;
 455      }
 456 }
 457
 458 /**
 459  * \brief Build Red/Green/Blue/Alpha LUT from File
 460  *         when (0028,0004),Photometric Interpretation = [PALETTE COLOR ]
 461  *          and (0028,1101),(0028,1102),(0028,1102)
 462  *            - xxx Palette Color Lookup Table Descriptor - are found
 463  *          and (0028,1201),(0028,1202),(0028,1202)
 464  *            - xxx Palette Color Lookup Table Data - are found
 465  * \warning does NOT deal with :
 466  *   0028 1100 Gray Lookup Table Descriptor (Retired)
 467  *   0028 1221 Segmented Red Palette Color Lookup Table Data
 468  *   0028 1222 Segmented Green Palette Color Lookup Table Data
 469  *   0028 1223 Segmented Blue Palette Color Lookup Table Data
 470  *   no known Dicom reader deals with them :-(
 471  * @return a RGBA Lookup Table
 472  */
 473 void PixelReadConvert::BuildLUTRGBA()
 474 {
 475    if ( LutRGBA )
 476    {
 477       return;
 478    }
 479    // Not so easy : see
 480    // http://www.barre.nom.fr/medical/dicom2/limitations.html#Color%20Lookup%20Tables
 481
 482    if ( ! IsPaletteColor )
 483    {
 484       return;
 485    }
 486
 487    if (   LutRedDescriptor   == GDCM_UNFOUND
 488        || LutGreenDescriptor == GDCM_UNFOUND
 489        || LutBlueDescriptor  == GDCM_UNFOUND )
 490    {
 491       gdcmWarningMacro( "(At least) a LUT Descriptor is missing" );
 492       return;
 493    }
 494
 495    ////////////////////////////////////////////
 496    // Extract the info from the LUT descriptors
 497    int lengthR;   // Red LUT length in Bytes
 498    int debR;      // Subscript of the first Lut Value
 499    int nbitsR;    // Lut item size (in Bits)
 500    int nbRead = sscanf( LutRedDescriptor.c_str(),
 501                         "%d\\%d\\%d",
 502                         &lengthR, &debR, &nbitsR );
 503    if( nbRead != 3 )
 504    {
 505       gdcmWarningMacro( "Wrong Red LUT descriptor" );
 506    }
 507
 508    int lengthG;  // Green LUT length in Bytes
 509    int debG;     // Subscript of the first Lut Value
 510    int nbitsG;   // Lut item size (in Bits)
 511    nbRead = sscanf( LutGreenDescriptor.c_str(),
 512                     "%d\\%d\\%d",
 513                     &lengthG, &debG, &nbitsG );
 514
 515    if( nbRead != 3 )
 516    {
 517       gdcmWarningMacro( "Wrong Green LUT descriptor" );
 518    }
 519
 520    int lengthB;  // Blue LUT length in Bytes
 521    int debB;     // Subscript of the first Lut Value
 522    int nbitsB;   // Lut item size (in Bits)
 523    nbRead = sscanf( LutRedDescriptor.c_str(),
 524                     "%d\\%d\\%d",
 525                     &lengthB, &debB, &nbitsB );
 526    gdcmWarningMacro(" lengthR " << lengthR << " debR "
 527                  << debR << " nbitsR " << nbitsR);
 528    gdcmWarningMacro(" lengthG " << lengthG << " debG "
 529                  << debG << " nbitsG " << nbitsG);
 530    gdcmWarningMacro(" lengthB " << lengthB << " debB "
 531                  << debB << " nbitsB " << nbitsB);
 532    if( nbRead != 3 )
 533    {
 534       gdcmWarningMacro( "Wrong Blue LUT descriptor" );
 535    }
 536
 537    if ( !lengthR ) // if = 2^16, this shall be 0 see : CP-143
 538       lengthR=65536;
 539    if( !lengthG ) // if = 2^16, this shall be 0
 540       lengthG=65536;
 541    if ( !lengthB ) // if = 2^16, this shall be 0
 542       lengthB=65536;
 543
 544    ////////////////////////////////////////////////////////
 545
 546    if ( ( ! LutRedData ) || ( ! LutGreenData ) || ( ! LutBlueData ) )
 547    {
 548       gdcmWarningMacro( "(At least) a LUT is missing" );
 549       return;
 550    }
 551
 552    // -------------------------------------------------------------
 553
 554    if ( BitsAllocated <= 8)
 555    {
 556
 557       // forge the 4 * 8 Bits Red/Green/Blue/Alpha LUT
 558       LutRGBA = new uint8_t[ 1024 ]; // 256 * 4 (R, G, B, Alpha)
 559       if ( !LutRGBA )
 560          return;
 561
 562       memset( LutRGBA, 0, 1024 );
 563
 564       int mult;
 565       if ( ( nbitsR == 16 ) && ( BitsAllocated == 8 ) )
 566       {
 567          // when LUT item size is different than pixel size
 568          mult = 2; // high byte must be = low byte
 569       }
 570       else
 571       {
 572          // See PS 3.3-2003 C.11.1.1.2 p 619
 573          mult = 1;
 574       }
 575
 576       // if we get a black image, let's just remove the '+1'
 577       // from 'i*mult+1' and check again
 578       // if it works, we shall have to check the 3 Palettes
 579       // to see which byte is ==0 (first one, or second one)
 580       // and fix the code
 581       // We give up the checking to avoid some (useless ?) overhead
 582       // (optimistic asumption)
 583       int i;
 584       uint8_t *a;
 585
 586       //take "Subscript of the first Lut Value" (debR,debG,debB) into account!
 587
 588       a = LutRGBA + 0 + debR;
 589       for( i=0; i < lengthR; ++i )
 590       {
 591          *a = LutRedData[i*mult+1];
 592          a += 4;
 593       }
 594
 595       a = LutRGBA + 1 + debG;
 596       for( i=0; i < lengthG; ++i)
 597       {
 598          *a = LutGreenData[i*mult+1];
 599          a += 4;
 600       }
 601
 602       a = LutRGBA + 2 + debB;
 603       for(i=0; i < lengthB; ++i)
 604       {
 605          *a = LutBlueData[i*mult+1];
 606          a += 4;
 607       }
 608
 609       a = LutRGBA + 3 ;
 610       for(i=0; i < 256; ++i)
 611       {
 612          *a = 1; // Alpha component
 613          a += 4;
 614       }
 615    }
 616    else
 617    {
 618       gdcmWarningMacro( "Sorry Palette Color Lookup Tables not yet dealt with"
 619                          << "for 16 Bits Per Pixel images" );
 620    }
 621 }
 622
 623 /**
 624  * \brief Swap the bytes, according to \ref SwapCode.
 625  */
 626 void PixelReadConvert::ConvertSwapZone()
 627 {
 628    unsigned int i;
 629
 630    if( BitsAllocated == 16 )
 631    {
 632       uint16_t *im16 = (uint16_t*)Raw;
 633       switch( SwapCode )
 634       {
 635          case 1234:
 636             break;
 637          case 3412:
 638          case 2143:
 639          case 4321:
 640             for( i = 0; i < RawSize / 2; i++ )
 641             {
 642                im16[i]= (im16[i] >> 8) | (im16[i] << 8 );
 643             }
 644             break;
 645          default:
 646             gdcmWarningMacro("SwapCode value (16 bits) not allowed.");
 647       }
 648    }
 649    else if( BitsAllocated == 32 )
 650    {
 651       uint32_t s32;
 652       uint16_t high;
 653       uint16_t low;
 654       uint32_t *im32 = (uint32_t*)Raw;
 655       switch ( SwapCode )
 656       {
 657          case 1234:
 658             break;
 659          case 4321:
 660             for( i = 0; i < RawSize / 4; i++ )
 661             {
 662                low     = im32[i] & 0x0000ffff;  // 4321
 663                high    = im32[i] >> 16;
 664                high    = ( high >> 8 ) | ( high << 8 );
 665                low     = ( low  >> 8 ) | ( low  << 8 );
 666                s32     = low;
 667                im32[i] = ( s32 << 16 ) | high;
 668             }
 669             break;
 670          case 2143:
 671             for( i = 0; i < RawSize / 4; i++ )
 672             {
 673                low     = im32[i] & 0x0000ffff;   // 2143
 674                high    = im32[i] >> 16;
 675                high    = ( high >> 8 ) | ( high << 8 );
 676                low     = ( low  >> 8 ) | ( low  << 8 );
 677                s32     = high;
 678                im32[i] = ( s32 << 16 ) | low;
 679             }
 680             break;
 681          case 3412:
 682             for( i = 0; i < RawSize / 4; i++ )
 683             {
 684                low     = im32[i] & 0x0000ffff; // 3412
 685                high    = im32[i] >> 16;
 686                s32     = low;
 687                im32[i] = ( s32 << 16 ) | high;
 688             }
 689             break;
 690          default:
 691             gdcmWarningMacro("SwapCode value (32 bits) not allowed." );
 692       }
 693    }
 694 }
 695
 696 /**
 697  * \brief Deal with endianness i.e. re-arange bytes inside the integer
 698  */
 699 void PixelReadConvert::ConvertReorderEndianity()
 700 {
 701    if ( BitsAllocated != 8 )
 702    {
 703       ConvertSwapZone();
 704    }
 705
 706    // Special kludge in order to deal with xmedcon broken images:
 707    if ( BitsAllocated == 16
 708      && BitsStored < BitsAllocated
 709      && !PixelSign )
 710    {
 711       int l = (int)( RawSize / ( BitsAllocated / 8 ) );
 712       uint16_t *deb = (uint16_t *)Raw;
 713       for(int i = 0; i<l; i++)
 714       {
 715          if( *deb == 0xffff )
 716          {
 717            *deb = 0;
 718          }
 719          deb++;
 720       }
 721    }
 722 }
 723
 724 /**
 725  * \brief Deal with Grey levels i.e. re-arange them
 726  *        to have low values = dark, high values = bright
 727  */
 728 void PixelReadConvert::ConvertFixGreyLevels()
 729 {
 730    if (!IsMonochrome1)
 731       return;
 732
 733    uint32_t i; // to please M$VC6
 734    int16_t j;
 735
 736    if (!PixelSign)
 737    {
 738       if ( BitsAllocated == 8 )
 739       {
 740          uint8_t *deb = (uint8_t *)Raw;
 741          for (i=0; i<RawSize; i++)
 742          {
 743             *deb = 255 - *deb;
 744             deb++;
 745          }
 746          return;
 747       }
 748
 749       if ( BitsAllocated == 16 )
 750       {
 751          uint16_t mask =1;
 752          for (j=0; j<BitsStored-1; j++)
 753          {
 754             mask = (mask << 1) +1; // will be fff when BitsStored=12
 755          }
 756
 757          uint16_t *deb = (uint16_t *)Raw;
 758          for (i=0; i<RawSize/2; i++)
 759          {
 760             *deb = mask - *deb;
 761             deb++;
 762          }
 763          return;
 764        }
 765    }
 766    else
 767    {
 768       if ( BitsAllocated == 8 )
 769       {
 770          uint8_t smask8 = 255;
 771          uint8_t *deb = (uint8_t *)Raw;
 772          for (i=0; i<RawSize; i++)
 773          {
 774             *deb = smask8 - *deb;
 775             deb++;
 776          }
 777          return;
 778       }
 779       if ( BitsAllocated == 16 )
 780       {
 781          uint16_t smask16 = 65535;
 782          uint16_t *deb = (uint16_t *)Raw;
 783          for (i=0; i<RawSize/2; i++)
 784          {
 785             *deb = smask16 - *deb;
 786             deb++;
 787          }
 788          return;
 789       }
 790    }
 791 }
 792
 793 /**
 794  * \brief  Re-arrange the bits within the bytes.
 795  * @return Boolean always true
 796  */
 797 bool PixelReadConvert::ConvertReArrangeBits() throw ( FormatError )
 798 {
 799    if ( BitsStored != BitsAllocated )
 800    {
 801       int l = (int)( RawSize / ( BitsAllocated / 8 ) );
 802       if ( BitsAllocated == 16 )
 803       {
 804          uint16_t mask = 0xffff;
 805          mask = mask >> ( BitsAllocated - BitsStored );
 806          uint16_t *deb = (uint16_t*)Raw;
 807          for(int i = 0; i<l; i++)
 808          {
 809             *deb = (*deb >> (BitsStored - HighBitPosition - 1)) & mask;
 810             deb++;
 811          }
 812       }
 813       else if ( BitsAllocated == 32 )
 814       {
 815          uint32_t mask = 0xffffffff;
 816          mask = mask >> ( BitsAllocated - BitsStored );
 817          uint32_t *deb = (uint32_t*)Raw;
 818          for(int i = 0; i<l; i++)
 819          {
 820             *deb = (*deb >> (BitsStored - HighBitPosition - 1)) & mask;
 821             deb++;
 822          }
 823       }
 824       else
 825       {
 826          gdcmWarningMacro("Weird image");
 827          throw FormatError( "Weird image !?" );
 828       }
 829    }
 830    return true;
 831 }
 832
 833 /**
 834  * \brief   Convert (Red plane, Green plane, Blue plane) to RGB pixels
 835  * \warning Works on all the frames at a time
 836  */
 837 void PixelReadConvert::ConvertRGBPlanesToRGBPixels()
 838 {
 839    uint8_t *localRaw = Raw;
 840    uint8_t *copyRaw = new uint8_t[ RawSize ];
 841    memmove( copyRaw, localRaw, RawSize );
 842
 843    int l = XSize * YSize * ZSize;
 844
 845    uint8_t *a = copyRaw;
 846    uint8_t *b = copyRaw + l;
 847    uint8_t *c = copyRaw + l + l;
 848
 849    for (int j = 0; j < l; j++)
 850    {
 851       *(localRaw++) = *(a++);
 852       *(localRaw++) = *(b++);
 853       *(localRaw++) = *(c++);
 854    }
 855    delete[] copyRaw;
 856 }
 857
 858 /**
 859  * \brief   Convert (cY plane, cB plane, cR plane) to RGB pixels
 860  * \warning Works on all the frames at a time
 861  */
 862 void PixelReadConvert::ConvertYcBcRPlanesToRGBPixels()
 863 {
 864    uint8_t *localRaw = Raw;
 865    uint8_t *copyRaw = new uint8_t[ RawSize ];
 866    memmove( copyRaw, localRaw, RawSize );
 867
 868    // to see the tricks about YBR_FULL, YBR_FULL_422,
 869    // YBR_PARTIAL_422, YBR_ICT, YBR_RCT have a look at :
 870    // ftp://medical.nema.org/medical/dicom/final/sup61_ft.pdf
 871    // and be *very* affraid
 872    //
 873    int l        = XSize * YSize;
 874    int nbFrames = ZSize;
 875
 876    uint8_t *a = copyRaw + 0;
 877    uint8_t *b = copyRaw + l;
 878    uint8_t *c = copyRaw + l+ l;
 879    int32_t R, G, B;
 880
 881    /// \todo : Replace by the 'well known' integer computation
 882    ///         counterpart. Refer to
 883    ///            http://lestourtereaux.free.fr/papers/data/yuvrgb.pdf
 884    ///         for code optimisation.
 885
 886    for ( int i = 0; i < nbFrames; i++ )
 887    {
 888       for ( int j = 0; j < l; j++ )
 889       {
 890          R = 38142 *(*a-16) + 52298 *(*c -128);
 891          G = 38142 *(*a-16) - 26640 *(*c -128) - 12845 *(*b -128);
 892          B = 38142 *(*a-16) + 66093 *(*b -128);
 893
 894          R = (R+16384)>>15;
 895          G = (G+16384)>>15;
 896          B = (B+16384)>>15;
 897
 898          if (R < 0)   R = 0;
 899          if (G < 0)   G = 0;
 900          if (B < 0)   B = 0;
 901          if (R > 255) R = 255;
 902          if (G > 255) G = 255;
 903          if (B > 255) B = 255;
 904
 905          *(localRaw++) = (uint8_t)R;
 906          *(localRaw++) = (uint8_t)G;
 907          *(localRaw++) = (uint8_t)B;
 908          a++;
 909          b++;
 910          c++;
 911       }
 912    }
 913    delete[] copyRaw;
 914 }
 915
 916 /// \brief Deals with the color decoding i.e. handle:
 917 ///   - R, G, B planes (as opposed to RGB pixels)
 918 ///   - YBR (various) encodings.
 919 ///   - LUT[s] (or "PALETTE COLOR").
 920
 921 void PixelReadConvert::ConvertHandleColor()
 922 {
 923    //////////////////////////////////
 924    // Deal with the color decoding i.e. handle:
 925    //   - R, G, B planes (as opposed to RGB pixels)
 926    //   - YBR (various) encodings.
 927    //   - LUT[s] (or "PALETTE COLOR").
 928    //
 929    // The classification in the color decoding schema is based on the blending
 930    // of two Dicom tags values:
 931    // * "Photometric Interpretation" for which we have the cases:
 932    //  - [Photo A] MONOCHROME[1|2] pictures,
 933    //  - [Photo B] RGB or YBR_FULL_422 (which acts as RGB),
 934    //  - [Photo C] YBR_* (with the above exception of YBR_FULL_422)
 935    //  - [Photo D] "PALETTE COLOR" which indicates the presence of LUT[s].
 936    // * "Planar Configuration" for which we have the cases:
 937    //  - [Planar 0] 0 then Pixels are already RGB
 938    //  - [Planar 1] 1 then we have 3 planes : R, G, B,
 939    //  - [Planar 2] 2 then we have 1 gray Plane and 3 LUTs
 940    //
 941    // Now in theory, one could expect some coherence when blending the above
 942    // cases. For example we should not encounter files belonging at the
 943    // time to case [Planar 0] and case [Photo D].
 944    // Alas, this was only theory ! Because in practice some odd (read ill
 945    // formated Dicom) files (e.g. gdcmData/US-PAL-8-10x-echo.dcm) we encounter:
 946    //     - "Planar Configuration" = 0,
 947    //     - "Photometric Interpretation" = "PALETTE COLOR".
 948    // Hence gdcm will use the folowing "heuristic" in order to be tolerant
 949    // towards Dicom-non-conformance files:
 950    //   << whatever the "Planar Configuration" value might be, a
 951    //      "Photometric Interpretation" set to "PALETTE COLOR" forces
 952    //      a LUT intervention >>
 953    //
 954    // Now we are left with the following handling of the cases:
 955    // - [Planar 0] OR  [Photo A] no color decoding (since respectively
 956    //       Pixels are already RGB and monochrome pictures have no color :),
 957    // - [Planar 1] AND [Photo B] handled with ConvertRGBPlanesToRGBPixels()
 958    // - [Planar 1] AND [Photo C] handled with ConvertYcBcRPlanesToRGBPixels()
 959    // - [Planar 2] OR  [Photo D] requires LUT intervention.
 960
 961    if ( ! IsRawRGB() )
 962    {
 963       // [Planar 2] OR  [Photo D]: LUT intervention done outside
 964       return;
 965    }
 966
 967    if ( PlanarConfiguration == 1 )
 968    {
 969       if ( IsYBRFull )
 970       {
 971          // [Planar 1] AND [Photo C] (remember YBR_FULL_422 acts as RGB)
 972          ConvertYcBcRPlanesToRGBPixels();
 973       }
 974       else
 975       {
 976          // [Planar 1] AND [Photo C]
 977          ConvertRGBPlanesToRGBPixels();
 978       }
 979       return;
 980    }
 981
 982    // When planarConf is 0, and RLELossless (forbidden by Dicom norm)
 983    // pixels need to be RGB-fied anyway
 984    if (IsRLELossless)
 985    {
 986      ConvertRGBPlanesToRGBPixels();
 987    }
 988    // In *normal *case, when planarConf is 0, pixels are already in RGB
 989 }
 990
 991 /// Computes the Pixels Size
 992 void PixelReadConvert::ComputeRawAndRGBSizes()
 993 {
 994    int bitsAllocated = BitsAllocated;
 995    // Number of "Bits Allocated" is fixed to 16 when it's 12, since
 996    // in this case we will expand the image to 16 bits (see
 997    //    \ref ReadAndDecompress12BitsTo16Bits() )
 998    if (  BitsAllocated == 12 )
 999    {
1000       bitsAllocated = 16;
1001    }
1002
1003    RawSize =  XSize * YSize * ZSize
1004                      * ( bitsAllocated / 8 )
1005                      * SamplesPerPixel;
1006    if ( HasLUT )
1007    {
1008       RGBSize = 3 * RawSize;
1009    }
1010    else
1011    {
1012       RGBSize = RawSize;
1013    }
1014 }
1015
1016 /// Allocates room for RGB Pixels
1017 void PixelReadConvert::AllocateRGB()
1018 {
1019   if ( RGB )
1020      delete [] RGB;
1021   RGB = new uint8_t[RGBSize];
1022 }
1023
1024 /// Allocates room for RAW Pixels
1025 void PixelReadConvert::AllocateRaw()
1026 {
1027   if ( Raw )
1028      delete [] Raw;
1029   Raw = new uint8_t[RawSize];
1030 }
1031
1032 //-----------------------------------------------------------------------------
1033 // Print
1034 /**
1035  * \brief        Print self.
1036  * @param indent Indentation string to be prepended during printing.
1037  * @param os     Stream to print to.
1038  */
1039 void PixelReadConvert::Print( std::ostream &os, std::string const &indent )
1040 {
1041    os << indent
1042       << "--- Pixel information -------------------------"
1043       << std::endl;
1044    os << indent
1045       << "Pixel Data: offset " << PixelOffset
1046       << " x(" << std::hex << PixelOffset << std::dec
1047       << ")   length " << PixelDataLength
1048       << " x(" << std::hex << PixelDataLength << std::dec
1049       << ")" << std::endl;
1050
1051    if ( IsRLELossless )
1052    {
1053       if ( RLEInfo )
1054       {
1055          RLEInfo->Print( os, indent );
1056       }
1057       else
1058       {
1059          gdcmWarningMacro("Set as RLE file but NO RLEinfo present.");
1060       }
1061    }
1062
1063    if ( IsJPEG2000 || IsJPEGLossless || IsJPEGLossy || IsJPEGLS )
1064    {
1065       if ( JPEGInfo )
1066       {
1067          JPEGInfo->Print( os, indent );
1068       }
1069       else
1070       {
1071          gdcmWarningMacro("Set as JPEG file but NO JPEGinfo present.");
1072       }
1073    }
1074 }
1075
1076 //-----------------------------------------------------------------------------
1077 } // end namespace gdcm
1078
1079 // NOTES on File internal calls
1080 // User
1081 // ---> GetImageData
1082 //     ---> GetImageDataIntoVector
1083 //        |---> GetImageDataIntoVectorRaw
1084 //        | lut intervention
1085 // User
1086 // ---> GetImageDataRaw
1087 //     ---> GetImageDataIntoVectorRaw
1088