JEMRIS  2.8.2
open-source MRI simulations
In depth tutorial - Part 2


Let us now consider extending the framework to achieve the same results. One may ask why this could be of any interest. First of all, the symbolic math toolbox requires additional computation power. The implementation of the nessecary code in the simulator framwork speeds up the simualtion process dramatically.

The workflow will involve the creation of two new classes one for the RF pulse and one for the gradient pulses.

Let us start with defining the template for an RF pulse, which inherits from the RF pulse functionality implemented in framework. This is done with a few lines for the object declaration.

Enter the directory in the [jemris-dir]/app/src/framework. Create a new header file for the STA RF pulse by the name STARFPulse.h. And type in the lines below.

#ifndef STARFPULSE_H_
#define STARFPULSE_H_
#include "RFPulse.h"
// Prototype for STA RF pulses
class STARFPulse : public RFPulse {
public:
// Default constructor
STARFPulse() {};
// Default copy constructor.
STARFPulse (const STARFPulse&) {};
// Default destructor.
~STARFPulse () {};
// Cloning mechanism as descibed by Gamma et al "Design patterns"
inline STARFPulse* Clone() const { return (new STARFPulse(*this)); };
// Prepare the pulse with nessecary parameters.
//
// param mode. Sets the preparation mode, one of enum PrepareMode
// {PREP_INIT,PREP_VERBOSE,PREP_UPDATE}.
// returns Success of the preparation.
virtual bool Prepare (PrepareMode mode);
// Returns a constant Magnitude for all times.
//
// param time. The magnitude is requested for "time". Offset is the start time of the pulse.
// returns Magnitude at time "time".
inline virtual double GetMagnitude (double time );
private:
// Constants as described in the paper.
double m_Gamma; // Parameter Gamma.
double m_Alpha; // Parameter Alpha.
double m_A; // Parameter A.
double m_Beta; // Parameter Beta.
int m_N; // Parameter N.
};
#endif /*STARFPULSE_H_*/

The interesting part of the class will be regarding Prepare and GetMagnitude functions. The Prepare does as name suggests. It will prepare the pulse at start up. While the GetMagnitude function implementation will do the calculation of the pulse during runtime.

The rest of the object declaration is similar for any newly added framework element inherting from the virtual implementation RFPulse.

Let us now have a look at the implementation. Add STARFPulse.cpp to the framework directory with the following content.

#include "STARFPulse.h"
/***********************************************************/
bool STARFPulse::Prepare (PrepareMode mode) {
bool btag = true;
// Every attribute, which needs to be set with values before the simulation starts
// needs to be declared to the framework's attribute decalaration and observation
// mechanism. With the macro ATTRIBUTE("ATTRIBUTE_NAME", &m_attribute_name);
ATTRIBUTE("Gamma", &m_Gamma);
ATTRIBUTE("Alpha", &m_Alpha);
ATTRIBUTE("A" , &m_A);
ATTRIBUTE("Beta" , &m_Beta);
ATTRIBUTE("N" , &m_N);
ATTRIBUTE("TPOIs", &m_more_tpois ); //number of TPOIs along the analytical expression
if ( mode == PREP_VERBOSE) {
// XML error checking
// This mechanism supplies the functionality for checking for proper
// decalartion of sequence modules at runtime.
if (!HasDOMattribute("Gamma")) {
btag = false;
cout << GetName() << "::Prepare() error: gamma required!!\n";
}
if (!HasDOMattribute("Alpha")) {
btag = false;
cout << GetName() << "::Prepare() error: Alpha required!!\n";
}
if (!HasDOMattribute("A")) {
btag = false;
cout << GetName() << "::Prepare() error: A required!!\n";
}
if (!HasDOMattribute("Beta")) {
btag = false;
cout << GetName() << "::Prepare() error: Beta required!!\n";
}
if (!HasDOMattribute("N")) {
btag = false;
cout << GetName() << "::Prepare() error: N required!!\n";
}
}
// Extremely vital step. Here the preparation is passed up the chain for
// inheriting framework members. This encapsulates common functionality
// and, thus, follows the object oriented paradigm.
btag = (RFPulse::Prepare(mode) && btag);
// Everything went well?
if (!btag && mode == PREP_VERBOSE)
cout << "\n warning in Prepare(1) of RFPULSE " << GetName() << endl;
// Then we are done and the pulse was successfully prepared.
return btag;
};
/***********************************************************/
double STARFPulse::GetMagnitude (double const time) {
// Strait forward: Return the value for this pulse at "time" in ms from the beginning
// of the pulse.
return m_Gamma * m_Alpha / GetDuration() * exp(-pow(m_Beta,2)*pow(1-time/GetDuration(),2)) *
sqrt( pow(2*PI*m_N*(1-time/GetDuration()),2) + 1);
}

Let us now implement the gradients. Since the two gradients only differ little we will implement both in the same class <>STASpiralGradient<>. Let us start with the object declaration STARGradPulse.h:

#ifndef _SPIRALGRADPULSE_H
#define _SPIRALGRADPULSE_H
#include <cmath>
#include "GradPulse.h"
// Prototype of a STA spiral gradient
class STASpiralGradPulse : public GradPulse {
public:
// Constructor
STASpiralGradPulse () {};
// Copy constructor.
STASpiralGradPulse (const SpiralGradPulse&) {};
// Destructor.
// Cloning mechanism
STASpiralGradPulse* Clone () const;
// Prepare the spiral gradient pulse.
// Please refer to above.
virtual bool Prepare (PrepareMode mode);
// Get gradient amplitude.
// Please refer to GetMagnitude above.
virtual double GetGradient (double const time);
protected:
// Get information for output.
virtual string GetInfo ();
// Constants as described in the paper.
double m_A;
double m_N;
};
#endif /*STASPIRALGRADPULSE_*/

The code of the implementation will now reflect the physics in STAGradPulse.cpp

#include "STASpiralGradPulse.h"
/***********************************************************/
double STASpiralGradPulse::GetGradient (double const time) {
double magn = 0.0;
double angle = 2*PI*m_N*time/GetDuration();
double ancos = cos(angle);
double ansin = sin(angle);
if (GetAxis() == AXIS_GX)
magn = - (m_A/GetDuration()) * ancos - 2*m_A*PI*m_N/GetDuration()*(1-time/GetDuration()) * ansin;
else if (GetAxis() == AXIS_GY)
magn = - (m_A/GetDuration()) * ansin + 2*m_A*PI*m_N/GetDuration()*(1-time/GetDuration()) * ancos;
return magn;
}
/***********************************************************/
bool STASpiralGradPulse::Prepare (PrepareMode mode) {
bool btag = true;
ATTRIBUTE("A" , &m_A);
ATTRIBUTE("N" , &m_N);
if ( mode == PREP_VERBOSE) {
if (!HasDOMattribute("Turns")) {
btag = false;
cout << GetName() << "::Prepare() error: Alpha required!!\n";
}
if (!HasDOMattribute("Pitch")) {
btag = false;
cout << GetName() << "::Prepare() error: Pitch required!!\n";
}
SetArea(0.0); // Inherently refocused pulse; no refocusing necessary!
}
btag = (GradPulse::Prepare(mode) && btag);
if (!btag && mode == PREP_VERBOSE)
cout << "\n warning in Prepare(1) of TRAPGRADPULSE " << GetName() << endl;
return btag;
}
/***********************************************************/
string STASpiralGradPulse::GetInfo() {
// Debug output
stringstream s;
s << GradPulse::GetInfo() << " , (A,N)= (" << m_A << "," << m_N << ")";
return s.str();
};

We now have to introduce the two new modules to the framework factory by adding them to the ModulePrototypeFactory. Just append the two include at the end of the module list.

...
#include "STASpiralGradPulse.h"
#include "STARFPulse.h"

Open ModulePrototypeFactory.cpp and add the headers of the new modules. And, further, add the object prototypes to the factory by appending the following two lines after the last entry in ModulePrototypeFactory::ModulePrototypeFactory:

...
m_Clonables.insert( pair<string,Module*>( "STASPIRALGRADPULSE", new STASpiralGradPulse () ));
m_Clonables.insert( pair<string,Module*>( "STARFPULSE", new STARFPulse () ));
}

We are done with the framework extension.


-- last change 24.05.2018 | Tony Stoecker | Imprint | Data Protection --