gao/analyseScript
1
0
Fork 0
Code Issues Pull Requests repo.project_board Releases Wiki Activity

Merge branch 'add_rotation' of https://git.physi.uni-heidelberg.de/gao/analyseScript into origin/add_rotation

Browse Source
This commit is contained in:
Jianshun Gao 2023-09-08 14:30:58 +02:00
commit 1b65934ede
2 changed files with 1624 additions and 4 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1599
20230630_Data_Analysis copy.ipynb Normal file
View File

File diff suppressed because one or more lines are too long

29
Analyser/FitAnalyser.py
View File

@ -26,7 +26,7 @@ from scipy.special import gamma as gamfcn
from scipy.special import wofz
from scipy.optimize import curve_fit
from scipy.interpolate import CubicSpline
from scipy.ndimage import gaussian_filter
from scipy.ndimage import gaussian_filter, rotate
import xarray as xr
@ -159,9 +159,23 @@ def polylog2_2d(x, y=0.0, centerx=0.0, centery=0.0, amplitude=1.0, sigmax=1.0, s
## Approximation of the polylog function with 2D gaussian as argument. -> discribes the thermal part of the cloud
return amplitude/1.643 * polylog_int(np.exp( -((x-centerx)**2/(2 * sigmax**2))-((y-centery)**2/( 2 * sigmay**2)) ))
def rotate_coord(x,y, rot_angle):
rot_angle *= 2*np.pi/360
x_og = x
x = x_og * np.cos(rot_angle) + y * np.sin(rot_angle)
y = -x_og * np.sin(rot_angle) + y * np.cos(rot_angle)
return x, y
def density_profile_BEC_2d(x, y=0.0, amp_bec=1.0, amp_th=1.0, x0_bec=0.0, y0_bec=0.0, x0_th=0.0, y0_th=0.0,
sigmax_bec=1.0, sigmay_bec=1.0, sigma_th=1.0):
sigmax_bec=1.0, sigmay_bec=1.0, sigma_th=1.0, rot_angle=None):
if rot_angle is not None:
x, y = rotate_coord(x,y, rot_angle)
x0_bec, y0_bec = rotate_coord(x0_bec, y0_bec, rot_angle)
x0_th, y0_th = rotate_coord(x0_th, y0_th, rot_angle)
return ThomasFermi_2d(x=x, y=y, centerx=x0_bec, centery=y0_bec,
amplitude=amp_bec, sigmax=sigmax_bec, sigmay=sigmay_bec
) + polylog2_2d(x=x, y=y, centerx=x0_th, centery=y0_th,
@ -402,12 +416,14 @@ class DensityProfileBEC2dModel(lmfit.Model):
self.set_param_hint('sigmax_bec', min=0)
self.set_param_hint('sigmay_bec', min=0)
self.set_param_hint('sigma_th', min=0)
self.set_param_hint('rot_angle', min=0, max=360)
self.set_param_hint('atom_number_bec', expr=f'{self.prefix}amp_bec / 5 * 2 * 3.14159265359 * {self.prefix}sigmax_bec * {self.prefix}sigmay_bec')
self.set_param_hint('atom_number_th', expr=f'{self.prefix}amp_th * 2 * 3.14159265359 * 1.20206 / 1.643 * {self.prefix}sigma_th * {self.prefix}sigma_th')
self.set_param_hint('condensate_fraction', expr=f'{self.prefix}atom_number_bec / ({self.prefix}atom_number_bec + {self.prefix}atom_number_th)')
def guess(self, data, x, y, pre_check=False, post_check=False, **kwargs):
def guess(self, data, x, y, rot_angle=0, vary_rot=False, pre_check=False, post_check=False, **kwargs):
"""Estimate and create initial model parameters for 2d bimodal fit, by doing a 1d bimodal fit along an integrated slice of the image
:param data: Flattened 2d array, in form [a_00, a_10, a_20, ..., a_01, a_02, .. ,a_XY] with a_xy, x_dim=X, y_dim=Y
@ -437,6 +453,9 @@ class DensityProfileBEC2dModel(lmfit.Model):
data = np.reshape(data, (y_width, x_width))
data = data.T
if rot_angle != 0:
data = rotate(data, rot_angle, reshape=False)
shape = np.shape(data)
if self.is_debug:
print(f'shape: {shape}')
@ -582,6 +601,8 @@ class DensityProfileBEC2dModel(lmfit.Model):
else:
print('Error in small width BEC recogintion, s_width_ind should be 0 or 1')
params[f'{self.prefix}rot_angle'].set(value=-rot_angle, min=-rot_angle-30, max=-rot_angle+30, vary=vary_rot)
if self.is_debug:
print('')
print('Init Params')