src/gdcmPixelReadConvert.cxx

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