362 lines
15 KiB
Python
362 lines
15 KiB
Python
import os
|
|
from math import pi
|
|
from io import BytesIO
|
|
import base64
|
|
|
|
from bokeh.models.axes import LinearAxis
|
|
from django.shortcuts import render
|
|
from django.utils.safestring import mark_safe
|
|
from django.conf import settings
|
|
import re
|
|
|
|
|
|
import numpy as np
|
|
from bokeh.plotting import figure
|
|
from bokeh.embed import components
|
|
from bokeh.models import Label, Node, MathML, Range1d, Span, ColumnDataSource, DataTable, TableColumn
|
|
|
|
from isotopapp.models import Isotope, FieldProfile
|
|
# Create your views here.
|
|
|
|
def home(request):
|
|
isotopes = [i for i in Isotope.objects.all() if (i.gamma != 0 or i.stable)]
|
|
return render(request, 'home.html', {'isotopes': [[f"{i.n_nucleons}{i.symbol}", mark_safe(f"<sup>{i.n_nucleons}</sup>{i.symbol}")] for i in isotopes],})
|
|
|
|
def sfg(request):
|
|
isotopes = [i for i in Isotope.objects.all() if (i.gamma != 0 or i.stable)]
|
|
magnets = FieldProfile.objects.all()
|
|
return render(request, 'sfg.html', {
|
|
'isotopes': [[f"{i.n_nucleons}{i.symbol}", mark_safe(f"<sup>{i.n_nucleons}</sup>{i.symbol}")] for i in isotopes],
|
|
'magnets': magnets,
|
|
})
|
|
|
|
def diffusion_form(request):
|
|
return render(request, 'diffusion_form.html', {})
|
|
|
|
|
|
def extract_isotope_parts(isotope_str):
|
|
"""Extracts the number and element from an isotope string (e.g., '23Na')."""
|
|
print(isotope_str)
|
|
match = re.match(r"(\d+)([A-Za-z]+)", isotope_str)
|
|
if not match:
|
|
raise ValueError("Invalid isotope format")
|
|
return int(match.group(1)), match.group(2)
|
|
|
|
def isotope_info(isotope, field):
|
|
f_larmor = field*isotope.gamma
|
|
if isotope.spin_quantum_number.as_integer_ratio()[1]==1:
|
|
spin = int(isotope.spin_quantum_number)
|
|
else:
|
|
spin = mark_safe(f"<sup>{isotope.spin_quantum_number.as_integer_ratio()[0]}</sup>⁄<sub>{isotope.spin_quantum_number.as_integer_ratio()[1]}</sub>")
|
|
return [isotope.n_nucleons,
|
|
mark_safe(f"<sup>{isotope.n_nucleons}</sup>{isotope.symbol}"),
|
|
isotope.name.capitalize(),
|
|
f"{f_larmor:.3f}",
|
|
spin,
|
|
f"{isotope.natural_abundance:.3f}",
|
|
f"{isotope.gamma:.5f}",
|
|
f"{2*pi*isotope.gamma*1e6:.5e}",
|
|
f"{relative_sensitivity(isotope):.4f}"]
|
|
|
|
def relative_sensitivity(iso):
|
|
riso = Isotope.objects.filter(symbol="H", n_nucleons="1").get()
|
|
i1 = riso.spin_quantum_number
|
|
g1 = riso.gamma
|
|
ab1 = riso.natural_abundance
|
|
i2 = float(iso.spin_quantum_number)
|
|
g2 = float(iso.gamma)
|
|
ab2 = float(iso.natural_abundance)
|
|
rel_sens = i2*(i2+1)*g2**3*ab2 / (i1*(i1+1)*g1**3*ab1)
|
|
print(rel_sens)
|
|
return rel_sens*100
|
|
|
|
def result(request):
|
|
n1, element1 = extract_isotope_parts(request.GET.get('isotope1'))
|
|
isotope1 = Isotope.objects.filter(symbol=element1, n_nucleons=n1).get()
|
|
|
|
freq = float(request.GET.get('freq'))
|
|
|
|
field_T = freq / isotope1.gamma
|
|
|
|
script = None
|
|
div = None
|
|
|
|
print(request.GET)
|
|
if request.GET.get('range_search') == "":
|
|
close_isotopes = []
|
|
freq_range = float(request.GET.get('freq_range'))
|
|
Isotope.objects.filter()
|
|
# calculate the frequency for all isotopes and compile a list of close by isotopes
|
|
for isotope in Isotope.objects.all():
|
|
if isotope.gamma == 0: continue
|
|
if not isotope.stable: continue
|
|
f_Larmor = field_T*isotope.gamma
|
|
if abs(f_Larmor - freq) <= freq_range:
|
|
close_isotopes.append(isotope_info(isotope, field_T))
|
|
ans = sorted(close_isotopes, key=lambda x: -float(x[3]))
|
|
div = f"Field B<sub>0</sub>: {field_T:.3g} T"
|
|
script = ""
|
|
|
|
elif request.GET.get('gradient_search') == "":
|
|
sample_diameter = 5e-3
|
|
gradient = float(request.GET.get('gradient'))
|
|
# create a plot (bokeh)
|
|
|
|
plot = figure(outer_width=500, outer_height=500, match_aspect=True)
|
|
plot.ellipse(x=[0], y=[0], width=5, height=5, color="#D5D9FF", alpha=0.8, line_width=1, line_color="black")
|
|
plot.ellipse(x=[0], y= [2.5], width=5, height=5, color="#D5D9FF", alpha=0.4, line_width=1, line_color="black")
|
|
plot.ellipse(x=[0], y=[-2.5], width=5, height=5, color="#D5D9FF", alpha=0.4, line_width=1, line_color="black")
|
|
plot.xaxis[0].axis_label = 'x / mm'
|
|
plot.yaxis[0].axis_label = 'z / mm'
|
|
|
|
close_isotopes = []
|
|
for isotope in Isotope.objects.all():
|
|
if isotope.gamma == 0: continue
|
|
if not isotope.stable: continue
|
|
z = (freq/isotope.gamma-field_T)/gradient
|
|
if abs(z) <= sample_diameter:
|
|
i_info = isotope_info(isotope, field_T)
|
|
#i_info[3] = f"{z*1e3:.1f} mm"
|
|
close_isotopes.append(i_info)
|
|
plot.rect(x=[0], y=[z*1e3], width=5, height=0.2, color="black", alpha=0.6)
|
|
label = Label(x=2.6, y=z*1e3, text=f"{isotope.n_nucleons}{isotope.symbol}", text_baseline="middle", text_align="left", text_font_size="16pt")
|
|
plot.add_layout(label)
|
|
|
|
frame_left = Node(target="frame", symbol="left", offset=5)
|
|
frame_bottom = Node(target="frame", symbol="bottom", offset=-5)
|
|
citation = Label(
|
|
x=frame_left,
|
|
y=frame_bottom,
|
|
anchor="bottom_left",
|
|
text=f"{isotope1.n_nucleons}{isotope1.symbol}: {freq:.1f} MHz\ng={gradient:.1f} T/m\n5 mm sample dia.",
|
|
padding=5,
|
|
border_radius=5,
|
|
border_line_color="#D5D9FF", background_fill_color="white",
|
|
)
|
|
plot.add_layout(citation)
|
|
|
|
# boke plot
|
|
script, div = components(plot)
|
|
ans = sorted(close_isotopes, key=lambda x: float(x[3]))
|
|
|
|
elif request.GET.get('transform') == "":
|
|
n2, element2 = extract_isotope_parts(request.GET.get('isotope2'))
|
|
isotope2 = Isotope.objects.filter(symbol=element2, n_nucleons=n2).get()
|
|
#isotope_info(isotope2, field_T)
|
|
ans = [isotope_info(isotope2, field_T), isotope_info(isotope1, field_T) ]
|
|
div = f"Field B<sub>0</sub>: {field_T:.3f} T"
|
|
script = ""
|
|
else:
|
|
ans = []
|
|
return render(request, 'result.html', {'ans': ans, 'script': script, 'div': div})
|
|
|
|
def position(request):
|
|
print(request)
|
|
n1, element1 = extract_isotope_parts(request.GET.get('isotope1'))
|
|
isotope1 = Isotope.objects.filter(symbol=element1, n_nucleons=n1).get()
|
|
probe_length = float(request.GET.get('probe_length'))
|
|
print(os.path.abspath("."))
|
|
magnet_id = request.GET.get('magnet')
|
|
try:
|
|
if magnet_id.isdigit():
|
|
magnet = FieldProfile.objects.get(pk=magnet_id)
|
|
else:
|
|
magnet = FieldProfile.objects.get(name=magnet_id)
|
|
|
|
step = magnet.step
|
|
cryo_length = magnet.cryo_length
|
|
if magnet.reverse_offset:
|
|
offset = -(cryo_length / step - probe_length / step)
|
|
else:
|
|
offset = cryo_length / step - probe_length / step
|
|
file_path = magnet.file.path
|
|
except (FieldProfile.DoesNotExist, ValueError):
|
|
if magnet_id == "oxford.dat":
|
|
step = 5e-3/6
|
|
cryo_length = 1113.0
|
|
offset = cryo_length / step - probe_length / step
|
|
file_path = os.path.join(settings.BASE_DIR, "magnet-profiles", "oxford.dat")
|
|
elif magnet_id == "magnex.dat":
|
|
step = 1e-3
|
|
cryo_length = 1262.2
|
|
offset = -(cryo_length / step - probe_length / step) # other direction
|
|
file_path = os.path.join(settings.BASE_DIR, "magnet-profiles", "magnex.dat")
|
|
else:
|
|
return render(request, 'home.html', {'error': f"Magnet profile {magnet_id} not found."})
|
|
|
|
data = np.loadtxt(file_path)
|
|
_z_coords = data[:,0]
|
|
_fields = data[:,1]
|
|
_gradients = data[:,2]
|
|
|
|
freq = float(request.GET.get('freq'))
|
|
#gradient = float(request.GET.get('gradient'))
|
|
|
|
field_T = freq / isotope1.gamma
|
|
|
|
sample_diameter = 5e-3
|
|
# create a plot (bokeh)
|
|
|
|
TOOLTIPS = [
|
|
#("index", "$index"),
|
|
("(x,y)", "($x, $y)"),
|
|
]
|
|
plot = figure(outer_width=500, outer_height=500, match_aspect=False, tooltips=TOOLTIPS)
|
|
frame_left = Node(target="frame", symbol="left", offset=5)
|
|
frame_bottom = Node(target="frame", symbol="bottom", offset=5)
|
|
frame_top = Node(target="frame", symbol="top", offset=5)
|
|
frame_right = Node(target="frame", symbol="right", offset=5)
|
|
#plot.ellipse(x=[0], y=[0], width=5, height=5, color="#D5D9FF", alpha=0.8, line_width=1, line_color="black")
|
|
#plot.ellipse(x=[0], y= [2.5], width=5, height=5, color="#D5D9FF", alpha=0.4, line_width=1, line_color="black")
|
|
#plot.ellipse(x=[0], y=[-2.5], width=5, height=5, color="#D5D9FF", alpha=0.4, line_width=1, line_color="black")
|
|
plot.xaxis[0].axis_label = 'steps'
|
|
plot.yaxis[0].axis_label = 'B0 in T'
|
|
plot.line(x=_z_coords, y=_fields, color="navy", line_width=2)
|
|
plot.yaxis[0].axis_label_text_color = "navy"
|
|
plot.y_range = Range1d(-10, 10)
|
|
|
|
|
|
plot.extra_y_ranges['gradient'] = Range1d(-200, 200)
|
|
plot.line(x=_z_coords, y=_gradients, color="crimson", line_width=2 , y_range_name="gradient")
|
|
ax2 = LinearAxis(y_range_name="gradient", axis_label="g in T/m")
|
|
ax2.axis_label_text_color = "crimson"
|
|
plot.add_layout(ax2, 'left')
|
|
|
|
# find fields where _fields == field:
|
|
pos_threshold_indices = np.where(np.diff(np.sign(_fields + field_T)) != 0)[0]
|
|
neg_threshold_indices = np.where(np.diff(np.sign(_fields - field_T)) != 0)[0]
|
|
indices = np.concatenate((pos_threshold_indices, neg_threshold_indices))
|
|
|
|
|
|
for p in indices:
|
|
#vline = Span(location=_z_coords[p]*1e3, dimension='height', line_color='gray', line_width=1)
|
|
#plot.renderers.extend([vline])
|
|
#zpos = _z_coords[p] + (1117 - probe_length)
|
|
label = Label(
|
|
x=_z_coords[p]+25,
|
|
y=-(np.sign(_z_coords[p])*30) + (_gradients[p]),
|
|
text=f"{_z_coords[p]:.0f}\n{_gradients[p]:.1f}T/m",
|
|
text_baseline="middle",
|
|
text_align="center", text_font_size="11pt",
|
|
border_line_color="#D5D9FF",
|
|
background_fill_color="white",
|
|
background_fill_alpha=0.7,
|
|
y_range_name="gradient",
|
|
)
|
|
plot.add_layout(label)
|
|
|
|
|
|
#plot.circle(x=_z_coords[indices]*1e3, y=_fields[indices], radius=2, color="gray")
|
|
plot.circle(x=_z_coords[indices], y=_gradients[indices], radius=2, color="gray",y_range_name="gradient" )
|
|
# plot the B0 field hlines
|
|
hline0 = Span(location=0, dimension='width', line_color='black', line_width=4)
|
|
hline1 = Span(location=field_T, dimension='width', line_color='gray', line_width=2)
|
|
label1 = Label(x=frame_right, y=field_T, text=f"{field_T:.1f}T/m",
|
|
text_baseline="middle")
|
|
hline2 = Span(location=-field_T, dimension='width', line_color='gray', line_width=2)
|
|
plot.renderers.extend([hline0, hline1, hline2, label1])
|
|
|
|
for i in indices:
|
|
plot.line(x=[_z_coords[i], _z_coords[i]], y=[_fields[i], _gradients[i]/20. ] , color="gray", line_dash="dashed" )
|
|
|
|
close_isotopes = []
|
|
for isotope in Isotope.objects.all():
|
|
if isotope.gamma == 0: continue
|
|
if not isotope.stable: continue
|
|
z = (freq/isotope.gamma-field_T)/_gradients[0]
|
|
if abs(z) <= sample_diameter:
|
|
i_info = isotope_info(isotope, field_T)
|
|
#i_info[3] = f"{z*1e3:.1f} mm"
|
|
close_isotopes.append(i_info)
|
|
#plot.rect(x=[0], y=[z*1e3], width=5, height=0.2, color="black", alpha=0.6)
|
|
#label = Label(x=2.6, y=z*1e3, text=f"{isotope.n_nucleons}{isotope.symbol}", text_baseline="middle", text_align="left", text_font_size="16pt")
|
|
#plot.add_layout(label)
|
|
|
|
text_label = f"{isotope1.n_nucleons}{isotope1.symbol}: {freq:.1f} MHz\nProbe position:\n"
|
|
for i in indices:
|
|
text_label += f"{_gradients[i]:.2f}T/m: {abs(round(offset + _z_coords[i]))}\n"
|
|
|
|
citation = Label(
|
|
x=frame_left,
|
|
y=frame_bottom,
|
|
anchor="bottom_left",
|
|
text=text_label.strip(),
|
|
padding=5,
|
|
border_radius=2,
|
|
border_line_color="#D5D9FF", background_fill_color="white",
|
|
)
|
|
plot.add_layout(citation)
|
|
|
|
bokeh_table_data = {"pos":[_z_coords[p] for p in indices] , "gradients":[_gradients[p] for p in indices]}
|
|
source = ColumnDataSource(bokeh_table_data)
|
|
|
|
columns = [
|
|
TableColumn(field="pos", title="Position"),
|
|
TableColumn(field="gradients", title="Gradient / T/m"),
|
|
]
|
|
data_table = DataTable(source=source, columns=columns, width=400, height=280)
|
|
|
|
script_plot, div_plot = components(plot)
|
|
|
|
script_table, div_table = components(data_table)
|
|
ans = sorted(close_isotopes, key=lambda x: float(x[3]))
|
|
print(div_plot)
|
|
return render(request, 'posiiton.html', {'ans': ans,
|
|
'script_plot': script_plot,
|
|
'div_plot': div_plot,
|
|
"script_table": script_table,
|
|
"div_table": div_table})
|
|
|
|
|
|
def diff(c1, c2, temp_Kelvin):
|
|
return 1e-9*np.exp(c1+c2*1000/temp_Kelvin)
|
|
|
|
def diffusion(request):
|
|
standard = request.GET.get("standard")
|
|
temperature = float(request.GET.get("temperature"))
|
|
if standard == "h2o":
|
|
diffusion_constant = 1.635e-8*(temperature/215.05-1)**2.063
|
|
standard_name = "H<sub>2</sub>O"
|
|
elif standard == "tetradecane":
|
|
diffusion_constant = diff(5.6286, -1.8726, temperature)
|
|
standard_name = "Tetradecane"
|
|
elif standard == "cyclohexane":
|
|
diffusion_constant = diff(5.8765, -1.6477, temperature)
|
|
standard_name = "Cyclohexane"
|
|
elif standard == "dioxane":
|
|
diffusion_constant = diff(5.7064, -1.6754, temperature)
|
|
standard_name = "Dioxane"
|
|
elif standard == "dodecane":
|
|
diffusion_constant = diff(5.3193 , -1.6483, temperature)
|
|
standard_name = "Dodecane"
|
|
elif standard == "dmso":
|
|
diffusion_constant = diff(5.6991 , -1.7927, temperature)
|
|
standard_name = "DMSO"
|
|
elif standard == "pentanol":
|
|
diffusion_constant = diff(8.4847, -2.9018, temperature)
|
|
standard_name = "Pentanol"
|
|
elif standard == "squalene":
|
|
diffusion_constant = 1e-12*np.exp(10.303 - 992.2/(temperature-155.09))
|
|
standard_name = "Squalene"
|
|
caption = "http://dx.doi.org/10.1021/acs.jced.5b00246"
|
|
elif standard == "ehb":
|
|
diffusion_constant = 1e-12*np.exp(9.8256 - 641.0/(temperature-163.00))
|
|
standard_name = "EHB"
|
|
caption = "http://dx.doi.org/10.1021/acs.jced.5b00246"
|
|
elif standard == "otp":
|
|
poly_otp = np.array([ -0.81665639, 7.85049752, -29.72962161, 51.26737628,
|
|
-42.33608564])
|
|
# from 273 to 407K
|
|
if 273 < temperature < 407:
|
|
diffusion_constant = 10**np.polyval(poly_otp, 1e3/temperature)
|
|
else:
|
|
diffusion_constant = 0
|
|
standard_name = "Orthoterphenyl"
|
|
|
|
|
|
|
|
return render(request, 'diffusion.html', {'standard_name': standard_name,
|
|
'temperature_kelvin': temperature,
|
|
'temperature_celsius': temperature-273.15,
|
|
'diffusion_constant': diffusion_constant})
|