#include "image.h"
#include "data.h"
#include <odinpara/sample.h>



Image::Image(const STD_string& label) : LDRblock(label) {
  magnitude.set_label("magnitude");
  magnitude.set_filemode(compressed);

  append_all_members();
}


Image& Image::operator = (const Image& i) {
  LDRblock::operator = (i);
  geo=i.geo;
  magnitude=i.magnitude;
  append_all_members();
  return *this;
}


unsigned int Image::size(axis ax) const {
  ndim nn=get_magnitude().get_extent();
  int index=nn.size()-1-ax;
  if(index>=0) return nn[index];
  return 1;
}



Image& Image::transpose_inplane(bool reverse_read, bool reverse_phase) {

  unsigned int d=magnitude.dim();
  if(d<2) return *this;

  geo.transpose_inplane(reverse_read,reverse_phase);

  // transpose image data
  farray magnitude_copy(magnitude);
  ndim shape=magnitude.get_extent();
  STD_swap(shape[d-1],shape[d-2]);
  magnitude.redim(shape);
  for(unsigned int i=0; i<magnitude.length(); i++) {
    ndim index(magnitude_copy.create_index(i));
    if(reverse_read) index[d-1]=shape[d-1]-1-index[d-1];
    if(reverse_phase) index[d-2]=shape[d-2]-1-index[d-2];
    STD_swap(index[d-1],index[d-2]);
    magnitude(index)=magnitude_copy[i];
  }

  return *this;
}




void  Image::append_all_members() {
  LDRblock::clear();
  LDRblock::merge(geo);
  append_member(magnitude);
}


/////////////////////////////////////////////////////////////////

ImageSet::ImageSet(const STD_string& label) : LDRblock(label) {
  Content.set_label("Content");

  append_all_members();
}


ImageSet::ImageSet(const Sample& smp) {
  Log<OdinData> odinlog(this,"ImageSet(Sample)");

  float min_fov=100.0;
  int min_size=64;


  farray sddata(smp.get_spinDensity());
  STD_string maplabel="Spin Density";
  sddata.normalize();
  ODINLOG(odinlog,normalDebug) << "sddata.get_extent()=" << sddata.get_extent() << STD_endl;


  float xFOV=smp.get_FOV(xAxis);
  float yFOV=smp.get_FOV(yAxis);
  float zFOV=smp.get_FOV(zAxis);
  float maxFOV=maxof3(xFOV,yFOV,zFOV);
  if(maxFOV<min_fov) maxFOV=min_fov;
  ODINLOG(odinlog,normalDebug) << "maxFOV/xFOV/yFOV/zFOV=" << maxFOV << "/" << xFOV << "/" << yFOV << "/" << zFOV << STD_endl;

  int nx=sddata.size(xDim);
  int ny=sddata.size(yDim);
  int nz=sddata.size(zDim);
  int nfreq=sddata.size(freqDim);
  int maxsize=int(maxof3(nx,ny,nz));
  if(maxsize<min_size) maxsize=min_size;
  ODINLOG(odinlog,normalDebug) << "maxsize/nx/ny/nz/nfreq=" << maxsize << "/" << nx << "/" << ny << "/" << nz << "/" << nfreq << STD_endl;

  STD_string dirlabel;

  // Sagittal & Coronal

  Geometry iasagcor;
  for(int ichan=0; ichan<n_directions; ichan++)  iasagcor.set_FOV(direction(ichan), maxFOV);

  farray data(maxsize,maxsize);

  dvector rpsvec(3);
  dvector xyzvec(3);

  for(int idir=sagittal; idir<=coronal; idir++) {

    data=0.0;

    iasagcor.set_orientation(sliceOrientation(idir));

    ODINLOG(odinlog,normalDebug) << "gradrotmatrix[" << idir << "]=" << iasagcor.get_gradrotmatrix().print() << STD_endl;

    for(int j=0; j<maxsize; j++) {
      for(int i=0; i<maxsize; i++) {

        rpsvec[0]=((i+0.5)/float(maxsize)-0.5)*maxFOV; // read
        rpsvec[1]=((j+0.5)/float(maxsize)-0.5)*maxFOV; // phase
        rpsvec[2]=0.0;
        xyzvec=iasagcor.transform(rpsvec);

        int ix=int(nx*(xyzvec[0]/xFOV+0.5)-0.5); // x
        int iy=int(ny*(xyzvec[1]/yFOV+0.5)-0.5); // y
        int iz=int(nz*(xyzvec[2]/zFOV+0.5)-0.5); // z

        if(ix>=0 && ix<nx && iy>=0 && iy<ny && iz>=0 && iz<nz) data(j,i)=sddata(0,0,iz,iy,ix);
      }
    }

    if(idir==sagittal) dirlabel="Sagittal";
    if(idir==coronal) dirlabel="Coronal";


    Image imgsagcor(maplabel+"("+dirlabel+")");
    iasagcor.set_label(dirlabel);
    imgsagcor.set_geometry(iasagcor);
    imgsagcor.set_magnitude(data);
    imgsagcor.transpose_inplane(); // for correct display in geoedit
    append_image(imgsagcor);
  }



  // Axial as multi slice
  dirlabel="Axial";
  Geometry ia(dirlabel);
  ia.set_FOV(readDirection, smp.get_FOV(xAxis));
  ia.set_FOV(phaseDirection,smp.get_FOV(yAxis));
  ODINLOG(odinlog,normalDebug) << "gradrotmatrix[axial]=" << ia.get_gradrotmatrix().print() << STD_endl;

  ia.set_nSlices(nz);
  float fovz=smp.get_FOV(zAxis);
  float dz=secureDivision(fovz,nz);
  ia.set_sliceThickness(dz);
  ia.set_sliceDistance(dz);

  Image img(maplabel+"("+dirlabel+")");
  img.set_geometry(ia);
  sddata.autosize();
  img.set_magnitude(sddata);
  append_image(img);
}




ImageSet& ImageSet::operator = (const ImageSet& is) {
  LDRblock::operator = (is);
  Content=is.Content;
  images=is.images;
  append_all_members();
  return *this;
}

ImageSet& ImageSet::append_image(const Image& img) {
  Log<OdinData> odinlog(this,"append_image");
  bool rename=false;
  if(img.get_label()=="" || LDRblock::parameter_exists(img.get_label())) rename=true;
  images.push_back(img);
  STD_list<Image>::iterator it=images.end();
  --it;
  if(rename) it->set_label("Image"+itos(images.size()-1));
  LDRblock::append(*it);

  unsigned int nimages=images.size();
  ODINLOG(odinlog,normalDebug) << "nimages=" << nimages << STD_endl;

  Content.resize(nimages);
  int index=0;
  for(it=images.begin(); it!=images.end(); ++it) {
    ODINLOG(odinlog,normalDebug) << "index/label=" << index << "/" << it->get_label() << STD_endl;
    Content[index]=it->get_label();
    index++;
  }
  return *this;
}


ImageSet& ImageSet::clear_images() {
  images.clear();
  Content.resize(0);
  return *this;
}


Image& ImageSet::get_image(unsigned int index) {
  Log<OdinData> odinlog(this,"get_image");
  unsigned int nimages=images.size();
  ODINLOG(odinlog,normalDebug) << "index/nimages=" << index << "/" << nimages << STD_endl;
  if(index>=nimages) return dummy;
  STD_list<Image>::iterator it=images.begin();
  for(unsigned int i=0; i<index; i++) ++it;
  return *it;
}


int ImageSet::load(const STD_string& filename, const LDRserBase& serializer) {
  Log<OdinData> odinlog(this,"load");

  clear_images();

  // check whether we have a whole set of images and determine the total num of images first
  int noiresult=Content.load(filename,serializer);
  svector contcopy=Content; // create copy because append_image will modify Content
  unsigned int nimages=contcopy.size();
  ODINLOG(odinlog,normalDebug) << "noiresult/nimages/Content=" << noiresult << "/" << nimages << "/" << contcopy.printbody() << STD_endl;

  int result=0;
  if(noiresult>0) {
    // create placeholders
    Image img;
    for(unsigned int i=0; i<nimages; i++) {
      img.set_label(contcopy[i]);
      append_image(img);
      ODINLOG(odinlog,normalDebug) << "Appended empty image for " << STD_string(contcopy[i]) << STD_endl;
    }
    ODINLOG(odinlog,normalDebug) << "Loading all images" << STD_endl;
    result=LDRblock::load(filename,serializer);
  } else {
    // try single image
    Image img;
    result=img.load(filename,serializer);
    if(result>0) {
      clear_images();
      append_image(img);
    }
  }
  return result;
}


void ImageSet::append_all_members() {
  LDRblock::clear();
  append_member(Content);
}