This user guide explains how to use spidyan for advanced pulse EPR simulations and to investigate spin dynamics.
It is easier to follow this guide, if you are familiar with part I and part II of the guide on saffron.
For the pulse definitions the documentation page of pulse and the corresponding guide are a great starting point.
This guide covers the following topics:
A spidyan script has the same layout as a saffron script. First a spin system is defined:
Sys.ZeemanFreq = 9.5; % in GHzHere, we already encounter the first difference between spidyan and saffron. saffron simulates powder averaged (or single crystal) EPR experiments, spidyan only works with one single spin system at a time. While spidyan still accepts
Sys.g
, it can be more convenient to define the Zeeman frequency Sys.ZeemanFreq
of the investigated spin packet.
After the spin system follows the experiment. Here we are creating a two-pulse sequence:
p90.Flip = pi/2; % flip angle in rad p90.tp = 0.02; % pulse length in μs p180.Flip = pi; % flip angle in rad p180.tp = 0.04; % pulse length in μs tau = 0.5; % delay in μs Exp.Sequence = {p90 tau p180 tau}; % pulse sequence Exp.mwFreq = 9.5; % carrierfrequency
Here we create two monochromatic rectangular pulses with flip angles π/2 and π that are separated by a free evolution period.
Running the above code with
spidyan(Sys,Exp)
will plot the evolution of the electron coherence during the entire sequence (this includes during the pulses).
You might be expecting an echo to form at some point, after all this is the common two-pulse echo sequence.
But keep in mind, an echo is a property of a distribution of spins with different resonance frequencies - the echo forms at a certain point of time, when all the spins are refocused.
Since we called spidyan with a single resonance frequency, no echo is to be expected.
This section is going to show three different inversion pulses and how to plot their effects on the spin. First, we start again with the spin system.
Sys.ZeemanFreq = 9.5; % in GHzNext, lets define some general experimental parameters, like the detection operator and the carrier frequency:
Exp.DetOperator = 'z1'; Exp.mwFreq = 9.5; % GHz
Since we are interested in the inversion of the spin we chose Ŝz as detection operator.
First, let us start with a commong monochromatic rectangular pulse:
Monochromatic.tp = 0.06; % pulse length in μs Monochromatic.Flip = pi; % flip angle in rad
This is the minimal required input for monochromatic pulse. If no Pulse.Type
is given, a rectangular pulse is assumed by default.
The same is true for the parameter Pulse.Frequency
, which is set to 0
by default, corresponding to a pulse with the frequency of the carrier frequency Exp.mwFreq
.
For completeness, here is an equivalent pulse definition:
Monochromatic.Type = 'rectangular'; Monochromatic.tp = 0.06; Monochromatic.Flip = pi; Monochromatic.Frequency = 0;
Next, let us build Exp.Sequence
, run the simulation and then plot it
Exp.Sequence = {Monochromatic}; [TimeAxis, Signal] = spidyan(Sys,Exp); % plotting: figure(1) clf hold on xlabel('t (ns)') ylabel('< S_z>') axis tight ylim([-1 1]) plot(TimeAxis*1000,real(Signal));The plot will show you the trajectory of the spin during the pulse, ending up at
<Ŝz> = 1
, as you would expect from a π pulse.
Next, let us look at what happens to a spin during adiabatic passage. We define the pulse and overwrite current Exp.Sequence
LinearChirp.Type = 'quartersin/linear'; LinearChirp.trise = 0.030; % rise time μs LinearChirp.tp = 0.2; % pulse length μs LinearChirp.Frequency = [-50 50]; % frequency band in MHz LinearChirp.Flip = pi; % flip angle in rad % Update sequence Exp.Sequence = {LinearChirp};
According to pulse, this gives a chirp pulse with a length of 200 ns and a bandwidth of of 100 MHz, centered around the carrier frequency Exp.mwFreq
.
To ensure a well-defined excitation profile in the frequency domain, the edges in the time domain are smoothed with a quarter sine, which has a length of 30 ns.
Run and plot this with
[TimeAxis, Signal] = spidyan(Sys,Exp); % Plotting plot(TimeAxis*1000,real(Signal));
and you will see how the spin follows the effective field during adiabatic inversion.
Once again, we define a new pulse and plot the simulation:
HS.Type = 'sech/tanh'; HS.beta = 10; % truncation parameter of 'sech' HS.n = 1; % exponent of the secant function argument HS.tp = 0.2; % pulse length HS.Frequency = [-50 50]; % excitation band, GHz HS.Flip = pi; % flip angle in rad % Update sequence Exp.Sequence = {HS}; [TimeAxis, Signal] = spidyan(Sys,Exp); % Plotting plot(TimeAxis*1000,real(Signal));
Although trajectory is slightly different, both adiabatic pulses end up in the same spin state.
The quality, of adiabatic pulses is in general described by a parameter called critical adiabaticity Qcrit
, which measures the ability of the pulse to invert the spin and depends on the amplitude and the sweep rate.
The flip angle and Qcrit
are related through the Landau-Zener-Stückelberg-Majorana equation:, e.g, a Qcrit
of 5 corresponds to a π pulse and a Qcrit
of about 0.7 to a π/2 pulse.
Hence, it is not only possible to define the pulses with Pulse.Flip
, but instead you can provide Pulse.Amplitude
or Pulse.Qcrit
for frequency-swept pulses.
Adiabatic pulses have the big advantage, that they invert all spins within their excitation band equally.
You can see this, by changing the resonance frequency of the spin Sys.ZeemanFreq
from 9.500
by to 25 MHz to Sys.ZeemanFreq = 9.525
and then comparing the three pulses.
A second advantage of adiabatic pulses is, that overturning of the spin is not possible, something that is common with rectangular pulses, see this example.
This means, no matter how strong the amplitude of your pulse gets, the flip angle can not became larger than π.
You can try this by removing Par.Flip
from the pulse definitions, and instead using Par.Amplitude
.
The aim of this example is to create a three-dimensional depiction of the spin trajectory during adiabatic passage. We start with the spin system:
Sys.ZeemanFreq = 9.500; % resonance frequency of the spin in GHz
Next, the pulse and experiment definition:
Pulse.Type = 'quartersin/linear'; Pulse.trise = 0.015; % rise time in μs Pulse.Qcrit = 5; % critical adiabaticity Pulse.tp = 0.2; % pulse length in μs Pulse.Frequency = [-100 100]; % frequency band in MHz Exp.Sequence = {Pulse}; % μs Exp.mwFreq = 9.500; % GHz
For the three-dimensional plot we need the expectation values of Ŝx, Ŝy and Ŝz.
Since Ŝ+ = Ŝx + i Ŝy
we can use this to our advantag and use for the detection:
Exp.DetOperator = {'z1','+1'}; % detection operators Exp.DetFreq = [0 9.5]; % downconversion frequency in GHz
Exp.DetFreq
contains the frequencies for downconversion in GHz.
Since Ŝz does not contain any rotating components, the down conversion has to be switched off for it, by setting it to 0
.
We can run the simulation and plot it: [TimeAxis, Signal] = spidyan(Sys,Exp);
% Plotting figure(1) clf plot(TimeAxis*1000,real(Signal)); xlabel('t (ns)') axis tight ylim([-1 1]) ylabel('< S_i>') legend(Exp.DetOperator)
In order to get a three-dimensional plot, we need use the real and imaginary parts of <Ŝy>
:
figure(2) clf plot3(real(Signal(:,2)),imag(Signal(:,2)),real(Signal(:,1))); xlabel('< S_x >') ylabel('< S_y>') zlabel('< S_z>') grid on
More on the theory of adiabatic passage can be found in
Gunnar Jeschke, Stephan Pribitzer, Andrin Doll
Coherence Transfer by Passage Pulses in Electron Paramagnetic Resonance Spectroscopy
The Journal of Physical Chemistry B 2015, DOI: 10.1021/acs.jpcb.5b02964
and
Philipp Spindler, et al.
Shaped Pulses in EPR
eMagRes 2016, DOI: 10.1002/9780470034590.emrstm1520
The use of the detection window Exp.DetWindow
was already discussed in part II of the guide on saffron.
With spidyan you get more control over what is detected and when.
The first topic concern detection operators.
By default, if Exp.DetOperator
is omitted, Ŝ+ of all electrons is detected.
But using Exp.DetOperator
it is not only possible to use several operators at once, you can even specify spin and transition selective operators, using the extended sop syntax.
The detection operators are specified in a cell array, e. g.
Exp.DetOperator = {'z1' 'x1' '+2'};
uses Ŝz and Ŝx of the first spin of the spin system and Ŝ+ of the second spin.
For each of the components of Exp.DetOperator
a frequency for downconversion in GHz must be provided.
Detection operators that correspond only to elements on the diagonal of the detection operator do not need to be down converted.
Such operators are any Ŝz-operators, written like 'z1'
, polarization operators between two levels, specified for example with z(1|2)1
and population operators selective operators, specified with like e(2)
.
Assuming that both spins are of the same type, we could write for the above example:
Exp.DetFreq = [0 9.5 9.5]; % down conversion frequency in GHz
This corresponds to no down conversion of 'z1'
and down conversion with 9.5 GHz for 'x1'
and 9.5 GHz for 'x2'
.
For counter-rotating operators, e.g. Ŝ-, the sign of the corresponding element in DetFreq
has to reflect this:
Exp.DetOperator = {'+1' '-1'} Exp.DetFreq = [9.5 -9.5]; % down conversion frequency in GHz
Instead of using Exp.DetWindow
, which only allows detection after the last pulse, you can also use Exp.DetSequence
, which allows you to detect during selected events.
If neither Exp.DetWindow
nor Exp.DetSequence
are defined, spidyan detects the entire sequence by default.
Detection is controlled through booleans for each element in Exp.Sequence
separately. Here are a few examples:
Exp.Sequence = {p90 tau p180 2*tau}; % the pulse experiment Exp.DetSequence = 1 % detect the entire sequence Exp.DetSequence = 0 % no detection at all Exp.DetSequence = [0 0 0 1] % only detect during the last event in Exp.Sequence Exp.DetSequence = [0 1 0 1] % detection during the second and fourth element in Exp.Sequence Exp.DetSequence = [0 0 1 0] % detection during the third event
Keep in mind, that detection not only requires a significant amount of computing time, it also keeps spidyan from using some of it speed-up tricks.
Hence, if you run big simulation with a lot indirect dimensions and data points, you might want to think about reducing detection to a minimum for the sake of speed.
This example is going over the enhancement of the central transition of a S = 3/2
system for a single spin packet.
For this we are going to use two frequency swept pulses that will transfer magnetization from the outer transitions to the central transition.
We can follow the polarization enhancement by using transition selective operators.
First, let us define a spin system with some zero-field splitting:
Sys.S = 3/2; Sys.ZeemanFreq = 33.500; % GHz Sys.D = 166; % MHz
Next we set up two chirp pulses that are being swept from the outside of the spectrum towards the center and put them right after each other in Exp.Sequence
:
Pulse1.Type = 'quartersin/linear'; % make it a chirp Pulse1.trise = 0.05; % smooth edges in time domain with a quarter sine, in mus Pulse1.Qcrit = 10; % use critical adiabaticity instead of Par.Flip or Par.Amplitude Pulse1.tp = 1; % pulse length in mus Pulse1.Frequency = [500 170]; % frequency band, relative to Exp.mwFreq Pulse2.Type = 'quartersin/linear'; % make it a chirp Pulse2.trise = 0.05; % smooth edges in time domain with a quarter sine, in mus Pulse2.Qcrit = 10; % use critical adiabaticity instead of Par.Flip or Par.Amplitude Pulse2.tp = 1; % pulse length in mus Pulse2.Frequency = [-500 -170]; % frequency band, relative to Exp.mwFreq % Sequence Exp.Sequence = {Pulse1 Pulse2};
We then set the carrier frequency and specify our detection operators:
Exp.mwFreq = 33.5; % GHz % The detection operators detect polarization between (1) levels 1 and 2, % levels 2 and (3) levels 3 and 4 Exp.DetOperator = {'z(1|2)' 'z(2|3)' 'z(3|4)'};
We make use of the transition-selective operators from sop:
With 'z(1|2)'
polarization between |3/2>
and |1/2>
is detected, with 'z(2|3)'
the central transition (|1/2>
and |-1/2>
) and with 'z(3|4)'
the polarization of the states |-1/2>
and |-3/2>
.
Running and plotting it with
[TimeAxis, Signal] = spidyan(Sys,Exp); % plotting figure(1) clf plot(TimeAxis,-real(Signal)); xlabel('t ({\mu}s)') axis tight ylim([-1 3]) ylabel('< S_i>') legend(Exp.DetOperator)
shows how the first pulse shuffles polarization from one of the outer transitions to the central transition, and increasing the polarization on the central transition from 1 to to 2.
The second pulse does the same with the other out transition, raising the polarization of the central transition to 3!
Try this with a S = 5/2
or S = 7/2
and see if you can achieve even larger enhancements - you might have to adapt the excitation bands of your pulses and add the correct detection operators!
Now, this simulation was done for a single spin packet (e.g. a single crystal) only. In reality, you encounter a distribution of spins (the spectrum) paired with a distribution of zero-field splitting values which make it harder to enhance the central transition. A real life example with Gd(III) can be found in
Andrin Doll, et al.
Sensitivity enhancement by population transfer in Gd(III) spin labels
Phys. Chem. Chem. Phys. 2015, DOI: 10.1039/C4CP05893C
This example explores how spin selective excitation and detection operators are used. First we start with the specification of a spin system. Since we are interested in seeing the pulses and not spin dynamics, we are going to set the spin-spin coupling to zero.
Sys.S = [1/2 1/2]; Sys.ZeemanFreq = [9.500 9.500]; % GHz Sys.J = 0; % spin-spin coupling, MHz
Next we are going to define pulses and then put them into a pulse sequence:
P90.tp = 0.016; % mus P90.Flip = pi/2; % rad P180.tp = 0.032; % mus P180.Flip = pi; % rad tau = 0.5; % mus % Sequence Exp.Sequence = {P90 tau P180 tau P180 tau}; Exp.mwFreq = 9.5; % GHz
To see how the pulses flip our spins, we are going to use selective Ŝz operators for both spins:
Exp.DetOperator = {'z1' 'z2'};
Next, let us use some selective excitation operators with Opt.ExcOperator
.
The indexing in Opt.ExcOperator
. corresponds to the pulses, e.g. the second element in Opt.ExcOperator
concerns only the second pulse.
For declaring special excitation operators, we can use the usual sop syntax.
Opt.ExcOperator = {'x1' [] 'x2'};
Now, the first pulse will use Ŝx for the first spin, the second pulse excites all spins in the system and the last pulse only excites the second spin. Run and plot it with
[TimeAxis, Signal] = spidyan(Sys,Exp,Opt); % plotting figure(1) clf plot(TimeAxis,real(Signal)); xlabel('t ({\mu}s)') axis tight ylim([-1 1]) ylabel('< S_i>') legend(Exp.DetOperator)
and you will see the spin flips.
Of course you can use transition selective excitation operators by using the sop syntax. If it is not possible to build the excitation operator with sop, it is even possible to manually design the excitation operator in matrix form, and then use that as input:
Sys.S = 3/2; % transition selective operator using sop syntax: Opt.ExcOperator = {'x(2|3)'}; % manually define x(2|3) for S = 3/2: Op = [0 0 0 0; 0 0 0.5 0; 0 0.5 0 0; 0 0 0 0]; Opt.ExcOperator = {Op};
This is example is going how to setup a rabi oscillation with spidyan using indirect dimensions. For more info on indirect dimension, please refer to the documentation. Once again, we start with the spin system: and then set up a very short 1 ns rectangular pulse.
% Spin System Sys.ZeemanFreq = 9.500; % GHz
In order to observe Rabi oscillations, we want a rectangular pulse with a constant amplitude.
We then are going to use the indirect dimension to increase the length of the pulse and detect the outcome of the experiment.
Let us start with the pulse and Exp.Sequence
% experiment setup Pulse.Type = 'rectangular'; Pulse.tp = 0.001; % pulse length in μs Pulse.Amplitude = 30; % pulse amplitude in MHz Exp.Sequence = {Pulse}; Exp.mwFreq = 9.5; % GHz
Next, detection. As detection operator we use
Exp.DetOperator = {'z1'};
and we are interested in the expectation value right after the pulse.
The easiest way is to use Exp.DetWindow
:
Exp.DetWindow = 0; % position of the detection window in μs
Since we are not providing an interval for the detection window, we are going to do single point detection.
The value 0
corresponds to the time in μs after the end of the last element in Exp.Sequence
.
For more information on how to use Exp.DetWindow
, see the documentation.
Alternatively, we could have adapted Exp.Sequence
and used Exp.DetSequence
(see here and here):
Exp.Sequence = {Pulse 0}; % append an event with length 0 μs Exp.DetSequence = [0 1]; % detect only during the the second element in Exp.Sequence
The indirect dimension, where we change the pulse length, can be added through:
Exp.nPoints = 99 % number of data points Exp.Dim1 = {'p1.tp' 0.001}; % increment length of the pulse by 1 ns per step
With this, we are going to record 99 data points, with the pulse length ranging from 1 ns to 100 ns:
spidyan(Sys,Exp); ylim([-1 1]) % to ensure right limits in the plot
Now you can investigate how the frequency of the Rabi oscillation depends on the frequency difference between the microwave pulse and the resonance frequency of the spin by changing Sys.ZeemanFreq
.
Or you can add some relaxation by adding:
Sys.T1 = 2; % longitudinal relaxation time in μs Sys.T2 = 1; % transverse relaxation time in μs Opt.Relaxation = 1; % then run spidyan(Sys,Exp,Opt);
Per call spidyan can process a single spin packet. If you want to simulate a pulse sequence with an entire spectrum, you would usually use saffron, unless you want to use special excitation operators, change the detection operator or detect between and during pulses. In such a case you can manually define an orientation loop which calls spidyan.
We start again with the spin system, but this time we are using Sys.g
instead of Sys.ZeemanFreq
.
Since we are going to manually rotate the interactions, we are going to set them up as tensors:
% Spin System Sys.S = 1/2; Sys.g = diag([2.00906 2.0062 2.0023]); Sys.A = diag([11.5 11.5 95]); Sys.Nucs = '14N';
We then determine the grid with the EasySpin functions hamsymm and sphgrid.
Symmetry = hamsymm(Sys); [phi,theta,Weights] = sphgrid(Symmetry,20);
Then we set up frequency-swept pulses, that cover the entire spectrum:
Chirp90.Type = 'quartersin/linear'; Chirp90.trise = 0.030; % edge smoothing, μs Chirp90.tp = 0.200; % pulse length, μs Chirp90.Flip = pi/2; % flip angle, rad Chirp90.Frequency = [-120 120]; % excitation band, GHz Chirp180.Type = 'quartersin/linear'; Chirp180.trise = 0.030; % smoothed edges, μs Chirp180.tp = 0.100; % pulse length, μs Chirp180.Flip = pi; % flip angle, rad Chirp180.Frequency = [-120 120]; % excitation band, GHzNext we set build the pulse sequence, set field and carrierfrequency and set up detection:
tau = 0.5; % μs Exp.Sequence = {Chirp90 tau Chirp180 tau+Chirp180.tp}; Exp.mwFreq = 9.1; % GHz Exp.Field = 324.9; % mT % Detection: Exp.DetWindow = [-0.05 0.05]; % μs Exp.DetOperator = {'+1'}; Exp.DetFreq = 9.1;
If you are not sure about where to place your pulses, keep in mind that you can call pepper with a Sys
and Exp.Field
, which is going to plot the spectrum for you.
For detection we are detecting the electron coherence in a 0.1 μs window centered around the end of the last element in Exp.Sequence
.
Since we are using frequency-swept pulses, the position of the echo is not at τ after the end of the π pulse, but we also need to take the pulse length into account.
With this we can loop over the orientations and accumulate the signal:
Signal = 0; % initialize Signal for iOrientation = 1 : numel(Weights) Sys_ = Sys; % create temporary copy of Sys R = erot(phi(iOrientation),theta(iOrientation),0); % rotation matrix Sys_.g = R'*Sys_.g*R; % rotate Sys.g Sys_.A = R'*Sys_.A*R; % rotate Sys.A [t, signal_] = spidyan(Sys_,Exp); Signal = Signal + signal_*Weights(iOrientation); % sum up signals end
After completion, the signal can be plotted:
% plotting figure(1) clf plot(t,real(Signal)) xlabel('Transient (\mus)') ylabel('Signal (arb.u.)')
If you are investigating spin dynamics, it is often not necessary to simulate entire spectra. Instead you can use any distribution of spin packets (a distribution of different resonance frequencies is needed for an echo to form).
Here we are going to cover how to loop spidyan over a Gaussian line. The resonance frequencies of the spins are centered around 9.5 GHz, with a standard deviation of 10 MHz. Spins are sampled from 9.450 to 9.550 GHz, with a sampling step 0.5 MHz, resulting in a total of 201 spin packets.
CenterFrequency = 9.5; % center frequency of Gaussian distribution, GHz GWidth = 0.01; % width of Gaussian distribution, GHz FreqStart = 9.45; % starting value for sampling FreqEnd = 9.55; % final value for sampling Sampling = 0.0005; % stepsize for sampling ZeemanFreqVec = FreqStart:Sampling:FreqEnd; % vector with resonance frequencies P = exp(-((CenterFrequency-ZeemanFreqVec)/GWidth).^2); % probabilities P = P/trapz(P); % normalization nSpinpackets = length(ZeemanFreqVec);
For each of the spin packets we now have a resonance frequency and a probability.
After setting up the pulse sequence
Pulse90.Type = 'rectangular'; Pulse90.tp = 0.025; % pulse length, mus Pulse90.Flip = pi/2; % flip angle, rad Pulse180.Type = 'rectangular'; Pulse180.tp = 0.025; % pulse length, mus Pulse180.Flip = pi/2; % flip angle, rad Exp.Sequence = {Pulse90 0.25 Pulse90 0.5}; Exp.mwFreq = 9.5; % GHz Exp.DetSequence = [0 0 0 1]; % detect the last event Exp.DetOperator = {'+1'}; % detect electron coherence Exp.DetFreq = 9.5; % GHz
Now, all that is left, is looping over the spins. In contrast to the example with the nitroxide, we do not need to rotate any tensors anymore, and can just load the resonance frequency and then plot the accumulated signal:
Signal = 0; % initialize the signal % Loop over the spin packets and sum up the traces ------------ for i = 1 : nSpinpackets Sys_.ZeemanFreq = ZeemanFreqVec(i); % Set Zeeman frequency [t, signal] = spidyan(Sys_,Exp); Signal = Signal + signal*P(i); end figure(1) clf plot(t,abs(Signal)); xlabel('t (\mus)') axis tight ylim([0 1])detect electron coherence Exp.DetFreq = 9.5; % GHz
Simultaneous irradiation (e.g. cross-polarization or decoupling) at two different frequencies is not implemented at this point, but can still be achieved:
The trick is to use Exp.mwPolarization
in combination with custom excitation operators.
This allows us to assign the real part of the wave to an excitation operator that acts on one spin or subset of spins (frequency) and the imaginary part to a second spin or subset of spins (frequency).
Look at the following example:
CombinedPulses.IQ = real(mwPulse) + imag(rfPulse); CombinedPulses.t = t_mwPulse); Opt.ExcOperator{1} = sop([1/2 1/2],'xe') + 1i*sop([1/2 1/2],'ex'); Exp.mwPolarization = 'circular'; Exp.Sequence = {CombinedPulses}; spidyan(Sys,Exp,Opt)Here, we are working with a two-spin system and (you can tell that by the way we are calling
sop
when creating the excitation operator).
We created two wave forms using pulse
(not shown) mwPulse
and rfPulse
that we then stored in the IQ
field of the CombindePulses
structure - the microwave pulse in the real part and the rf pulse as the imaginary part.
The corresponding time axis is saved to CombinedPulses.t
(both pulses must have the same length and time step).
In the next step we assign spin selective operators: The electron spin is to get excited only real part of the pulse and the nucleus by the imaginary part.
With Exp.mwPolarization = 'circular'
we are letting spidyan
(or saffron
) know that it should expect a wave form that contains a real and an imaginary part.
If the above code is now run, the simulation will use the real part of the wave form to propagate the electron spin and the imaginary part to propagate the nucleus, as if we were irradiating the spin system with two pulses at different frequencies at the same time.
If you are interested in the expectation values of Ŝx, it is usually beneficial to use Ŝ+ as the actual detection operator and then take the real part of the obtained signal. This removes artifacts at the beginning and end of the time traces that are introduced when translating a purely real signal during the signal processing.
Ŝ+ and Ŝ- are not Hermitian operators. Hence, if you use Ŝ+ as detection operator, the signal that is returned will in fact correspond to ⟨Ŝ-⟩. The real part of the signal (Ŝx) is nof affected by this. But if your data processing involves the imaginary part of the signal and you encounter frequencies with a wrong sign, you might want to consider using Ŝ- as detection operator instead.
In some simulations you may observe echoes in your timetrace that should not be there. Besides physical reasons (incomplete phase cycling, different refocusing conditions,...), such artifact echoes can arise from an aliasing effect if the spectrum (e.g. the nitroxide spectrum in the example above) is insufficiently sampled. This can easily be checked: By increasing the number of samples from your distribution the artifact echoes should move, spread further apart or vanish, while physical echoes stay in place. You can also compute the required increment for the spin packet distribution in frequency domain by an inverse Nyquist criterion. If you want to observe up to a time tmax after the first excitation event (first pulse), the frequency increment should not be larger than 1/(2 tmax). Otherwise use Monte Carlo sampling of spin packet frequencies, which will convert the artifacts to 'Monte Carlo noise' that averages with increasing number of trials.
If the signal processing fails, spidyan
returns the signal in the simulation frame.
This can happen if you try to translate Ŝz or provide a wrong frequency.
Instead of adapting down conversion and rerunning the entire simulation, you can call signalprocessing
with the returned signal and the correct down conversion frequencies.
More info here.