src/gdcmSegmentedPalette.h

   1 /*=========================================================================
   2
   3   Program:   gdcm
   4   Module:    $RCSfile: gdcmSegmentedPalette.h,v $
   5   Language:  C++
   6   Date:      $Date: 2007/10/08 15:20:17 $
   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 #ifndef _GDCMSEGMENTEDPALETTE_H_
  20 #define _GDCMSEGMENTEDPALETTE_H_
  21
  22 #include "gdcmFile.h"
  23 #include "gdcmTagKey.h"
  24 #include "gdcmDataEntry.h"
  25
  26 // Ref:
  27 // http://blog.goo.ne.jp/satomi_takeo/e/3643e5249b2a9650f9e10ef1c830e8b8
  28 // I bet the code was compiled on VS6. Make it compile on other platform:
  29 // * typedef are not inherited
  30 // * need to explicitely add typename keyword
  31 // * Uint8 / Uint16 are neither C nor C++
  32 // * replace all dcmtk code with equivalent gdcm code
  33 // * Add extra parameter to ReadPaletteInto to save extracted info
  34
  35 // MM: Over 50000 DICOM files I only found two that had a Segmented Palette LUT
  36 // And both were coming from a ALOKA SSD-4000 station
  37 #include <assert.h>
  38 #include <algorithm>
  39 #include <deque>
  40 #include <map>
  41 #include <vector>
  42 #include <iterator>
  43
  44 // Hack for VS6
  45 #if defined(_MSC_VER) && (_MSC_VER < 1310)
  46 #define GDCM_TYPENAME class
  47 #else
  48 #define GDCM_TYPENAME typename
  49 #endif
  50
  51 namespace GDCM_NAME_SPACE
  52 {
  53     // abstract class for segment.
  54     template <GDCM_TYPENAME EntryType>
  55     class Segment {
  56     public:
  57         typedef std::map<const EntryType*, const Segment*> SegmentMap;
  58         virtual bool Expand(const SegmentMap& instances,
  59             std::vector<EntryType>& expanded) const = 0;
  60         const EntryType* First() const { return _first; }
  61         const EntryType* Last() const { return _last; }
  62         struct ToMap {
  63             std::pair<
  64                 typename SegmentMap::key_type,
  65                 typename SegmentMap::mapped_type
  66             >
  67                 operator()(const Segment* segment) const
  68             { return std::make_pair(segment->First(), segment); }
  69         };
  70     protected:
  71         Segment(const EntryType* first, const EntryType* last) {
  72             _first = first; _last = last;
  73         }
  74         const EntryType* _first;
  75         const EntryType* _last;
  76     };
  77
  78     // discrete segment (opcode = 0)
  79     template <GDCM_TYPENAME EntryType>
  80     class DiscreteSegment : public Segment<EntryType> {
  81     public:
  82         typedef typename Segment<EntryType>::SegmentMap SegmentMap;
  83         DiscreteSegment(const EntryType* first)
  84             : Segment<EntryType>(first, first+2+*(first+1)) {}
  85         virtual bool Expand(const SegmentMap&,
  86             std::vector<EntryType>& expanded) const
  87         {
  88             std::copy(this->_first + 2, this->_last, std::back_inserter(expanded));
  89             return true;
  90         }
  91     };
  92
  93     // linear segment (opcode = 1)
  94     template <GDCM_TYPENAME EntryType>
  95     class LinearSegment : public Segment<EntryType> {
  96     public:
  97         typedef typename Segment<EntryType>::SegmentMap SegmentMap;
  98         LinearSegment(const EntryType* first)
  99             : Segment<EntryType>(first, first+3) {}
 100         virtual bool Expand(const SegmentMap&,
 101             std::vector<EntryType>& expanded) const
 102         {
 103             if ( expanded.empty() ) {
 104                 // linear segment can't be the first segment.
 105                 return false;
 106             }
 107             EntryType length = *(this->_first + 1);
 108             EntryType y0 = expanded.back();
 109             EntryType y1 = *(this->_first + 2);
 110             double y01 = y1 - y0;
 111             for ( EntryType i = 0; i <length; ++i ) {
 112                 double value_float
 113                     = static_cast<double>(y0)
 114                     + (static_cast<double>(i)/static_cast<double>(length)) * y01;
 115                 EntryType value_int = static_cast<EntryType>(value_float + 0.5);
 116                 expanded.push_back(value_int);
 117             }
 118             return true;
 119         }
 120     };
 121
 122     // indirect segment (opcode = 2)
 123     template <GDCM_TYPENAME EntryType>
 124     class IndirectSegment : public Segment<EntryType> {
 125     public:
 126         typedef typename Segment<EntryType>::SegmentMap SegmentMap;
 127         IndirectSegment(const EntryType* first)
 128             : Segment<EntryType>(first, first+2+4/sizeof(EntryType)) {}
 129         virtual bool Expand(const SegmentMap& instances,
 130             std::vector<EntryType>& expanded) const
 131         {
 132             if ( instances.empty() ) {
 133                 // some other segments are required as references.
 134                 return false;
 135             }
 136             const EntryType* first_segment = instances.begin()->first;
 137             const unsigned short* pOffset
 138                 = reinterpret_cast<const unsigned short*>(this->_first + 2);
 139             unsigned long offsetBytes
 140                 = (*pOffset) | (static_cast<unsigned long>(*(pOffset + 1)) << 16);
 141             const EntryType* copied_part_head
 142                 = first_segment + offsetBytes / sizeof(EntryType);
 143             typename SegmentMap::const_iterator ppHeadSeg = instances.find(copied_part_head);
 144             if ( ppHeadSeg == instances.end() ) {
 145                 // referred segment not found
 146                 return false;
 147             }
 148             EntryType nNumCopies = *(this->_first + 1);
 149             typename SegmentMap::const_iterator ppSeg = ppHeadSeg;
 150             while ( std::distance(ppHeadSeg, ppSeg) < nNumCopies ) {
 151                 assert( ppSeg != instances.end() );
 152                 ppSeg->second->Expand(instances, expanded);
 153                 ++ppSeg;
 154             }
 155             return true;
 156         }
 157     };
 158
 159     template <GDCM_TYPENAME EntryType>
 160     void ExpandPalette(const EntryType* raw_values, uint32_t length,
 161         std::vector<EntryType>& palette)
 162     {
 163         typedef std::deque<Segment<EntryType>*> SegmentList;
 164         SegmentList segments;
 165         const EntryType* raw_seg = raw_values;
 166         while ( (std::distance(raw_values, raw_seg) * sizeof(EntryType)) <length ) {
 167             Segment<EntryType>* segment = NULL;
 168             if ( *raw_seg == 0 ) {
 169                 segment = new DiscreteSegment<EntryType>(raw_seg);
 170             } else if ( *raw_seg == 1 ) {
 171                 segment = new LinearSegment<EntryType>(raw_seg);
 172             } else if ( *raw_seg == 2 ) {
 173                 segment = new IndirectSegment<EntryType>(raw_seg);
 174             }
 175             if ( segment ) {
 176                 segments.push_back(segment);
 177                 raw_seg = segment->Last();
 178             } else {
 179                 // invalid opcode
 180                 break;
 181             }
 182         }
 183         typename Segment<EntryType>::SegmentMap instances;
 184         std::transform(segments.begin(), segments.end(),
 185             std::inserter(instances, instances.end()), typename Segment<EntryType>::ToMap());
 186         typename SegmentList::iterator ppSeg = segments.begin();
 187         typename SegmentList::iterator endOfSegments = segments.end();
 188         for ( ; ppSeg != endOfSegments; ++ppSeg ) {
 189             (*ppSeg)->Expand(instances, palette);
 190         }
 191         ppSeg = segments.begin();
 192         for ( ; ppSeg != endOfSegments; ++ppSeg ) {
 193             delete *ppSeg;
 194         }
 195     }
 196
 197     void ReadPaletteInto(GDCM_NAME_SPACE::File* pds, const GDCM_NAME_SPACE::TagKey& descriptor,
 198       const GDCM_NAME_SPACE::TagKey& segment, uint8_t* lut)
 199       {
 200       unsigned int desc_values[3] = {};
 201       unsigned long count = 0;
 202       //if ( pds->findAndGetUint16Array(descriptor, desc_values, &count).good() )
 203       std::string desc_values_str = pds->GetEntryString(descriptor.GetGroup(), descriptor.GetElement() );
 204       count = sscanf( desc_values_str.c_str(), "%u\\%u\\%u", desc_values, desc_values+1, desc_values+2 );
 205         {
 206         assert( count == 3 );
 207         unsigned int num_entries = desc_values[0];
 208         if ( num_entries == 0 ) {
 209           num_entries = 0x10000;
 210         }
 211         unsigned int min_pixel_value = desc_values[1];
 212         assert( min_pixel_value == 0 ); // FIXME
 213         unsigned int entry_size = desc_values[2];
 214         assert( entry_size == 8 || entry_size == 16 );
 215         //DcmElement* pe = NULL;
 216         GDCM_NAME_SPACE::DataEntry* pe = NULL;
 217
 218         pe = pds->GetDataEntry(segment.GetGroup(), segment.GetElement() );
 219           {
 220           //if ( pds->findAndGetElement(segment, pe).good() )
 221           unsigned long length = pe->GetLength();
 222           if ( entry_size == 8 )
 223             {
 224             uint8_t* segment_values = NULL;
 225             //if ( pe->getUint8Array(segment_values).good() )
 226             segment_values = (uint8_t*)pe->GetBinArea();
 227               {
 228               std::vector<uint8_t> palette;
 229               palette.reserve(num_entries);
 230               ExpandPalette(segment_values, length, palette);
 231               memcpy(lut, &palette[0], palette.size() );
 232               }
 233             }
 234           else if ( entry_size == 16 )
 235             {
 236             uint16_t* segment_values = NULL;
 237             segment_values = (uint16_t*)pe->GetBinArea();
 238               {
 239               //if ( pe->getUint16Array(segment_values).good() )
 240               std::vector<uint16_t> palette;
 241               palette.reserve(num_entries);
 242               ExpandPalette(segment_values, length, palette);
 243               memcpy(lut, &palette[0], palette.size()*2 );
 244
 245               }
 246             }
 247           }
 248         }
 249       }
 250 } // end namespace gdcm
 251
 252
 253 // do not pollute namespace:
 254 #undef GDCM_TYPENAME
 255
 256 #endif