filter_range.h

/***************************************************************************
                          filter_range.h  -  description
                             -------------------
    begin                : Tue Jan 29 2008
    copyright            : (C) 2001 by Thies Jochimsen
    email                : jochimse@cns.mpg.de
 ***************************************************************************/

/***************************************************************************
 *                                                                         *
 *   This program is free software; you can redistribute it and/or modify  *
 *   it under the terms of the GNU General Public License as published by  *
 *   the Free Software Foundation; either version 2 of the License, or     *
 *   (at your option) any later version.                                   *
 *                                                                         *
 ***************************************************************************/

#ifndef FILTER_RANGE_H
#define FILTER_RANGE_H

#include <odindata/filter_step.h>


bool str2range(const STD_string& str, Range& range, int srcsize);


template<int Dir>
class FilterRange : public FilterStep {

  JDXstring range;

  STD_string label() const {return STD_string(1,STD_string(dataDimLabel[Dir])[0])+"range";} // Use only 1st letter
  STD_string description() const {return "Select range in "+STD_string(dataDimLabel[Dir])+" direction";}

  bool process(Data<float,4>& data, Protocol& prot) const {
    Range all=Range::all();

    Range rng[4];
    for(int i=0; i<4; i++) rng[i]=Range::all();

    if(!str2range(range,rng[Dir],data.extent(Dir))) return false;

    TinyVector<int,4> newshape=data.shape();
    newshape(Dir)=rng[Dir].length();

    float fovfraction=secureDivision(rng[Dir].last()-rng[Dir].first()+1, data.shape()[Dir]);

    float offsetfactor=secureDivision(0.5*(rng[Dir].first()+rng[Dir].last()), data.shape()[Dir])-0.5;

    Data<float,4> data_copy(data.copy());
    data.resize(newshape);
    data(all,all,all,all)=data_copy(rng[timeDim],rng[sliceDim],rng[phaseDim],rng[readDim]);

    // adapt protocol
    int stride=rng[Dir].stride();
    if(Dir==timeDim) {
      prot.seqpars.set_NumOfRepetitions(newshape(timeDim));
      if(stride>1) prot.seqpars.set_RepetitionTime(prot.seqpars.get_RepetitionTime()*stride);

    } else {

      direction geodir=direction(3-Dir);
      prot.geometry.set_offset(geodir,prot.geometry.get_offset(geodir)+offsetfactor*prot.geometry.get_FOV(geodir));
      prot.geometry.set_FOV(geodir,fovfraction*prot.geometry.get_FOV(geodir));
      prot.seqpars.set_MatrixSize(geodir,newshape(Dir));

      if(Dir==sliceDim) {
        if(prot.geometry.get_Mode()==slicepack) {
          prot.geometry.set_nSlices(newshape(sliceDim));
          if(stride>1) prot.geometry.set_sliceDistance(prot.geometry.get_sliceDistance()*stride);
          prot.seqpars.set_MatrixSize(sliceDirection,1);
        }
      }
    }

    return true;
  }

  FilterStep*  allocate() const {return new FilterRange<Dir>();}

  void init() {
    range.set_description("Single value or range, optionally with increment (e.g. 1-10:3)");
    append_arg(range,"range");
  }
};

#endif