mirror of
https://github.com/nqrduck/quackseq-simulator.git
synced 2024-12-22 05:00:24 +00:00
196 lines
6.2 KiB
Python
196 lines
6.2 KiB
Python
import unittest
|
|
import logging
|
|
import matplotlib.pyplot as plt
|
|
from quackseq.phase_table import PhaseTable
|
|
from quackseq.pulsesequence import QuackSequence
|
|
from quackseq.event import Event
|
|
from quackseq.functions import RectFunction
|
|
from quackseq_simulator.simulator import Simulator
|
|
|
|
logging.basicConfig(level=logging.INFO)
|
|
|
|
logger = logging.getLogger(__name__)
|
|
|
|
|
|
class TestQuackSequence(unittest.TestCase):
|
|
|
|
def test_event_creation(self):
|
|
seq = QuackSequence("test - event creation")
|
|
seq.add_pulse_event("tx", "10u", 100, 90.0, RectFunction())
|
|
seq.add_blank_event("blank", "3u")
|
|
seq.add_readout_event("rx", "100u")
|
|
seq.add_blank_event("TR", "1m")
|
|
|
|
sim = Simulator()
|
|
sim.set_averages(100)
|
|
|
|
sim.settings.noise = 0
|
|
|
|
result = sim.run_sequence(seq)
|
|
self.assertIsNotNone(result)
|
|
self.assertTrue(hasattr(result, "tdx"))
|
|
self.assertTrue(hasattr(result, "tdy"))
|
|
self.assertGreater(len(result.tdx[0]), 0)
|
|
self.assertGreater(len(result.tdy[0]), 0)
|
|
|
|
logger.info("Plotting imaginary part")
|
|
plt.plot(result.tdx[0], result.tdy[0].imag, label="imaginary")
|
|
logger.info("Plotting real part")
|
|
plt.plot(result.tdx[0], result.tdy[0].real, label="real")
|
|
plt.plot(result.tdx[0], abs(result.tdy[0]), label="abs")
|
|
plt.legend()
|
|
plt.show()
|
|
|
|
def test_simulation_run_sequence(self):
|
|
seq = QuackSequence("test - simulation run sequence")
|
|
|
|
tx = Event("tx", "10u", seq)
|
|
seq.add_event(tx)
|
|
seq.set_tx_amplitude(tx, 100)
|
|
seq.set_tx_phase(tx, 0)
|
|
|
|
rect = RectFunction()
|
|
seq.set_tx_shape(tx, rect)
|
|
|
|
blank = Event("blank", "3u", seq)
|
|
seq.add_event(blank)
|
|
|
|
rx = Event("rx", "100u", seq)
|
|
seq.set_rx(rx, True)
|
|
seq.add_event(rx)
|
|
|
|
TR = Event("TR", "1m", seq)
|
|
seq.add_event(TR)
|
|
|
|
sim = Simulator()
|
|
sim.set_averages(100)
|
|
|
|
result = sim.run_sequence(seq)
|
|
self.assertIsNotNone(result)
|
|
self.assertTrue(hasattr(result, "tdx"))
|
|
self.assertTrue(hasattr(result, "tdy"))
|
|
self.assertGreater(len(result.tdx[0]), 0)
|
|
self.assertGreater(len(result.tdy[0]), 0)
|
|
|
|
plt.plot(result.tdx[0], abs(result.tdy[0]))
|
|
plt.show()
|
|
|
|
def test_phase_array_generation(self):
|
|
seq = QuackSequence("test - phase array generation")
|
|
|
|
tx = Event("tx", "10u", seq)
|
|
seq.add_event(tx)
|
|
seq.set_tx_amplitude(tx, 100)
|
|
seq.set_tx_phase(tx, 0)
|
|
seq.set_tx_n_phase_cycles(tx, 2)
|
|
seq.set_tx_phase_cycle_group(tx, 0)
|
|
|
|
rect = RectFunction()
|
|
seq.set_tx_shape(tx, rect)
|
|
|
|
tx2 = Event("tx2", "10u", seq)
|
|
seq.add_event(tx2)
|
|
seq.set_tx_amplitude(tx2, 100)
|
|
seq.set_tx_phase(tx2, 1)
|
|
seq.set_tx_n_phase_cycles(tx2, 2)
|
|
seq.set_tx_phase_cycle_group(tx2, 1)
|
|
|
|
tx3 = Event("tx3", "10u", seq)
|
|
seq.add_event(tx3)
|
|
seq.set_tx_amplitude(tx3, 100)
|
|
seq.set_tx_phase(tx3, 2)
|
|
seq.set_tx_n_phase_cycles(tx3, 2)
|
|
seq.set_tx_phase_cycle_group(tx3, 3)
|
|
|
|
sim = Simulator()
|
|
sim.set_averages(100)
|
|
|
|
result = sim.run_sequence(seq)
|
|
|
|
plt.plot(result.tdx[0], abs(result.tdy[0]))
|
|
plt.plot(result.tdx[1], abs(result.tdy[1]))
|
|
plt.show()
|
|
|
|
phase_table = PhaseTable(seq)
|
|
|
|
logger.info(phase_table.n_phase_cycles)
|
|
self.assertEqual(phase_table.n_phase_cycles, 8)
|
|
self.assertEqual(phase_table.n_parameters, 3)
|
|
|
|
def test_phase_cycling(self):
|
|
|
|
seq = QuackSequence("test - phase cycling")
|
|
|
|
# Simple spin echo sequence with phase cycling.
|
|
# Create the first 90 degree pulse
|
|
pi_half = Event("pi_half", "6u", seq)
|
|
seq.add_event(pi_half)
|
|
seq.set_tx_amplitude(pi_half, 100)
|
|
seq.set_tx_phase(pi_half, 0)
|
|
seq.set_tx_n_phase_cycles(pi_half, 4) # 0 90 180 270
|
|
seq.set_tx_phase_cycle_group(pi_half, 0)
|
|
|
|
# Wait for TE/2 (approx)
|
|
seq.add_blank_event("te/2", "150u")
|
|
|
|
# Create the 180 degree pulse
|
|
pi = Event("pi", "12u", seq)
|
|
seq.add_event(pi)
|
|
seq.set_tx_amplitude(pi, 100)
|
|
seq.set_tx_phase(pi, 90)
|
|
seq.set_tx_n_phase_cycles(pi, 1)
|
|
seq.set_tx_phase_cycle_group(pi, 0)
|
|
|
|
# Wait another blanking time
|
|
seq.add_blank_event("blank", "100u")
|
|
|
|
result = Simulator().run_sequence(seq)
|
|
|
|
# Plotting the results
|
|
plt.title("Phase cycling")
|
|
logger.info(f"Number of data sets {len(result.tdy)}")
|
|
plt.plot(result.tdx[0], result.tdy[0].real, label="pc 1")
|
|
plt.plot(result.tdx[1], result.tdy[1].real, label="pc 2")
|
|
plt.plot(result.tdx[2], result.tdy[2].real, label="pc 3")
|
|
plt.plot(result.tdx[3], result.tdy[3].real, label="pc 4")
|
|
plt.legend()
|
|
plt.show()
|
|
|
|
rx = Event("rx", "100u", seq)
|
|
seq.add_event(rx)
|
|
seq.set_rx(rx, True)
|
|
|
|
# Readout scheme for phase cycling TX pulses have the scheme 0 90 180 270 for the first pulse and 180 always for the second pulse
|
|
readout_scheme = [[1, 0], [1, 90], [1, 180], [1, 270]]
|
|
|
|
seq.set_rx_readout_scheme(rx, readout_scheme)
|
|
|
|
result = Simulator().run_sequence(seq)
|
|
|
|
# Plotting the results
|
|
plt.title("Phase cycling")
|
|
logger.info(f"Number of data sets {len(result.tdy)}")
|
|
plt.plot(result.tdx[0], result.tdy[0].real, label="pc 1 real")
|
|
plt.plot(result.tdx[0], result.tdy[0].imag, label="pc 2 imag")
|
|
plt.legend()
|
|
plt.show()
|
|
plt.plot(result.tdx[1], result.tdy[1].real, label="pc 2")
|
|
plt.plot(result.tdx[1], result.tdy[1].imag, label="pc 2")
|
|
plt.legend()
|
|
plt.show()
|
|
plt.plot(result.tdx[2], result.tdy[2].real, label="pc 3")
|
|
plt.plot(result.tdx[2], result.tdy[2].imag, label="pc 3")
|
|
plt.legend()
|
|
plt.show()
|
|
plt.plot(result.tdx[3], result.tdy[3].real, label="pc 4")
|
|
plt.plot(result.tdx[3], result.tdy[3].imag, label="pc 4")
|
|
plt.legend()
|
|
plt.show()
|
|
plt.plot(result.tdx[4], abs(result.tdy[4]), label="Phase Cycling")
|
|
plt.legend()
|
|
plt.show()
|
|
# seq.add_event(rx)
|
|
|
|
|
|
if __name__ == "__main__":
|
|
unittest.main()
|