#include "itkGDCMImporter.h"

#include <itkGDCMImageIOFactory.h>
#include <itkGDCMImageIO.h>

#include <itkImageFileWriter.h>

#include <sstream>


namespace itk
{
  
  GDCMImporter::GDCMImporter()
  {
    itk::GDCMImageIOFactory::RegisterOneFactory();
    
    this->DICOMScanner = itk::GDCMSeriesFileNames::New();
    this->DICOMScanner->SetUseSeriesDetails(true);
    
    this->InputDirectory = "";
    this->m_Recursive = true;
  }


  
  GDCMImporter::~GDCMImporter()
  {

    // this->Reset(); // no need to reset when object is deleted
    
  }

  
  uint16_t GDCMImporter::GetDICOMTagGroup(unsigned long tag)
  {
    switch(tag)
    {
	case D_MODALITY :
	  return 0x0008;
	case D_SERIESUID :
	  return 0x0020;
	case D_STUDYDESCRIPTION :
	  return 0x0008;
	case D_SERIESDESCRIPTION :
	  return 0x0008;
	case D_ACQUISITIONTIME :
	  return 0x0008;
	default :
	  return 0x0020;
    }
  }

  
  uint16_t GDCMImporter::GetDICOMTagElement(unsigned long tag)
  {
    switch(tag)
    {
	case D_MODALITY :
	  return 0x0060;
	case D_SERIESUID :
	  return 0x000e;
	case D_STUDYDESCRIPTION :
	  return 0x1030;
	case D_SERIESDESCRIPTION :
	  return 0x103e;
	case D_ACQUISITIONTIME :
	  return 0x0032;
	default :
	  return 0x0020;
    }
  }
  

  bool GDCMImporter::CanReadFile (const char* filename)
  {
    ReaderType::Pointer reader = ReaderType::New();
    reader->SetFileName(filename);
    try
    {
      reader->Update();
      return true;
    }
    catch (itk::ExceptionObject & e)
    {
      std::cerr << e;
      return false;
    }
  }


  void GDCMImporter::Reset(void)
  {
    this->DICOMVolumeList.clear();
    this->OutputVolumeList.clear();
  }


  void GDCMImporter::Scan (void)
  {
    this->UpdateProgress(0.05);
    
    if (!this->InputDirectory.size())
    {
      itkWarningMacro(<<"Input directory not set !"<<std::endl<<"No scan performed."<<std::endl);
      return;
    }

    this->Reset();
    
    this->UpdateProgress(0.1);
    
    this->DICOMScanner->SetRecursive( this->GetRecursive() );
	
    this->DICOMScanner->SetDirectory( this->InputDirectory );

    std::vector<std::string> UIDs = this->DICOMScanner->GetSeriesUIDs();
    
    this->UpdateProgress(0.2);

    std::vector<std::string> totalfilenames = this->DICOMScanner->GetInputFileNames();
    this->UpdateProgress(0.3);

    unsigned int total_size = totalfilenames.size();
    unsigned int total_iterator = 0;

    this->UpdateProgress(0.0);

    for (unsigned int i=0; i<UIDs.size(); i++)
    {

      std::string uid = UIDs[i];

      
      std::vector<std::string> filenamelist  = this->DICOMScanner->GetFileNames(uid);
      std::vector<gdcm::File*>* gdcmfilelist = this->DICOMScanner->GetSeriesHelper()->GetSingleSerieUIDFileSet(uid);

        
      std::ostringstream os;
      
      if( gdcmfilelist->size()>0 )
      {
	std::string date = (*gdcmfilelist)[0]->GetEntryValue(0x0008,0x0021);
	if (date == gdcm::GDCM_UNFOUND || date=="")
        {
          date = "nodate";
        }

        std::string sequenceName = (*gdcmfilelist)[0]->GetEntryValue(0x0018, 0x0024);
        if( sequenceName == gdcm::GDCM_UNFOUND || sequenceName == "")
        {
          char number[256];
          sprintf (number, "%d", i);
          sequenceName = number;
        }
        
        
        os << date << "_" << (*gdcmfilelist)[0]->GetEntryValue(0x0008,0x0060) << "_" << sequenceName ;
        
      }
      
      bool Is2DFiles = false;
      for (unsigned int file_it=0; file_it<gdcmfilelist->size(); file_it++)
      {
	int x = (*gdcmfilelist)[file_it]->GetXSize();
	int y = (*gdcmfilelist)[file_it]->GetYSize();
	int z = (*gdcmfilelist)[file_it]->GetZSize();

	Is2DFiles = (x == 1 || y == 1 || z == 1);
	if (!Is2DFiles)
	  break;
      }
      if (!Is2DFiles && gdcmfilelist->size())
      {
	// cannot construct from 3D files
	itkWarningMacro(<<std::endl<<"file : "<<(*gdcmfilelist)[0]->GetFileName().c_str()
			<<std::endl<<"cannot import 3D (or 2D+t) DICOM files !"
			<<std::endl<<"volume "<<os.str()<<" will not be imported"<<std::endl);

	continue;
      }

      this->AddVolume(os.str(), filenamelist, gdcmfilelist);
      this->UpdateProgress((double)i/(double)UIDs.size());
      
    }    
    
  }


  void GDCMImporter::SplitVolumeByTag (DICOMVolume::Pointer inputvolume, unsigned long tag)
  {
    
    uint16_t group   = GDCMImporter::GetDICOMTagGroup(tag);
    uint16_t element = GDCMImporter::GetDICOMTagElement(tag);

    gdcm::XCoherentFileSetmap filesetmap = this->DICOMScanner->GetSeriesHelper()->SplitOnTagValue(inputvolume->GetgdcmFileList(), group, element);
    gdcm::XCoherentFileSetmap::iterator it;

    std::vector<itk::DICOMVolume::Pointer> listToInsert;
    
    for (it = filesetmap.begin(); it != filesetmap.end(); it++)
    {
      std::string splittag = (*it).first;
      if (splittag == gdcm::GDCM_UNFOUND)
	splittag = "notag";
      
      gdcm::FileList* gdcmfilelist = (*it).second;

      if (!gdcmfilelist->size())
      {
        continue; // What is it supposed to do here?
	          // Sometime when the input volume doesn't have to be splitted,
	          // we end up with a empty secondary file list (size = 0),
	          // if so then we don't want insert an empty-file-list volume, do we ?
	          // this is what it is supposed to do
      }
      
      itk::DICOMVolume::Pointer dicomvolume = itk::DICOMVolume::New();
      std::string description = inputvolume->GetDescription()+"_"+splittag;
      dicomvolume->SetDescription (description);
      dicomvolume->SetgdcmFileList (gdcmfilelist);
      dicomvolume->SetSerieHelper(this->DICOMScanner->GetSeriesHelper());

      std::vector<std::string> filelist;
      gdcm::FileList::iterator it2;
      for(it2 = gdcmfilelist->begin(); it2 != gdcmfilelist->end(); it2++ )
      {
        filelist.push_back((*it2)->GetFileName());
      }
      	
      dicomvolume->SetFileList (filelist);

      listToInsert.push_back (dicomvolume);
    }


    // now insert all elements of listToInsert into this->DicomVolume
    std::vector<itk::DICOMVolume::Pointer> testList (1);
    testList[0] = inputvolume;

    std::vector<itk::DICOMVolume::Pointer>::iterator toInsert = search (this->DICOMVolumeList.begin(), this->DICOMVolumeList.end(),
                                                                        testList.begin(), testList.end());
    if( toInsert == this->DICOMVolumeList.end() )
    {
      std::cerr << "Error: Cannot find dicom volume in input list." << std::endl;
      return;
    }
    
    this->DICOMVolumeList.insert (++toInsert, listToInsert.begin(), listToInsert.end());
    
    this->RemoveDICOMVolume(inputvolume);
  }



  

  void GDCMImporter::SplitVolumeByPositionConsistency (DICOMVolume::Pointer inputvolume)
  {
    gdcm::XCoherentFileSetmap filesetmap = this->SplitOnPosition(inputvolume->GetgdcmFileList());
    gdcm::XCoherentFileSetmap newfilesetmap;
    
    bool multiple = false;

    gdcm::XCoherentFileSetmap::iterator it = filesetmap.begin();
    while( it != filesetmap.end() )
    {
      gdcm::FileList* gdcmfilelist = (*it).second;
      
      if (!gdcmfilelist->size())
      {
        continue; // What is it supposed to do here?
	          // Same as above : 
	          // Sometime when the input volume doesn't have to be splitted,
	          // we end up with a empty secondary file list (size = 0),
	          // if so then we don't want insert an empty-file-list volume, do we ?
	          // this is what it is supposed to do
      }
      
      unsigned int fileid = 0;
      gdcm::FileList::iterator it2 = gdcmfilelist->begin();
      while( it2 != gdcmfilelist->end() )
      {
        
	std::ostringstream strid;
	strid << inputvolume->GetDescription().c_str()<<"_item_"<<fileid;
        
	if (!newfilesetmap.count (strid.str()))
	{
	  newfilesetmap[strid.str()] = new gdcm::FileList;
	}
	
	newfilesetmap[strid.str()]->push_back( (*it2) );
	
	fileid++;

	if (fileid > 1)
        {
          multiple = true;
        }

        ++it2;
      }

      ++it;
    }


    
    std::vector<itk::DICOMVolume::Pointer> listToInsert;
    
    for (it = newfilesetmap.begin(); it != newfilesetmap.end(); it++)
    {
      gdcm::FileList* gdcmfilelist = (*it).second;
      
      if (!gdcmfilelist->size())
	continue;
      
      
      itk::DICOMVolume::Pointer dicomvolume = itk::DICOMVolume::New();
      std::string description = (*it).first;
      
      dicomvolume->SetDescription (description);
      dicomvolume->SetgdcmFileList (gdcmfilelist);
      dicomvolume->SetSerieHelper(this->DICOMScanner->GetSeriesHelper());

      gdcm::FileList::iterator it2;

      this->DICOMScanner->GetSeriesHelper()->OrderFileList( gdcmfilelist );

      std::vector<std::string> filelist;
      for(it2 = gdcmfilelist->begin(); it2 != gdcmfilelist->end(); it2++ )
      {
        filelist.push_back((*it2)->GetFileName());
      }
      
      dicomvolume->SetFileList (filelist);

      listToInsert.push_back (dicomvolume);
    }


    std::vector<itk::DICOMVolume::Pointer> testList (1);
    testList[0] = inputvolume;
    
    std::vector<itk::DICOMVolume::Pointer>::iterator toInsert = search (this->DICOMVolumeList.begin(), this->DICOMVolumeList.end(),
                                                                        testList.begin(), testList.end());

    if( toInsert == this->DICOMVolumeList.end() )
    {
      std::cerr << "Error: Cannot find dicom volume in input list." << std::endl;
      return;
    }
    
    this->DICOMVolumeList.insert (++toInsert, listToInsert.begin(), listToInsert.end());
    
    this->RemoveDICOMVolume(inputvolume);
  }




  void GDCMImporter::RemoveDICOMVolume (DICOMVolume::Pointer dcmvolume)
  {
    std::vector<itk::DICOMVolume::Pointer> testList (1);
    testList[0] = dcmvolume;
    
    std::vector<itk::DICOMVolume::Pointer>::iterator toRemove = search (this->DICOMVolumeList.begin(), this->DICOMVolumeList.end(),
                                                                        testList.begin(), testList.end());

    if( toRemove == this->DICOMVolumeList.end() )
    {
      std::cerr << "Error: Cannot find dicom volume to be removed from input list." << std::endl;
      return;
    }
    
    this->DICOMVolumeList.erase (toRemove);
  }


  
  
  void GDCMImporter::AddVolume (std::string description, std::vector<std::string> filelist, std::vector<gdcm::File*> *gdcmfilelist)
  {

    itk::DICOMVolume::Pointer dicomvolume = itk::DICOMVolume::New();
    dicomvolume->SetDescription (description);
    dicomvolume->SetFileList (filelist);
    dicomvolume->SetgdcmFileList (gdcmfilelist);
    dicomvolume->SetSerieHelper(this->DICOMScanner->GetSeriesHelper());
    
    try
    {
      this->DICOMVolumeList.push_back(dicomvolume);
    }
    catch (itk::ExceptionObject & e)
    {
      std::cerr << e;
      throw itk::ExceptionObject (__FILE__,__LINE__,"Error in GDCMImporter::AddVolume()");
    }
  }


  

  void GDCMImporter::BuildAllVolumes (void)
  {

    this->UpdateProgress(0.0);
    
    for (unsigned int i=0; i<this->DICOMVolumeList.size(); i++)
    {
      try
      {
      	this->DICOMVolumeList[i]->Build();
	this->UpdateProgress((double)i/(double)this->DICOMVolumeList.size());
      }
      catch (itk::ExceptionObject & e)
      {
	std::cerr << e;
        //throw itk::ExceptionObject (__FILE__,__LINE__,"Error in GDCMImporter::BuildAllVolumes()");
        // We don't throw exception here as we want to reconstruct as many volumes
        // as possible, i.e. even if one volume fails does not mean that all volumes
        // will fail.
      }
    }
  }

  
  DICOMVolume::FloatImageType::Pointer DICOMVolume::TemporaryBuild ()
  {

    if ( this->GetFileList().size() < 2)
    {
      return 0;
    }
    
    SeriesReaderType::Pointer dicomreader = SeriesReaderType::New();
    itk::GDCMImageIO::Pointer io = itk::GDCMImageIO::New();
    
    dicomreader->SetFileNames( this->GetFileList() );
    dicomreader->SetImageIO( io );
    
    try
    {
      dicomreader->Update();

      FloatImageType::Pointer image = dicomreader->GetOutput();

      /**
         The GDCM reader seems to flip in the A-P direction the image.
         We flip it back until this is resolved internally in ITK.
       */
      /*
      FlipFilterType::Pointer flipper = FlipFilterType::New();
      flipper->SetInput (dicomreader->GetOutput());
      FlipFilterType::FlipAxesArrayType axis;
      axis[0] = true;
      axis[1] = true;
      axis[2] = false;
      flipper->SetFlipAxes (axis);

      flipper->Update();

      image = flipper->GetOutput();
      */

      return image;
      
    }
    catch (itk::ExceptionObject & e)
    {
      std::cerr << e;
      throw itk::ExceptionObject (__FILE__,__LINE__,"Error in DICOMVolume::Build()");
    }
  }
  
  
  void DICOMVolume::Build (void)
  {

    /*
    if ( this->GetFileList().size() < 2)
    {
      return;
    }

//     if (this->gdcmFileList->size())
//     {
//       std::cout<<this->Description.c_str()<<", position : "
// 	       <<(*this->gdcmFileList)[0]->GetEntryValue(0x0018,0x5100)<<std::endl
// 	       <<(*this->gdcmFileList)[0]->GetEntryValue(0x0020,0x0020)<<std::endl
// 	       <<(*this->gdcmFileList)[0]->GetEntryValue(0x0020,0x0035)<<std::endl
// 	       <<(*this->gdcmFileList)[0]->GetEntryValue(0x0020,0x0037)<<std::endl
// 	       <<(*this->gdcmFileList)[0]->GetEntryValue(0x0054,0x0410)<<std::endl
//  	       <<(*this->gdcmFileList)[0]->GetEntryValue(0x0018,0x9302)<<std::endl
//  	       <<(*this->gdcmFileList)[0]->GetEntryValue(0x0018,0x1310)<<std::endl
//  	       <<(*this->gdcmFileList)[0]->GetEntryValue(0x0054,0x0410)<<std::endl
//  	       <<(*this->gdcmFileList)[0]->GetEntryValue(0x0018,0x1314)<<std::endl
//  	       <<(*this->gdcmFileList)[0]->GetEntryValue(0x0070,0x0041)<<std::endl<<std::endl;
//     }
    
    SeriesReaderType::Pointer dicomreader = SeriesReaderType::New();
    itk::GDCMImageIO::Pointer io = itk::GDCMImageIO::New();
    
    dicomreader->SetFileNames( this->GetFileList() );
    dicomreader->SetImageIO( io );
    
    try
    {
      dicomreader->Update();

      FloatImageType::Pointer image = dicomreader->GetOutput();

      
      //  The GDCM reader seems to flip in the A-P direction the image.
      //  We flip it back until this is resolved internally in ITK.       
      FlipFilterType::Pointer flipper = FlipFilterType::New();
      flipper->SetInput (dicomreader->GetOutput());
      FlipFilterType::FlipAxesArrayType axis;
      axis[0] = true;
      axis[1] = true;
      axis[2] = false;
      flipper->SetFlipAxes (axis);

      flipper->Update();

      image = flipper->GetOutput();
      
      this->SetImage (image);
    }
    catch (itk::ExceptionObject & e)
    {
      std::cerr << e;
      throw itk::ExceptionObject (__FILE__,__LINE__,"Error in DICOMVolume::Build()");
    }
    */

    try
    {
      this->SetImage ( this->TemporaryBuild() );
    }
    catch (itk::ExceptionObject & e)
    {
      std::cerr << e;
      throw itk::ExceptionObject (__FILE__,__LINE__,"Error in DICOMVolume::Build()");
    }
    
    
  }



  void DICOMVolume::Save (const char* filename)
  {

    if ( this->GetFileList().size() < 2)
    {
      return;
    }

    FloatImageType::Pointer image = this->GetImage();

    if( image.IsNull() )
    {
      
      SeriesReaderType::Pointer dicomreader = SeriesReaderType::New();
      itk::GDCMImageIO::Pointer io = itk::GDCMImageIO::New();
      
      dicomreader->SetFileNames( this->GetFileList() );
      dicomreader->SetImageIO( io );

      /**
         The GDCM reader seems to flip in the A-P direction the image.
         We flip it back until this is resolved internally in ITK.
      */
      /*
      FlipFilterType::Pointer flipper = FlipFilterType::New();
      flipper->SetInput (dicomreader->GetOutput());
      FlipFilterType::FlipAxesArrayType axis;
      axis[0] = true;
      axis[1] = true;
      axis[2] = false;
      flipper->SetFlipAxes (axis);
      
      try
      {
        flipper->Update();
        image = flipper->GetOutput();
      }
      catch (itk::ExceptionObject & e)
      {
        std::cerr << e;
        throw itk::ExceptionObject (__FILE__,__LINE__,"Error in GDCMImporter::Save()");
      }
      */

      try
      {
        dicomreader->Update();
        image = dicomreader->GetOutput();
      }
      catch (itk::ExceptionObject & e)
      {
        std::cerr << e;
        throw itk::ExceptionObject (__FILE__,__LINE__,"Error in GDCMImporter::Save()");
      }
      
    }

    
    itk::Matrix<double,3,3> cosines;
    cosines[0][0]= 1;
    cosines[0][1]= 0;
    cosines[0][2]= 0;
    cosines[1][0]= 0;
    cosines[1][1]=-1;
    cosines[1][2]= 0;
    cosines[2][0]= 0;
    cosines[2][1]= 0;
    cosines[2][2]= 1;
    
    image->SetDirection(cosines);
    
    
    itk::ImageFileWriter<FloatImageType>::Pointer myWriter = itk::ImageFileWriter<FloatImageType>::New();
    myWriter->SetInput( image );
    myWriter->SetFileName( filename );


    try
    {
      myWriter->Update();
    }
    catch (itk::ExceptionObject &e)
    {
      std::cerr << e;
      throw itk::ExceptionObject (__FILE__,__LINE__,"Error while writing the image.");
    }
    
  }
    

  

  void GDCMImporter::Save(int n, const char* filename)
  {
    // save a volume without keeping it in memory
    if( n<0 || n>=this->DICOMVolumeList.size() )
    {
      throw itk::ExceptionObject(__FILE__,__LINE__, "Error in GDCMImporter::Save(): Index is out of range.");
    }


    try
    {
      this->DICOMVolumeList[n]->Save (filename);
    }
    catch (itk::ExceptionObject & e)
    {
      std::cerr << e;
      throw itk::ExceptionObject (__FILE__,__LINE__,"Error in GDCMImporter::Save()");
    }
  
  }
  
  
  void GDCMImporter::SaveAll(const char* directory)
  {
    for( unsigned int i=0; i<this->DICOMVolumeList.size(); i++)
    {
      try
      {
        std::string filename = directory;
        filename += "/" + this->DICOMVolumeList[i]->GetDescription() + ".hdr";
        this->Save(i, filename.c_str());
      }
      catch(itk::ExceptionObject &e)
      {
        std::cerr << e;
      }
      this->UpdateProgress( double(i)/double(this->DICOMVolumeList.size()) );
    }
  }




  gdcm::XCoherentFileSetmap GDCMImporter::SplitOnPosition(gdcm::FileList *fileSet)
  {
    gdcm::XCoherentFileSetmap CoherentFileSet;
    
    size_t nb = fileSet->size();
    if (nb == 0 )
      return CoherentFileSet;
    float pos[3];
    std::string strImPos;  // read on disc
    itksys_ios::ostringstream ossPosition;
    std::string strPosition; // re computed
    gdcm::FileList::const_iterator it = fileSet->begin();
    for ( ; it != fileSet->end(); it++)
    {     
      // Information is in :      
      // 0020,0032 : Image Position Patient
      // 0020,0030 : Image Position (RET)
      
      strImPos = (*it)->GetEntryValue(0x0020,0x0032);
      if ( strImPos == gdcm::GDCM_UNFOUND)
      {
        std::cerr << "Unfound Image Position Patient (0020,0032)\n";
        strImPos = (*it)->GetEntryValue(0x0020,0x0030); // For ACR-NEMA images
        if ( strImPos == gdcm::GDCM_UNFOUND )
        {
          std::cerr << "Unfound Image Position (RET) (0020,0030)\n";
          // User wants to split on the 'Position'
          // No 'Position' info found.
          // We return an empty Htable !
          return CoherentFileSet;
        }  
      }
      
      if ( sscanf( strImPos.c_str(), "%f \\%f \\%f ", 
                   &pos[0], &pos[1], &pos[2]) != 3 )
      {
		  std::cerr << "Wrong number for Position : ["
                       << strImPos << "]\n";
             return CoherentFileSet;
      }

      // Let's build again the 'position' string, to be sure of it's format      

      ossPosition << pos[0];      
      for (int i = 1; i < 3; i++)
      {
        ossPosition << "\\";
        ossPosition << pos[i]; 
      }      
      strPosition = ossPosition.str();
      ossPosition.str("");

      if ( CoherentFileSet.count(strPosition) == 0 )
      {
        //std::cerr << " New Position :[" << strPosition << "]\n";
        // create a File set in 'position' position
        CoherentFileSet[strPosition] = new gdcm::FileList;
      }
      // Current Position and DICOM header match; add the file:
      CoherentFileSet[strPosition]->push_back( (*it) );
    }   
    return CoherentFileSet;
  }
  
  unsigned int* DICOMVolume::GetDimensions(void)
  {
    
    unsigned int* ret = new unsigned int[3];
    ret[0] = 0;
    ret[1] = 0;
    ret[2] = 0;
    
    
    if (this->gdcmFileList->size())
    {
      ret[0] = (*this->gdcmFileList)[0]->GetXSize();
      ret[1] = (*this->gdcmFileList)[0]->GetXSize();
    }
    ret[2] = this->gdcmFileList->size();
    
    return ret;
  }
  
  
}