Example/exExtractCSATag.cxx

   1 /*=========================================================================
   2
   3   Program:   gdcm
   4   Module:    $RCSfile: exExtractCSATag.cxx,v $
   5   Language:  C++
   6   Date:      $Date: 2007/10/19 10:18:25 $
   7   Version:   $Revision: 1.3 $
   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 #include "gdcmFile.h"
  19 #include "gdcmFileHelper.h"
  20 #include "gdcmCommon.h"
  21 #include "gdcmDebug.h"
  22 #include "gdcmDocEntry.h"
  23 #include "gdcmDataEntry.h"
  24 #include "gdcmSeqEntry.h"
  25 #include "gdcmSQItem.h"
  26
  27 #include <iomanip>
  28
  29 namespace gdcm
  30 {
  31
  32    static const char *lookupTable1[] = {
  33       "UsedPatientWeight",
  34       "NumberOfPrescans",
  35       "TransmitterCalibration",
  36       "PhaseGradientAmplitude",
  37       "ReadoutGradientAmplitude",
  38       "SelectionGradientAmplitude",
  39       "GradientDelayTime",
  40       "RfWatchdogMask",
  41       "RfPowerErrorIndicator",
  42       "SarWholeBody",
  43       "Sed",
  44       "SequenceFileOwner",
  45       "Stim_mon_mode",
  46       "Operation_mode_flag",
  47       "dBdt_max",
  48       "t_puls_max",
  49       "dBdt_thresh",
  50       "dBdt_limit",
  51       "SW_korr_faktor",
  52       "Stim_max_online",
  53       "Stim_max_ges_norm_online",
  54       "Stim_lim",
  55       "Stim_faktor",
  56       "CoilForGradient",
  57       "CoilTuningReflection",
  58       "CoilId",
  59       "MiscSequenceParam",
  60       "MrProtocolVersion",
  61       "MrProtocol",
  62       "DataFileName",
  63       "RepresentativeImage",
  64       "PositivePCSDirections",
  65       "RelTablePosition",
  66       "ReadoutOS",
  67       "LongModelName",
  68       "SliceArrayConcatenations",
  69       "SliceResolution",
  70       "MrEvaProtocol",
  71       "AbsTablePosition",
  72       "AutoAlignMatrix",
  73       "MeasurementIndex",
  74       "CoilString",
  75       "PATModeText",
  76       "PatReinPattern",
  77       NULL
  78    };
  79
  80    static const char *lookupTable2[] = {
  81       "EchoLinePosition",
  82       "EchoColumnPosition",
  83       "EchoPartitionPosition",
  84       "UsedChannelMask",
  85       "Actual3DImaPartNumber",
  86       "ICE_Dims",
  87       "B_value",
  88       "Filter1",
  89       "Filter2",
  90       "ProtocolSliceNumber",
  91       "RealDwellTime",
  92       "PixelFile",
  93       "PixelFileName",
  94       "SliceMeasurementDuration",
  95       "SequenceMask",
  96       "AcquisitionMatrixText",
  97       "MeasuredFourierLines",
  98       "FlowEncodingDirection",
  99       "FlowVenc",
 100       "PhaseEncodingDirectionPositive",
 101       "NumberOfImagesInMosaic",
 102       "DiffusionGradientDirection",
 103       "ImageGroup",
 104       "SliceNormalVector",
 105       "DiffusionDirectionality",
 106       "TimeAfterStart",
 107       "FlipAngle",
 108       "SequenceName",
 109       "RepetitionTime",
 110       "EchoTime",
 111       "NumberOfAverages",
 112       "VoxelThickness",
 113       "VoxelPhaseFOV",
 114       "VoxelReadoutFOV",
 115       "VoxelPositionSag",
 116       "VoxelPositionCor",
 117       "VoxelPositionTra",
 118       "VoxelNormalSag",
 119       "VoxelNormalCor",
 120       "VoxelNormalTra",
 121       "VoxelInPlaneRot",
 122       "ImagePositionPatient",
 123       "ImageOrientationPatient",
 124       "PixelSpacing",
 125       "SliceLocation",
 126       "SliceThickness",
 127       "SpectrumTextRegionLabel",
 128       "Comp_Algorithm",
 129       "Comp_Blended",
 130       "Comp_ManualAdjusted",
 131       "Comp_AutoParam",
 132       "Comp_AdjustedParam",
 133       "Comp_JobID",
 134       "FMRIStimulInfo",
 135       "FlowEncodingDirectionString",
 136       "RepetitionTimeEffective",
 137     NULL
 138    };
 139
 140    static const char *lookupTable3[] = {
 141       "ImageNumber",
 142       "ImageComments",
 143       "ReferencedImageSequence",
 144       "PatientOrientation",
 145       "ScanningSequence",
 146       "SequenceName",
 147       "RepetitionTime",
 148       "EchoTime",
 149       "InversionTime",
 150       "NumberOfAverages",
 151       "ImagingFrequency",
 152       "ImagedNucleus",
 153       "EchoNumbers",
 154       "MagneticFieldStrength",
 155       "NumberOfPhaseEncodingSteps",
 156       "EchoTrainLength",
 157       "PercentSampling",
 158       "PercentPhaseFieldOfView",
 159       "TriggerTime",
 160       "ReceivingCoil",
 161       "TransmittingCoil",
 162       "AcqusitionMatrix",
 163       "PhaseEncodingDirection",
 164       "FlipAngle",
 165       "VariableFlipAngleFlag",
 166       "SAR",
 167       "dBdt",
 168       "Rows",
 169       "Columns",
 170       "SliceThickness",
 171       "ImagePositionPatient",
 172       "ImageOrientationPatient",
 173       "SliceLocation",
 174       "EchoLinePosition",
 175       "EchoColumnPosition",
 176       "EchoPartitionPosition",
 177       "Actual3DImaPartNumber",
 178       "RealDwellTime",
 179       "ProtocolSliceNumber",
 180       "DataFile",
 181       "DataFileName",
 182       "ICE_Dims",
 183       "PixelSpacing",
 184       "SourceImageSequence",
 185       "PixelBandwidth",
 186       "SliceMeasurementDuration",
 187       "SequenceMask",
 188       "AcquisitionMatrixText",
 189       "MeasuredFourierLines",
 190       "CsiGridshiftVector",
 191       NULL
 192    };
 193
 194    /*
 195     * What if SIEMENS decide to add another entry in this table, all the offsets are completely off ...
 196     * Also we do not respect the ordering anymore.
 197     * TODO: Need to check if elements always comes in the same order
 198     */
 199    unsigned int GetLookupCSAIndex(const char *csa_name, const char **lookuptable = lookupTable1)
 200       {
 201      const char **p = lookuptable;
 202       while( *p && strcmp(*p, csa_name) != 0 )
 203          {
 204       ++p;
 205          }
 206     assert( strcmp(*p, csa_name) == 0 );
 207     return p - lookuptable + 1;  // Start at 1 to avoid being on the 0000 position
 208       }
 209
 210       // Looks like there is mapping in between syngodt and vr...
 211       // O <=> UN
 212       // 3 <=> DS
 213       // 4 <=> FD
 214       // 5 <=> FL
 215       // 6 <=> IS
 216       // 9 <=> UL
 217       // 10 <=> US
 218       // 16 <=> CS
 219       // 19 <=> LO
 220       // 20 <=> LT
 221       // 22 <=> SH
 222       // 25 <=> UI
 223
 224 struct equ
 225 {
 226    uint32_t syngodt;
 227    char vr[2+1];
 228 };
 229
 230 static equ mapping[] = {
 231    {     0 , "UN" },
 232    {     3 , "DS" },
 233    {     4 , "FD" },
 234    {     5 , "FL" },
 235    {     6 , "IS" },
 236    {     7 , "SL" },
 237    {     8 , "SS" },
 238    {     9 , "UL" },
 239    {    10 , "US" },
 240    {    16 , "CS" },
 241    {    19 , "LO" },
 242    {    20 , "LT" },
 243    {    22 , "SH" },
 244    {    23 , "ST" },
 245    {    25 , "UI" },
 246    {    27 , "UT" }
 247 };
 248
 249 bool check_mapping(uint32_t syngodt, const char *vr)
 250 {
 251   static const unsigned int max = sizeof(mapping) / sizeof(equ);
 252   unsigned int s = 0;
 253   const equ *p = mapping;
 254   assert( syngodt <= mapping[max-1].syngodt );
 255    while(p->syngodt < syngodt )
 256       {
 257     ++p;
 258       }
 259   assert( p->syngodt == syngodt ); // or else need to update mapping
 260   const char* lvr = p->vr;
 261   int check = strcmp(vr, lvr) == 0;
 262   assert( check );
 263   return true;
 264 }
 265
 266 uint32_t readCSAHeader(std::istream &is)
 267 {
 268    char dummy[4+1];
 269    dummy[4] = 0;
 270    is.read(dummy, 4);
 271    std::cout << dummy << std::endl;
 272    if( strcmp( dummy, "SV10" )  )
 273       {
 274       std::cerr << "Either not a SV10 header or filled with 0..." << std::endl;
 275       return 1;
 276       }
 277    // wotsit ?
 278    is.read(dummy, 4);
 279    if( strcmp( dummy, "\4\3\2\1" )  )
 280       {
 281       std::cerr << "Either not a SV10 header or filled with 0..." << std::endl;
 282       return 1;
 283       }
 284    std::cout << dummy << std::endl;
 285    uint32_t n;
 286    is.read((char*)&n, sizeof(n));
 287    std::cout << n << std::endl;
 288    uint32_t unused;
 289    is.read((char*)&unused, sizeof(unused));
 290    std::cout << unused << std::endl;
 291    assert( unused == 77 ); // 'M' character...
 292
 293   return n;
 294 }
 295
 296 DataEntry *readCSAElement(std::istream &is)
 297 {
 298       char name[64+1];
 299       name[64] = 0; // security
 300       //std::cout << "Pos 0x" << std::hex << is.tellg() << std::dec << std::endl;
 301       is.read(name, 64);
 302       std::cout << "Name=" << name << std::endl;
 303      unsigned int element = GetLookupCSAIndex(name,lookupTable2);
 304       uint32_t vm;
 305       is.read((char*)&vm, sizeof(vm));
 306       std::cout << "vm=" << vm << std::endl;
 307       char vr[4];
 308       is.read(vr, 4);
 309     assert( vr[2] == vr[3] && vr[2] == 0 );
 310       std::cout << "vr=" << vr << std::endl;
 311
 312       DataEntry *de = DataEntry::New(0x0029, element, vr);
 313
 314       uint32_t syngodt;
 315       is.read((char*)&syngodt, sizeof(syngodt));
 316       check_mapping(syngodt, vr);
 317
 318       std::cout << "syngodt=" << syngodt << std::endl;
 319       uint32_t nitems;
 320       is.read((char*)&nitems, sizeof(nitems));
 321       std::cout << "nitems=" << nitems<< std::endl;
 322       uint32_t xx;
 323       is.read((char*)&xx, sizeof(xx));
 324       //std::cout << "xx=" << xx<< std::endl;
 325       assert( xx == 77 || xx == 205 );
 326       for( uint32_t j = 0; j < nitems; ++j)
 327          {
 328          uint32_t item_xx[4];
 329          is.read((char*)&item_xx, 4*sizeof(uint32_t));
 330          std::cout << "item_xx=" << item_xx[0] << " " << item_xx[1] << " " << item_xx[2] << " " << item_xx[3] << std::endl;
 331          //std::cout << "0x" << std::hex << is.tellg() << std::dec << std::endl;
 332          assert( item_xx[2] == 77 || item_xx[2] == 205 );
 333          uint32_t len = item_xx[1]; // 2nd element
 334          std::cout << "len=" << len << std::endl;
 335          assert( item_xx[0] == item_xx[1] && item_xx[1] == item_xx[3] );
 336          char *val = new char[len+1];
 337          val[len] = 0; // security
 338          is.read(val,len);
 339          // WARNING vr does not means anything AFAIK, simply print the value as if it was IS/DS or LO (ASCII)
 340          std::cout << "val=" << val << std::endl;
 341       if( !j )
 342       de->SetString( std::string(val,len) );
 343          delete[] val;
 344          char dummy[4];
 345          uint32_t dummy_len = (4 - len % 4) % 4;
 346          is.read(dummy, dummy_len );
 347          for(uint32_t d= 0; d < dummy_len; ++d)
 348             {
 349             // I think dummy should always be 0
 350             if( dummy[d] )
 351                {
 352                std::cout << "dummy=" << (int)dummy[d] << std::endl;
 353                }
 354             }
 355          }
 356
 357
 358   return de;
 359 }
 360
 361 int convertCSA(std::istream &is, File *f)
 362 {
 363    f->RemoveEntry( f->GetDataEntry(0X0029,0x1010) );
 364
 365     SeqEntry *sq = SeqEntry::New(0x0029,0x1010);
 366     SQItem *sqi = SQItem::New(1);
 367 //    DataEntry *e_0008_1150 = DataEntry::New(0x0008, 0x1150, "UI");
 368 //    e_0008_1150->SetString( "coucou" );
 369 //    sqi->AddEntry(e_0008_1150);
 370 //    e_0008_1150->Delete();
 371 //
 372 //    DataEntry *e_0008_1155 = DataEntry::New(0x0008, 0x1155, "UI");
 373 //    e_0008_1155->SetString( "mathieu" );
 374 //    sqi->AddEntry(e_0008_1155);
 375 //    e_0008_1155->Delete();
 376
 377
 378
 379   uint32_t n = readCSAHeader(is);
 380
 381    for(uint32_t i = 0; i < n; ++i)
 382       {
 383     DataEntry *de = readCSAElement(is);
 384     sqi->AddEntry(de);
 385     de->Delete();
 386       }
 387
 388     sq->AddSQItem(sqi,1);
 389     sqi->Delete();
 390     //sq->Print( std::cout );
 391    f->AddEntry( sq );
 392     sq->Delete();
 393    f->Print( std::cout );
 394
 395    GDCM_NAME_SPACE::FileHelper *fh = GDCM_NAME_SPACE::FileHelper::New(f);
 396          fh->SetWriteTypeToDcmExplVR();
 397          fh->Write("/tmp/csa2.dcm");
 398
 399    return 0;
 400 }
 401 } // end namespace gdcm
 402
 403 int main(int argc, char *argv[])
 404 {
 405    GDCM_NAME_SPACE::File *f;
 406
 407    if( argc < 5 )
 408    {
 409       std::cerr << "Usage :" << argv[0] << " input.dcm  group element outputfile" << std::endl;
 410       std::cerr << "  Ex: " << argv[0] << " /tmp/bla.dcm 0029 2110 /tmp/out.raw" << std::endl;
 411       return 1;
 412    }
 413    std::string fileName = argv[1];
 414
 415    std::cout << fileName << std::endl;
 416 // ============================================================
 417 //   Read the input image.
 418 // ============================================================
 419
 420    f = GDCM_NAME_SPACE::File::New( );
 421
 422    //f->SetLoadMode(GDCM_NAME_SPACE::LD_NOSEQ | GDCM_NAME_SPACE::LD_NOSHADOW);
 423    f->SetFileName( fileName );
 424    f->SetMaxSizeLoadEntry(0xffff);
 425    bool res = f->Load();
 426
 427    if( GDCM_NAME_SPACE::Debug::GetDebugFlag())
 428    {
 429       std::cout << "---------------------------------------------" << std::endl;
 430       f->Print();
 431       std::cout << "---------------------------------------------" << std::endl;
 432    }
 433    if (!res) {
 434        std::cerr << "Sorry, " << fileName << " not a gdcm-readable "
 435            << "DICOM / ACR File"
 436            << std::endl;
 437       f->Delete();
 438       return 1;
 439    }
 440    std::cout << " ... is readable " << std::endl;
 441
 442    // Find the dicom tag, and extract the string
 443    uint16_t group, elem;
 444    std::istringstream convert;
 445    convert.str( argv[2] );
 446    convert >> std::hex >> group;
 447    convert.clear(); //important
 448    convert.str( argv[3] );
 449    convert >> std::hex >> elem;
 450    std::cout << "Extracting tag: (0x" << std::hex << std::setw(4) << std::setfill('0')
 451      << group << ",0x" << std::setw(4) << std::setfill('0') << elem << ")" << std::endl;
 452    std::string dicom_tag_value = f->GetEntryString(group, elem);
 453    if (dicom_tag_value == GDCM_NAME_SPACE::GDCM_UNFOUND)
 454    {
 455      GDCM_NAME_SPACE::DictEntry *dictEntry = f->GetPubDict()->GetEntry( group, elem);
 456      std::cerr << "Image doesn't contain any tag: " << dictEntry->GetName() << std::endl;
 457      f->Delete();
 458      return 1;
 459    }
 460
 461    GDCM_NAME_SPACE::DocEntry *dicom_tag_doc = f->GetDocEntry(group, elem);
 462    GDCM_NAME_SPACE::DataEntry *dicom_tag = dynamic_cast<GDCM_NAME_SPACE::DataEntry *>(dicom_tag_doc);
 463    if( !dicom_tag )
 464    {
 465       std::cerr << "Sorry DataEntry only please" << std::endl;
 466       f->Delete();
 467       return 1;
 468    }
 469
 470    // Write out the data as a file:
 471    std::ofstream o(argv[4]);
 472    if( !o )
 473    {
 474       std::cerr << "Problem opening file: " << argv[4] << std::endl;
 475       f->Delete();
 476       return 1;
 477    }
 478    o.write((char*)dicom_tag->GetBinArea(), dicom_tag->GetLength());
 479    o.close();
 480
 481    std::istringstream is;
 482    is.str( std::string( (char*)dicom_tag->GetBinArea(), dicom_tag->GetLength()) );
 483
 484    GDCM_NAME_SPACE::convertCSA(is, f);
 485
 486    f->Delete();
 487
 488
 489
 490    return 0;
 491 }
 492
 493