odin 1.8.8-2ubuntu1 / usr / share / odin / sequences / odinpilot.cpp

//////////////////////////////////////////////////////////////////////////////
//
// A navigator sequence for the ODIN framework
//
//////////////////////////////////////////////////////////////////////////////



// includes for the sequence objects
#include <odinseq/seqall.h>



// The method itself is a class that is derived from
// an abstract base class 'SeqMethod' that contains all
// routines common to all methods:

class METHOD_CLASS : public SeqMethod {

 public:

  // Constructor that takes the methods identifier (unique string) as its argument
  METHOD_CLASS(const STD_string& label);

  // virtual functions that are overwritten in this method to build the sequence,
  // calculate parameter relations, etc.:
  void method_pars_init();
  void method_seq_init();
  void method_rels();
  void method_pars_set();

 private:

  // Parameters for this method:
  JDXfloat T1Ernst;
  JDXint   DummyCycles;

  JDXint NumSagSlices;
  JDXint NumCorSlices;
  JDXint NumAxiSlices;
  JDXfloat PilotSliceDist;

  // Sequence objects for this method:
  SeqPulsar exc;

  SeqDelay relaxdelay;

  SeqObjLoop peloop;
  SeqObjLoop sliceloop;
  SeqObjLoop reploop;
  SeqObjLoop dummyloop;

  SeqObjList scanpart;
  SeqObjList dummypart;

  SeqDelay dummypad;

  SeqGradEcho grech;

  SeqGradTrapezParallel crusher;
  SeqDelay   crusherdelay;

  SeqRotMatrixVector pilotrot;


  unsigned int nslices;
  ivector nslices_dir;
  dvector dirvec;
  Geometry area[n_orientations];

};

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

METHOD_CLASS::METHOD_CLASS (const STD_string& label)
                                   : SeqMethod(label) {
  // Put in here all stuff that will once be initialised and
  // never changed

  // Specify a short description of the method
  set_description("Navigator which acquires three orthogonal sets of slices." );
}


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

void METHOD_CLASS::method_pars_init() {
  // In this function the methods parameters will be initialised


  // Assign default values:
  commonPars->set_RepetitionTime(300.0);
  commonPars->set_AcqSweepWidth(25.0);
  commonPars->set_MatrixSize(readDirection,128);


  // set default values and a short a description (which
  // will appear as a tooltip in the UI) for each parameter
  T1Ernst=1300.0;
  T1Ernst.set_minmaxval(0.0,5000.0).set_description("For optimum SNR, the flip angle will be set to the Ernst angle using this T1");

  DummyCycles=3;
  DummyCycles.set_description("Number of dummy shots before actual acquisition");


  NumSagSlices=3;
  NumSagSlices.set_description("Number of sagittal slices in Pilot");

  NumCorSlices=3;
  NumCorSlices.set_description("Number of coronal slices in Pilot");

  NumAxiSlices=3;
  NumAxiSlices.set_description("Number of axial slices in Pilot");

  PilotSliceDist=25.0;
  PilotSliceDist.set_description("Slice distance of Pilot");


  // Make method specific parameters visible in the user interface
  append_parameter(T1Ernst,"T1Ernst");
  append_parameter(DummyCycles,"DummyCycles");
  append_parameter(NumSagSlices,"NumSagSlices");
  append_parameter(NumCorSlices,"NumCorSlices");
  append_parameter(NumAxiSlices,"NumAxiSlices");
  append_parameter(PilotSliceDist,"PilotSliceDist");

}


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


void METHOD_CLASS::method_seq_init() {
  Log<Seq> odinlog(this,"method_seq_init");

  // Put in here all stuff to create the layout of the sequence


  // number of slices in each direction
  nslices_dir.resize(n_orientations);
  nslices_dir[sagittal]=NumSagSlices;
  nslices_dir[coronal]=NumCorSlices;
  nslices_dir[axial]=NumAxiSlices;


  nslices=nslices_dir.sum();

  geometryInfo->reset(); // Use default orientations



  ///////////////// Excitation Pulse: /////////////////////

  // Get the slice thickness from the global geometry handler 'geometryInfo':

  float slicethick=4.0;

  // Create the sinc shaped excitation pulse.
  exc=SeqPulsarSinc("exc",slicethick,false,2.0,commonPars->get_FlipAngle());

  // This is useful for simulating/visualization of the sequence
  exc.set_pulse_type(excitation);

  // Do not use maxDistEncoding
  exc.set_encoding_scheme(linearEncoding);

  // Pilot rotation and slice offset

  dvector offset(nslices);
  offset=0.0;

  pilotrot=SeqRotMatrixVector("pilotrot");

  dirvec.resize(nslices);

  int ioffset=0;
  for(int idir=0; idir<n_orientations; idir++) {

    area[idir].set_orientation(sliceOrientation(idir));

    area[idir].set_nSlices(nslices_dir[idir]);
    area[idir].set_sliceDistance(PilotSliceDist);
    dvector slicevec=area[idir].get_sliceOffsetVector();

    for(int islice=0; islice<nslices_dir[idir]; islice++) {
      offset[ioffset]=slicevec[islice];
      dirvec[ioffset]=idir;
      pilotrot.append(area[idir].get_gradrotmatrix());
      ioffset++;
    }
  }

  ODINLOG(odinlog,significantDebug) << "offset=" << offset << STD_endl;

  // Set the frequency list of the excitation pulse so that we will excite all slices in the slicepack
  exc.set_freqlist( systemInfo->get_gamma("") * exc.get_strength() / (2.0*PII) *  offset );


  ivector iv(nslices);
  iv=0;
  exc.set_indexvec(iv); // Always set zero slice index, indexing will be done by userdef dimension


  // This loop object is used to loop over the slices
  sliceloop=SeqObjLoop("sliceloop");



  ////////////////// Geometry: /////////////////////////////////


  // calculate the resolution in the read Channel and set the number of phase encoding
  // steps so that we will obtain a uniform resolution in read and phase Channel:
  float resolution=secureDivision(geometryInfo->get_FOV(readDirection),commonPars->get_MatrixSize(readDirection));
  commonPars->set_MatrixSize(phaseDirection,int(secureDivision(geometryInfo->get_FOV(phaseDirection),resolution)+0.5),noedit);


  //////////////// Delays: //////////////////////////////

  // relaxation delay after each readout
  relaxdelay=SeqDelay("relaxdelay");

  //////////////// Gradient Echo Module: //////////////////////////////

  float os_factor=2.0; // avoid cut-off in read direction for off-center FOV

  // This is the gradient-recalled echo kernel of the sequence.
  // Calculations of read/phase gradient strengts and durations are
  // performed by this module itself so we just have to pass in the size and FOV
  // in each direction. This module has a tight timing scheme so that short TE is possible.
  grech=SeqGradEcho("grech", exc, commonPars->get_AcqSweepWidth(),
              commonPars->get_MatrixSize(readDirection), geometryInfo->get_FOV(readDirection),
              commonPars->get_MatrixSize(phaseDirection), geometryInfo->get_FOV(phaseDirection),
              linearEncoding, noReorder, 1, commonPars->get_ReductionFactor(), DEFAULT_ACL_BANDS,
              false, commonPars->get_PartialFourier(),0.0,false,os_factor);

  grech.set_gradrotmatrixvector(pilotrot); // rotate grech


  //////////////// Crusher Gradient: //////////////////////////////

  float crusher_strength=0.5*systemInfo->get_max_grad();
  float crusher_integral=4.0*fabs(grech.get_gradintegral().sum());

  crusher=SeqGradTrapezParallel("crusher",crusher_integral,crusher_integral,crusher_integral, crusher_strength);
  crusherdelay=SeqDelay("crusherdelay",0.1); // Small delay to avoid gradient-induced stimulation


  //////////////// Loops: //////////////////////////////

  // Construct a loop object to iterate through the phase encoding steps
  peloop=SeqObjLoop("peloop");

  // Construct a loop object to perform repetitions of the experiment
  reploop=SeqObjLoop("reploop");

  // Construct a loop object to perform dummy cylces
  dummyloop=SeqObjLoop("dummyloop");


  //////////////// Constructing the Sequence: ///////////////


  // This sequence container will hold the objects for one readout:
  scanpart=SeqObjList("scanpart");
  scanpart = grech + crusherdelay + crusher;


  // Contains kernel of one dummy cycle
  dummypart=SeqObjList("dummypart");
  dummypad=SeqDelay("dummypad",scanpart.get_duration()-exc.get_duration()); // Padding delay for dummy cycle
  dummypart = exc + dummypad + relaxdelay;



  // Finally, build the whole sequence
  set_sequence(
    dummyloop(
      sliceloop(
        dummypart
      )[grech.get_exc_vector()][pilotrot]
    )[DummyCycles]
    + peloop(
        sliceloop(
          scanpart + relaxdelay   // the sequence kernel
        )[grech.get_exc_vector()][pilotrot]
      )[grech.get_pe_vector()]
  );

}


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


void METHOD_CLASS::method_rels() {
  // Put in here all stuff that has to be performed whenever one of the sequence parameters
  // has been changed by the user


  // calculate relaxdelay to get the desired repetition time
  float scandur=scanpart.get_duration()*float(nslices);
  if(scandur>commonPars->get_RepetitionTime()) commonPars->set_RepetitionTime(scandur);
  relaxdelay.set_duration( (commonPars->get_RepetitionTime()-scandur)/float(nslices) );


  // calculate Ernst angle accordng to TR
  float flipangle=180.0/PII * acos( exp ( -secureDivision ( commonPars->get_RepetitionTime(), T1Ernst) ) );
  commonPars->set_FlipAngle( flipangle, noedit );
  exc.set_flipangle( flipangle );


  // retrieve echo time from gradient echo module to display it in the user interface
  commonPars->set_EchoTime(grech.get_echo_time(), noedit );
}


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

void METHOD_CLASS::method_pars_set() {
  // Put in here all stuff that has to be performed after the parameters have been edited by the user
  // and before the sequence is played out

  // for odinreco
  grech.set_reco_vector(userdef,pilotrot,dirvec); // distinguish slices by their orientation
  recoInfo->set_Recipe("kspace | fft | offset | oversampling | usercoll | pilot("+ftos(PilotSliceDist)+") | store");

}


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

// entry point for the sequence module
ODINMETHOD_ENTRY_POINT