diff --git a/iexcode/instruments/VLS_PGM.py b/iexcode/instruments/VLS_PGM.py index 4aae5e76c4331f820c8eb345b0383ab1e1b73db1..cb1f34624b4e2548b09f2bebe3a4fdba0fb4abdc 100644 --- a/iexcode/instruments/VLS_PGM.py +++ b/iexcode/instruments/VLS_PGM.py @@ -13,6 +13,8 @@ from iexcode.instruments.IEX_BL_config import * from iexcode.instruments.utilities import print_warning_message,read_dict from iexcode.instruments.scanRecord import * + + ############################################################################################################## ################################ mono ############################## ############################################################################################################## @@ -84,21 +86,29 @@ def _mono_mirror_num(): """returns the current mirror num""" return caget('29idmonoMIR_TYPE_MON') -def _mono_pvs(grt_num=None,mir_num=None): +def mono_pvs(): + """ + returns a dictionary with the pvs for energy, energy_sp and grt_density + """ + d={'energy':"29idmono:ENERGY_MON", + 'energy_sp':"29idmono:ENERGY_SP", + 'grt_density':"29idmono:GRT_DENSITY", + } + return d + +def _mono_pvs_extended(grt_num=None,mir_num=None): """ returns a dictionary with all the mono pvs """ + if grt_num is None: grt_num = _mono_grating_num() if mir_num is None: mir_num = _mono_mirror_num() ext=_mono_extensions() - - d={ - 'energy':"29idmono:ENERGY_MON", - 'energy_sp':"29idmono:ENERGY_SP", - 'grt_density':"29idmono:GRT_DENSITY", + d = mono_pvs() + d.update({ 'grt_offset':'29idmonoGRT:P_OFFSETS.'+ext[grt_num], 'grt_b2':'29idmonoGRT:B2_CALC.'+ext[grt_num], 'grt_pos':'29idmonoGRT:X_DEF_POS.'+ext[grt_num], @@ -136,10 +146,8 @@ def _mono_pvs(grt_num=None,mir_num=None): 'mir_pos':'29idmonoMIR:X_DEF_POS.'+ext[mir_num], 'mono_stop':"29idmono:STOP_CMD.PROC", 'cff':'29idmono:CC_MON', - 'arm':'29idmono:PARAMETER.G' - - - } + 'arm':'29idmono:PARAMETER.G', + }) return d def mono_get_all_extended(verbose=False): @@ -156,7 +164,7 @@ def mono_get_all_extended(verbose=False): ext=_mono_extensions() - pvs = _mono_pvs() + pvs = _mono_pvs_extended() vals={ 'grt_num':grt_num, 'mir_num':mir_num, @@ -220,7 +228,7 @@ def mono_energy_set(hv_eV,verbose=True): Previously: SetMono """ hv_min, hv_max = mono_energy_range() - pv = _mono_pvs()['energy_sp'] + pv = _mono_pvs_extended()['energy_sp'] if hv_min <= hv_eV <= hv_max: caput(pv,hv_eV,wait=True,timeout=60) time.sleep(2.5) @@ -241,8 +249,8 @@ def mono_scan_pvs(): """ returns the pvs to be used in scanning """ - val_pv = _mono_pvs()['energy_sp'] - rbv_pv = _mono_pvs()['energy'] + val_pv = _mono_pvs_extended()['energy_sp'] + rbv_pv = _mono_pvs_extended()['energy'] return val_pv, rbv_pv def mono_status_get(): @@ -255,7 +263,7 @@ def mono_status_get(): Previously Mono_Status """ - pvs = _mono_pvs() + pvs = _mono_pvs_extended() mir_P_status = int(caget(pvs['mir_P_status'])) grt_P_status = int(caget(pvs['grt_P_status'])) mir_X_status = int(caget(pvs['mir_X_status'])) @@ -301,7 +309,7 @@ def _mono_optic_change(optic,name, verbose=True): print_warning_message(name+' not a valid '+optic) return - pvs = _mono_pvs () + pvs = _mono_pvs_extended () if optic == 'grating': type_get = mono_grating_get() type_sp = pvs['grt_type_sp'] @@ -339,7 +347,7 @@ def mono_resest_pitch(): resets the MIR and GRT pitch interlocks after a following error or limit Previoulsy part of SetMono """ - pv=_mono_pvs() + pv=_mono_pvs_extended() caput(pv['grt_P_stop'],1) time.sleep(1) caput(pv['mir_P_stop'],1) @@ -352,7 +360,7 @@ def mono_kill(): Previously: Kill_Mono """ - pv=_mono_pvs() + pv=_mono_pvs_extended() caput(pv['grt_P_kill'],1) caput(pv['mir_P_kill'],1) caput(pv['grt_X_kill'],1) @@ -368,7 +376,7 @@ def mono_stop(): Previously: Stop_Mono """ - pv=_mono_pvs() + pv=_mono_pvs_extended() caput(pv['mono_stop'],1) time.sleep(5) @@ -378,7 +386,7 @@ def mono_enable(): Previously: Enable_Mono """ - pv=_mono_pvs() + pv=_mono_pvs_extended() caput(pv['grt_P_enable'],1) caput(pv['mir_P_enable'],1) caput(pv['grt_X_enable'],1) @@ -392,7 +400,7 @@ def mono_limits_reset(): Previously: Reset_Mono_Limits """ - pv=_mono_pvs()['energy_sp'] + pv=_mono_pvs_extended()['energy_sp'] caput(pv+".DRVL",200) caput(pv+".LOW",200) caput(pv+".LOLO",200) @@ -404,7 +412,7 @@ def mono_parameters_pv(grt_num=None,mir_num=None): """ returns dictionary with mono_parameter for current grating/mirror """ - pvs=_mono_pvs(grt_num,mir_num) + pvs=_mono_pvs_extended(grt_num,mir_num) d={ @@ -617,7 +625,7 @@ def mono_angles_set(alpha,beta): #JM modified to monitor the ready, moving seque """ alpha=alpha*1.0 beta=beta*1.0 - pvs = _mono_pvs() + pvs = _mono_pvs_extended() #don't wait during move so they move in parallel (less likely to crash) caput(pvs['mir_P_sp'],alpha) @@ -636,7 +644,7 @@ def mono_pink_beam(): Previously: Mono_pinkbeam """ mono_angles_set(0,0) - pvs=_mono_pvs() + pvs=_mono_pvs_extended() mono_motor_move('MIR:X',-52,wait=False) mono_motor_move('GRT:X',210,wait=True) @@ -654,7 +662,7 @@ def mono_CC(): caput("29idmono:ENERGY_SP", 440) def mono_grating_offset_get(grt_num=None): - pvs=_mono_pvs(grt_num,None) + pvs=_mono_pvs_extended(grt_num,None) return caget(pvs['grt_offset']) def mono_grating_offset_set(val,grt_num=None): @@ -663,12 +671,12 @@ def mono_grating_offset_set(val,grt_num=None): Previously: Mono_Set_GRT0 """ - pvs=_mono_pvs(grt_num,None) + pvs=_mono_pvs_extended(grt_num,None) caput(pvs['grt_offset'],val) mono_get_all_extended(verbose=True) def mono_mirror_offset_get(mir_num=None): - pvs=_mono_pvs(None,mir_num) + pvs=_mono_pvs_extended(None,mir_num) return caget(pvs['mir_offset']) def mono_mirror_offset_set(val,mir_num=None): @@ -677,7 +685,7 @@ def mono_mirror_offset_set(val,mir_num=None): Previously: Mono_Set_MIR0 """ - pvs=_mono_pvs(None,mir_num) + pvs=_mono_pvs_extended(None,mir_num) caput(pvs['mir_offset'],val) mono_get_all_extended(verbose=True) @@ -713,7 +721,7 @@ def mono_grating_mirror_offsets_set(mirror_offset): mono_get_all_extended(verbose=True) def mono_grating_b2_get(grt_num=None): - pvs = _mono_pvs(grt_num,None) + pvs = _mono_pvs_extended(grt_num,None) return caget(pvs['grt_b2']) def mono_grating_b2_set(val,grt_num=None): @@ -722,7 +730,7 @@ def mono_grating_b2_set(val,grt_num=None): Previously: Mono_Set_b2 """ hv = mono_energy_get() - pvs = _mono_pvs(grt_num,None) + pvs = _mono_pvs_extended(grt_num,None) caput(pvs['grt_b2t'],val) time.sleep(1) @@ -755,7 +763,7 @@ def mono_cff_set(val,tune_num): Previously: Mono_Set_cff """ hv = mono_energy_get() - pvs = _mono_pvs() + pvs = _mono_pvs_extended() #get current info current_vals = mono_get_all_extended(verbose=False) @@ -780,7 +788,7 @@ def mono_cff_set(val,tune_num): print('Differences : ',grt_dif,mir_dif) def mono_arm_get(): - pvs = _mono_pvs() + pvs = _mono_pvs_extended() return caget(pvs(['arm'])) def mono_arm_set(distance_mm): @@ -790,7 +798,7 @@ def mono_arm_set(distance_mm): Previously: Mono_Set_ExitArm """ hv = mono_energy_get() - pvs = _mono_pvs() + pvs = _mono_pvs_extended() #get current info current_vals = mono_get_all_extended(verbose=False) diff --git a/iexcode/instruments/cfg.py b/iexcode/instruments/cfg.py index e5767b1f89575b95bc1f02a1289cbd67d843febe..c44dc83474f83e1d93f3b602f0521ff1b134a056 100644 --- a/iexcode/instruments/cfg.py +++ b/iexcode/instruments/cfg.py @@ -5,4 +5,5 @@ BL=None mpa=None kappa_detectors_dictionary={} -EA=None \ No newline at end of file +EA=None +vortex = None \ No newline at end of file diff --git a/iexcode/macros/A_hutch/A_hutch.py b/iexcode/macros/A_hutch/A_hutch.py index f9b41d07d2ad1823f608fc249ba1f4e56b129406..90d49ec60bb150580e63509b51d21b6cd1485e30 100644 --- a/iexcode/macros/A_hutch/A_hutch.py +++ b/iexcode/macros/A_hutch/A_hutch.py @@ -199,7 +199,8 @@ def _Ahutch_detector_dictionary(): 12:"29idb:ca12:read", 13:"29idb:ca13:read", 14:"29idb:ca14:read", - 15:"29idb:ca15:read" + 15:"29idb:ca15:read", + 68:"29idb:ca14:read", } Kappa_scalers={ diff --git a/iexcode/macros/Kappa/__init__.py b/iexcode/macros/Kappa/__init__.py index 1580436f76c2e3aa81a033e81f5377d9677f1977..f3db549d2cd0c4cc3e7661080c0742401009a1ca 100644 --- a/iexcode/macros/Kappa/__init__.py +++ b/iexcode/macros/Kappa/__init__.py @@ -7,7 +7,7 @@ from iexcode.instruments.kappa_angle_calcs import * from iexcode.instruments.m3r_align import * from iexcode.instruments.m3r import * from iexcode.instruments.MPA import * -from iexcode.instruments.vortex import * +from iexcode.instruments.Vortex import * from iexcode.instruments.SRS_current_amplifiers import * from iexcode.macros.Kappa.Kappa_mca import * diff --git a/iexcode/macros/Kappa/hkl_iex.py b/iexcode/macros/Kappa/hkl_iex.py index 9e9038b2bcabaf8ba4bc86a7241b37fa674905e2..852befa63809b8b0859ebedf5731715e74ed0011 100644 --- a/iexcode/macros/Kappa/hkl_iex.py +++ b/iexcode/macros/Kappa/hkl_iex.py @@ -55,6 +55,7 @@ class FourCircleDiffractometer(SimulatedE4CV): # create the 4-circle diffractometer object fourc = FourCircleDiffractometer('', name="fourc") select_diffractometer(fourc) +az_ref = [0, 0, 1] # What the available modes with this diffractometer? print(f"\n{fourc.name} modes: {fourc.engine.modes}") @@ -381,7 +382,12 @@ def sample_new(name,A,B,C,Alpha,Beta,Gamma): lattice=Lattice( a=A, b=B, c=C, alpha=Alpha, beta=Beta, gamma=Gamma)) - print(mysample) + # print(mysample) + print(f'Sample = {mysample.name}') + print(f'Lattice = {mysample.lattice}') + print(f'U = {mysample.U}') + print(f'UB = {mysample.UB}') + print(f'Reflections: {mysample.reflections}') _UB_put('UBsample',fourc.calc.sample.name) _UB_put('UBlattice',fourc.calc.sample.lattice) except ValueError: @@ -846,7 +852,7 @@ def pa(diffractometer=fourc): print(f'\treal space = {a:.3f} {b:.3f} {c:.3f} / {alpha:.3f} {beta:.3f} {gamma:.3f}') print(f'\treciprocal space = {r_a:.3f} {r_b:.3f} {r_c:.3f} / {r_alpha:.3f} {r_beta:.3f} {r_gamma:.3f}\n') print(f'Azimuthal Reference:') - # print(f'\tH K L = {fourc.az_ref[0]} {fourc.az_ref[1]} {fourc.az_ref[2]}') + print(f'\tH K L = {az_ref[0]} {az_ref[1]} {az_ref[2]}') print(f'[U]: \n{np.round(sample.U, 3)}\n') print(f'[UB]: \n{np.round(sample.UB, 3)}') @@ -917,12 +923,60 @@ def set_az(g_haz, g_kaz, g_laz): ''' set azimutal orientation ''' - pass + az_ref = [g_haz, g_kaz, g_laz] + + +def angle_between(v1, v2): + angle = np.arccos(np.dot(v1, v2)/(np.linalg.norm(v1)*np.linalg.norm(v2))) + return np.rad2deg(angle) + + +def rot(v, axis, deg): + rad = np.deg2rad(deg) + if axis.lower()=='x': + r = np.array([[1, 0, 0],[0, np.cos(rad), -np.sin(rad)],[0, np.sin(rad), np.cos(rad)]]) + elif axis.lower()=='y': + r = np.array([[np.cos(rad), 0, np.sin(rad)],[0, 1, 0],[-np.sin(rad), 0, np.cos(rad)]]) + elif axis.lower()=='z': + r = np.array([[np.cos(rad), -np.sin(rad), 0],[np.sin(rad), np.cos(rad), 0],[0, 0, 1]]) + else: + print('axis should be either x, y, or z') + return np.dot(r, v) + + +def calc_psi_alpha_beta(az_ref=az_ref): + h = getattr(fourc, 'h').position + k = getattr(fourc, 'k').position + l = getattr(fourc, 'l').position + omega_now = omega_rbv.get() + phi_now = phi_rbv.get() + chi_now = chi_rbv.get() + tth_now = tth_rbv.get() + + incident = np.array((1, 0, 0)) + diffracted = np.array((np.cos(np.deg2rad(tth_now)), np.sin(np.deg2rad(tth_now)), 0)) + q_hkl = diffracted - incident + q_hkl /= np.linalg.norm(q_hkl) + scattering_plane = np.array((0, 0, 1)) + az_now = np.dot(fourc.calc.sample.UB, az_ref) + az_now = rot(az_now, 'y', phi_now) + az_now = rot(az_now, 'x', 90-chi_now) + az_now = rot(az_now, 'z', omega_now) + az_now /= np.linalg.norm(az_now) + az_projected = az_now - np.dot(az_now, q_hkl)*q_hkl + incident_projected = incident - np.dot(incident, q_hkl)*q_hkl + incident_projected /= np.linalg.norm(incident_projected) + az_projected_x = np.dot(az_projected, incident_projected) + az_projected_y = np.dot(az_projected, (0, 0, -1)) + psi = np.rad2deg(np.arctan2(az_projected_y, az_projected_x)) + alpha = angle_between(incident, az_now) - 90 + beta = 90-angle_between(az_now, diffracted) + return psi, alpha, beta def fixed_alpha(alpha): ''' - fix alpha mode, which is not avaiable in hklpy as of 202408. Need to use external variables + fix alpha mode, which is not avaiable in hklpy as of 202408. Need to use external variables az_ref SPEC modes: g_mode = diff --git a/iexcode/macros/Octupole/__init__.py b/iexcode/macros/Octupole/__init__.py index fa65535336ff20838febdbfc8eb13ba8d6a12417..3027ecbd4e41485b37991be793429771835eae3e 100644 --- a/iexcode/macros/Octupole/__init__.py +++ b/iexcode/macros/Octupole/__init__.py @@ -5,6 +5,6 @@ from iexcode.macros.Octupole_macros.Octupole_mca import * from iexcode.instruments.m3r_align import * from iexcode.instruments.m3r import * -from iexcode.instruments.vortex import * +from iexcode.instruments.Vortex import * from iexcode.instruments.SRS_current_amplifiers import * \ No newline at end of file diff --git a/iexcode/macros/resolution.py b/iexcode/macros/resolution.py index b1aaa8f9c2b88a22618cf33e01d7aa4942b5af94..851f93640cf317cf465bb9971aa4f4ea6de555e9 100644 --- a/iexcode/macros/resolution.py +++ b/iexcode/macros/resolution.py @@ -3,14 +3,19 @@ import numpy as np from epics import caget -#import iexcode.instruments.cfg as iex +import iexcode.instruments.cfg as iex from iexcode.instruments.IEX_BL_config import BLconfig_endstation_name from iexcode.instruments.utilities import take_closest_value from iexcode.macros.xrays.xrays import getE, slit_get from iexcode.instruments.VLS_PGM import mono_grating_get -from iexcode.instruments.electron_analyzer import EA, resolution_EA, SES_slit_get + +from iexcode.instruments.electron_analyzer import resolution_EA, SES_slit_get + + from iexcode.instruments.ARPES import ARPES_extra_pvs from iexcode.instruments.scanRecord import * + + ############################################################################################################# ############################## Resolution ############################## ############################################################################################################# @@ -26,7 +31,7 @@ def resolution(): if BLconfig_endstation_name() == "ARPES": slit_SES = SES_slit_get() - PE = int(EA.PassEnergy) + PE = int(iex.EA.PassEnergy) Tsample = caget(ARPES_extra_pvs()['TA']) resolution_ARPES(grt,hv_eV,slit_size,PE,slit_SES,Tsample,verbose=True) print("Calculate: resolution_ARPES(grt,hv_eV,slit_size,PE,slit_SES,Tsample)") @@ -97,6 +102,8 @@ def resolution_ARPES(grt,hv_eV,slit_size,PE,slit_SES,T,verbose=True): """ resolution of the ARPES """ + + bl = resolution_beamline(grt,hv_eV,slit_size,verbose=False) KbT = T*25/300 SES = resolution_EA(PE,slit_SES) diff --git a/iexcode/macros/staff/ID_energy_calibration.py b/iexcode/macros/staff/ID_energy_calibration.py index d2aaa738f1c893440a575a3a47b4ffd7298dca6e..57d5a24d3313e70a96ab17370d5ce109e75497bf 100644 --- a/iexcode/macros/staff/ID_energy_calibration.py +++ b/iexcode/macros/staff/ID_energy_calibration.py @@ -16,9 +16,9 @@ import iexcode.instruments.cfg as iex #from iexcode.macros.xrays.xrays import * from iexcode.macros.xrays.xrays import _energy_range_min_max, polarization, apertures_set from iexcode.instruments.AUX100_diode import current2flux -from iexcode.instruments.IEX_VPU import ID_QP_ratio_get, ID_set_eV +from iexcode.instruments.IEX_VPU import ID_QP_ratio_get, ID_set_eV, ID_pvs from iexcode.instruments.VLS_PGM import mono_grating_density_get, mono_grating_get, mono_energy_set, mono_scan - +from iexcode.instruments.VLS_PGM import mono_pvs from iexcode.instruments.diagnostics import diagnostics_all_out from iexcode.instruments.m3r import m3r_branch from iexcode.instruments.m3r_align import m3r_align @@ -26,8 +26,17 @@ from iexcode.instruments.m3r_align import m3r_align from iexcode.instruments.utilities import range_up, read_dict, today from iexcode.instruments.Logfile import log_name_set,log_print -from iexplot.mda_quick_plot import * +from iexplot.mda_quick_plot import * +from iexplot.utilities import make_num_list +from iexplot.pynData.nmda import nmda +from iexplot.IEX_funcs import IEX_scanNum_fpath + +#mono max and min energies +global plotMin,plotMax,ID_dict_path +plotMin = 100 +plotMax = 3000 +ID_dict_path = '/home/beams22/29IDUSER/Documents/User_Macros/Macros_29id/IEX_Dictionaries/Dict_IDCal.txt' ############################################################################################################## ################################ ID calibration scripts ############################## ############################################################################################################## @@ -128,276 +137,201 @@ def ID_calibration_scan(ID_start,ID_stop,ID_step,bandwidth=10,QP=None,Harm=1,sca ######################################### IDCalibration fitting ################################### ############################################################################################### +#======================= fitting and plotting ================================== +def calibration_metadata_pvs(): + d = { + 'id_sp':ID_pvs()['energy_sp'], + 'id_rbv':ID_pvs()['energy_rbv'], + 'id_mode':ID_pvs()['mode_rbv'], + 'mono_grt_den':mono_pvs()['grt_density'], + } + return d -def id_calibration_fit(first_scan,last_scan,det,PolyRank,**kwargs): - + + +def calibration_read_mda(*scans,detNum,**kwargs): + """ + *scans + **kwargs: - #id_calibration_fit(FirstLast[0],FirstLast[1],det,poly_order,mytitle=mytitle,countby=countby,plotMin=plt_min,plotMax=plt_max,plotType=plotType,filepath=filepath,prefix=prefix) """ - Fits calibration curves fpr each GRTs and modes included in the data set. Creates 3 dictionnaries: + kwargs.setdefault('ext','mda') + + modeDict={'CW, RCP':'RCP','CCW, LCP':'LCP','H':'H','V':'V'} + grtDict={1200:'MEG',2400:'HEG'} + + scanNum_list = make_num_list(*scans) - coefs={[GRT][mode][[breakpoint1,[coefs]],[breakpoint2,[coefs]...} - xdata= {[GRT][mode][[breakpoint1,[flux curve x axis]],[breakpoint2,[flux curve x axis]...} - fdata= {[GRT][mode][[breakpoint1,[flux curve y axis]],[breakpoint2,[flux curve y axis]...} + data_detNum = [] + data_mono_rbv = [] + id_sp_eV = [] + id_mode =[] + grt_name = [] + mono_max = [] + current_max = [] + flux_max =[] + mono_res = [] + current_res=[] + flux_res = [] + + for scanNum in scanNum_list: + #read the data + fpath = IEX_scanNum_fpath(scanNum,fullpath=True,**kwargs) + mda = nmda(fpath) + + x = mda.posx[0].data + data_mono_rbv.append(x) + y = mda.det[detNum].data + data_detNum.append(y) - kwargs: - countby = 1 by default - mytitle = '' by default - plotMin = min x range for plotting (default 250) - plotMax = max x range for plotting (default 3000) - plotType = ['dif,'fit,'flux'], full set by default - filepath = None, uses current mda filepath unless specified - e.g. user : filepath='/net/s29data/export/data_29idc/2018_2/UserName/mda/' - e.g. staff: filepath='/net/s29data/export/data_29idb/2018_2/mda/' - prefix = None, uses current mda prefix unless specified - scanIOC = None, uses BL_ioc() unless specified + #extra pvs + id_sp_eV.append(mda.header.all[meta_data_pvs()['id_sp']][2][0]*1000) + key = mda.header.all[meta_data_pvs()['id_mode']][2] + id_mode.append(modeDict[key]) + key = mda.header.all[meta_data_pvs()['mono_grt_den']][2][0] + grt_name.append(grtDict[key]) + #find max current + hv,current,flux = calibration_find_max(x,y) + mono_max.append(hv) + current_max.append(current) + flux_max.append(flux) + + return data_mono_rbv,data_detNum, id_sp_eV, id_mode,grt_name,mono_max, current_max, flux_max + + +def calibration_find_max(mono_rbv,intensity): """ + Finds the max intensity from a series of mono scans and + returns mono_energy_maxI,id_sp,current,flux - kwargs.setdefault('countby',1) - kwargs.setdefault('mytitle','') - kwargs.setdefault('plotMin',225) - kwargs.setdefault('plotMax',3000) - kwargs.setdefault('plotType',['fit','dif','flux']) - kwargs.setdefault('filepath',iex.BL.mda.filepath()) - kwargs.setdefault('prefix',iex.BL.mda.prefix()) - - countby=kwargs['countby'] - mytitle=kwargs['mytitle'] - plotMin=kwargs['plotMin'] - plotMax=kwargs['plotMax'] - filepath=kwargs['filepath'] - prefix=kwargs['prefix'] - + intensity = data_detNum + mono_rbv = data_mono_rbv - def calibration_curve(first,last,det,countby,filepath,prefix): - mono_list=[] - ID_list=[] - max_list=[] - flux_list=[] - print("filepath = ",filepath) - print("prefix = ",prefix) - print("First scan = ",first) - print("Last scan = ",last) - for i in range(first,last,countby): - x,y,x_name,y_name=mda_1D(i,det,1,0,filepath,prefix) #mda_1D(ScanNum,DetectorNum,coeff=1,bckg=0,filepath=None,prefix=None,scanIOC=None): - v,w,v_name,w_name=mda_1D(i, 4 ,1,0,filepath,prefix) - if y != []: - n=y.index(max(y)) # finds the peak max intensity index - e=round(x[n],2) # finds the corresponding mono energy - sp=round(w[2]*1000,0) # reads the ID set point - mono_list.append(e) - max_list.append(max(y)) - ID_list.append(sp) - flux_list.append(current2flux(max(y)),e,p=None) - return mono_list,ID_list,max_list,flux_list - - Mono_max,ID_SP,int_max,flux=calibration_curve(first_scan,last_scan+1,det,countby,filepath,prefix) - nPt_gaus=200 - x_HR = np.linspace(Mono_max[0], Mono_max[-1], nPt_gaus) + """ - # Data - xdata = np.array(Mono_max) - ydata = np.array(ID_SP) - zdata = np.array(int_max) - fdata = np.array(flux) - # Fitting - coefs = poly.polyfit(xdata, ydata, PolyRank) - ffit_HR = poly.polyval(x_HR, coefs) - ffit_Coarse = poly.polyval(xdata, coefs) - # Residual - Dif=np.array(ID_SP)-np.array(ffit_Coarse) - ratio=(np.array(ID_SP)-np.array(ffit_Coarse))/(np.array(ID_SP)/100) - # Plot differences - #print('plotMin = ',str(plotMin)) - if 'dif' in kwargs['plotType']: - fig = plt.figure(figsize=(12,6)) - plt.plot(Mono_max,ID_SP,marker='x',markersize=5,color='g',linewidth=0,label='data') - plt.plot(x_HR,ffit_HR,color='b',label='SP-fit') - plt.plot(xdata,Dif*100,marker='x',color='r',label='Dif x100') - plt.plot(xdata,ratio*1000,marker='x',color='g',label='Difx1000/ID') - plt.ylabel('ID SP') - plt.xlabel('Mono') - plt.xlim(plotMin,plotMax) - plt.legend(ncol=2, shadow=True, title=mytitle, fancybox=True) - plt.grid(linestyle='-', linewidth='0.5', color='grey') - plt.show() - # Plot raw data + fit - if 'fit' in kwargs['plotType']: - fig = plt.figure(figsize=(12,6)) - a1 = fig.add_axes([0,0,1,1]) - a1.plot(xdata+Dif,zdata,marker='*',markersize=10,color='r',linewidth=0,label='Interpolated ID SP') - for i in range(first_scan,last_scan+1): - x,y,x_name,y_name=mda_1D(i,det,1,0,filepath,prefix) - a1.plot(x,y,color='b') - a1.set_xlim(plotMin,plotMax) - a1.set(xlabel='Mono') - a1.set(ylabel='ID SP') - a1.ticklabel_format(style='sci', axis='y', scilimits=(0,0)) - a1.grid(linestyle='-', linewidth='0.5', color='grey') - plt.legend(ncol=2, shadow=True, title=mytitle, fancybox=True) - plt.show() - # Plot flux curves: - fdata_x=xdata+Dif - if 'flux' in kwargs['plotType']: - fig = plt.figure(figsize=(12,6)) - a1 = fig.add_axes([0,0,1,1]) - a1.plot(fdata_x,fdata,color='r',linewidth=1,label='Flux curves') - a1.set_xlim(plotMin,plotMax) - a1.set(xlabel='Mono') - a1.set(ylabel='ID SP') - a1.set_yscale('log') - #a1.ticklabel_format(style='sci', axis='y', scilimits=(0,0)) - a1.grid(linestyle='-', linewidth='0.5', color='grey') - plt.legend(ncol=2, shadow=True, title=mytitle, fancybox=True) - plt.show() - return coefs, fdata_x,fdata - - -def read_id_files(first,last,filepath=None,prefix=None,q=True): + # finds the peak max intensity index + intensity_max = np.max(intensity) + index_max = list(intensity).index(intensity_max) + mono_max = round(mono_rbv[index_max],2) + flux_max = current2flux(intensity_max,mono_max,verbose=False) + return mono_max,intensity_max,flux_max + +def Au_bkpt(id_sp_eV,Au=2100): """ - Reads extra PVs - Return a list of [[ScanNum,ID_SP,grt,mode],...] for mda files between first and last. + returns the index for the id setpoint closest to Au """ - if filepath == None: - filepath=iex.BL.mda.filepath() - if prefix == None: - prefix=iex.BL.mda.prefix()[:-1] - elif prefix[-1] == "_": - prefix=prefix[:-1] - mydata=mdaFile(first,last,filepath=filepath,prefix=prefix,q=q) - value=[] - modeDict={'CW, RCP':'RCP','CCW, LCP':'LCP','H':'H','V':'V'} - grtDict={1200:'MEG',2400:'HEG'} - for mdaNum in mydata.scanList: - if first<=mdaNum<=last: - extraPVs=mydata.header[mdaNum].all - try: - ID=round(extraPVs['ID29:EnergySet.VAL'][2][0]*1000,2) - mode=modeDict[extraPVs['ID29:ActualMode'][2]] # extraPVs return 'CW, RCP' - grt=grtDict[extraPVs['29idmono:GRT_DENSITY'][2][0]] # extraPVs return 1200 - except KeyError: - ID=round(extraPVs[b'ID29:EnergySet.VAL'][2][0]*1000,2) - mode=modeDict[str(extraPVs[b'ID29:ActualMode'][2])[2:-1]] - grt=grtDict[extraPVs[b'29idmono:GRT_DENSITY'][2][0]] # extraPVs return 1200 - if len(value)>0 and value[-1][1:] == [ID,grt,mode]: - pass - else: - value.append([mdaNum,ID,grt,mode]) - return value - - - -def id2num(ID,grt,mode,first=0,last=inf,ignore=[],filepath=None,prefix=None,q=True): # Not efficient - requires to read all 600 files everytime + Au=2100 + sp_Au = take_closest_value(id_sp_eV,Au) + n = id_sp_eV.index(sp_Au) + if sp_Au <= Au: + return n + else: + return n-1 + +def calibration_bkpts(id_sp_eV,**kwargs): """ - Return ScanNum corresponding to a given ID_SP from ExtraPVs (using mdaFile) """ - if filepath == None: - filepath=iex.BL.mda.filepath() - if prefix == None: - prefix=iex.BL.mda.prefix()[:-1] - elif prefix[-1] == "_": - prefix=prefix[:-1] - ScanNum = 0 - data_list = read_id_files(first,last,filepath,prefix,q) - data_short=[x for x in data_list if x[2:]==[grt,mode] and x[0] not in ignore] - step=data_short[1][1]-data_short[0][1] - precision=int(step/2) - ID1=ID-precision - ID2=ID+precision - ScanNum=[x[0] for x in data_short if ID1<= x[1]<= ID2] - if len(ScanNum)==0: result=None - else: result=ScanNum[0] - return result - - - -def num2id(ScanNum,grt,mode,first=0,last=inf,ignore=[],filepath=None,prefix=None,q=True): # Not efficient - requires to read all 600 files everytime + kwargs.setdefault('bkpts',[0,Au_bkpt(id_sp_eV)]) + n1,n2 = kwargs['bkpts'] + return n1,n2 + + +def calibration_do_fitting(mono_vals,id_sp_eV,**kwargs): """ - Return ID SP corresponding to a given ScanNum from ExtraPVs (using mdaFile) + fits the ID set point vs the mono values to a polynomial of rank PolyRand + returns id_sp_eV_fit,coefs + **kwargs + poly_rank (default: 1) """ - if filepath == None: - filepath=iex.BL.mda.filepath() - if prefix == None: - prefix=iex.BL.mda.prefix()[:-1] - elif prefix[-1] == "_": - prefix=prefix[:-1] - ID = 0 - data_short=[] - data_list = read_id_files(first,last,filepath,prefix,q) - data_short=[x for x in data_list if x[2:]==[grt,mode] and x[0] not in ignore] - ID=[x[1] for x in data_short if x[0] == ScanNum] - return ID[0] - - - -def extract_id(first,last,ignore=[],breakpts={'RCP':[600],'H':[400,600],'V':[600],'LCP':[600],'HEG':[],'MEG':[2200,2475]},filepath=None,prefix=None,q=True): + kwargs.setdefault('poly_rank',1) + + xdata = np.array(mono_vals) + ydata = np.array(id_sp_eV) + + # Fitting + n1,n2 = calibration_bkpts(id_sp_eV,**kwargs) + coefs = poly.polyfit(xdata[n1:n2], ydata[n1:n2], 4) + ft = poly.polyval(np.array(mono_vals), coefs) + + return ft,coefs + +def calibration_interp_fit(mono_vals,coefs,**kwargs): + """ + interpolates the mono_values with kwarg[npts] (default: 200) + and calculates id_setpoint based on the fit + returns mono_new, id_new + """ + kwargs.setdefault('npts',200) + mono_new = np.linspace(mono_vals[0], mono_vals[-1], kwargs['npts']) + id_new = poly.polyval(mono_new, coefs) + return mono_new, id_new + +def calibration_calc_delta(id_sp_eV, id_sp_eV_fit): + """ + returns the difference between the ID setpint and the interpolated value + """ + id_delta = np.array(id_sp_eV) - id_sp_eV_fit + return id_delta + +def calibration_plot_fit(hv_max,id_sp_eV,id_sp_eV_fit,**kwargs): """ - Breaksdown the info from a calibration files into grt & mode with the corresponding breakpoints (hardcoded): - [[first, last, 'HEG', 'RCP', [600]], - [first, last, 'HEG', 'H', [400, 600]]...] + plots the fit to t """ - if filepath == None: - filepath=iex.BL.mda.filepath() - if prefix == None: - prefix=iex.BL.mda.prefix()[:-1] - elif prefix[-1] == "_": - prefix=prefix[:-1] - IDlog=read_id_files(first,last,filepath,prefix,q) - #breakpts={'RCP':[600],'H':[400,600],'V':[600],'LCP':[600],'HEG':[],'MEG':[2200]} - #breakpts={'RCP':[600],'H':[400,600],'V':[600],'LCP':[600],'HEG':[],'MEG':[2200]} - #breakpts={'RCP':[600],'H':[600],'V':[600],'LCP':[600],'HEG':[],'MEG':[2200,2475]} # FR changed H breakpoints, missing low energy scans 06/14/2021 - data=[] - for grt in ['HEG','MEG']: - for mode in ['RCP','H','V','LCP']: - tmp=[x for x in IDlog if x[2:]==[grt,mode] and x[0] not in ignore] - #print(tmp) - if len(tmp)>0: - tmp=[tmp[0][0],tmp[-1][0],grt,mode,breakpts[mode]+breakpts[grt]] - data.append(tmp) - return data - - - -def update_id_dict(first,last,det,update_file,ignore,**kwargs): + fig = plt.figure(figsize=(12,6)) + #plotting data and difference between fit and data + plt.subplot(2,1,1) + plt.plot(hv_max, id_sp_eV,color='r',marker='o',label='data') + plt.plot(hv_max, id_sp_eV_fit,color='k',label='fit') + n1,n2 = calibration_bkpts(id_sp_eV,**kwargs) + plt.vlines(id_sp_eV[n2],plotMin,plotMax,color='g',label='breakpoint') + plt.legend() + plt.grid(linestyle='-', linewidth='0.5', color='grey') + + #difference plot + plt.subplot(2,1,2) + plt.plot(hv_max,id_delta,marker='.',label='difference') + plt.legend() + n1,n2 = calibration_bkpts(id_sp_eV) + ylims=np.min(id_delta[n1:n2]),np.max(id_delta[n1:n2]) + plt.ylim(ylims) + plt.grid(linestyle='-', linewidth='0.5', color='grey') + +def calibration_plot_data(data_mono_rbv,data_detNum,hv_max,current_max): + fig = plt.figure(figsize=(12,6)) + a1 = fig.add_axes([0,0,1,1]) + #plotting interpolated ID setpoint + a1.plot(np.array(hv_max)+id_delta,np.array(current_max),marker='*',markersize=10,color='r',linewidth=0,label='Interpolated ID SP') + #plotting data curves + for i,intensity in enumerate(data_detNum): + a1.plot(data_mono_rbv[i],data_detNum[i],color='b') + #plt.xlim(kwargs['plotMin'],['plotMax']) + a1.set(xlabel='Mono') + a1.set(ylabel='Intensity') + a1.ticklabel_format(style='sci', axis='y', scilimits=(0,0)) + a1.grid(linestyle='-', linewidth='0.5', color='grey') + plt.legend(ncol=2, shadow=True, title='', fancybox=True) + +def ID_calibration_analysis(*scans,detNum=15,**kwargs): """ Calculate new calibration curve for a full set of data.' But what if the set is not complete? To be addressed by Future Me If update_file == True (ONLY): - \tupdate the ID dictionary '/home/beams22/29IDUSER/Documents/User_Macros/Macros_29id/IEX_Dictionaries/Dict_IDCal.txt' + \tupdate the ID dictionary '' """ - - if update_file == True: - foo=input('\nAre you sure you want to update the ID calibration dictionary?\n>') - if foo == 'Y' or foo == 'y' or foo == 'yes'or foo == 'YES' or foo == 'Yes': - foo = 'yes' - print('\nWill save new dictionary to: \'/home/beams22/29IDUSER/Documents/User_Macros/Macros_29id/IEX_Dictionaries/Dict_IDCal.txt\'\n') - else: - print('\nCalibration curves will not be saved.\n') - else: - print('\nCalibration curves will not be saved.\n') - - kwargs.setdefault('countby',1) kwargs.setdefault('mytitle','') - kwargs.setdefault('plotMin',225) - kwargs.setdefault('plotMax',3000) + kwargs.setdefault('plotMin',plotMin) + kwargs.setdefault('plotMax',plotMax) kwargs.setdefault('plotType',['fit','dif','flux']) - kwargs.setdefault('filepath',iex.BL.mda.filepath()) - kwargs.setdefault('prefix',iex.BL.mda.prefix()) - kwargs.setdefault('q',True) - kwargs.setdefault('breakpts',{'RCP':[600],'H':[400,600],'V':[600],'LCP':[600],'HEG':[],'MEG':[2200,2475]}) + kwargs.setdefault('filepath',mda_filepath()) + kwargs.setdefault('prefix',mda_prefix()) + kwargs.setdefault('verbose',False) + kwargs.setdefault('breakpts',calibration_breakpts()) - - breakpts=kwargs['breakpts'] - plotType=kwargs['plotType'] - filepath=kwargs['filepath'] - prefix=kwargs['prefix'] - q=kwargs['q'] - - ID_data=extract_id(first,last,ignore,breakpts,filepath,prefix,q) + #extract_id returns data + id_sp, id_mode,grt_name,hv_max, current_max, flux_max = calibration_read_mda(*scans,detNum,**kwargs) ##### Extract parameters: - mode_dict ={'RCP':0,'LCP':1,'V':2,'H':3} id_coef={} new_id_function={} @@ -454,16 +388,15 @@ def update_id_dict(first,last,det,update_file,ignore,**kwargs): if grt=='MEG' and ID[1]>2900: poly_order=5 plt_min=ID[0]-200 - print('det =',det) + print('det =',detNum) print('poly_order =',poly_order) subkwargs=kwargs # kwargs contains plotType, filepath, prefix and q newkwargs={'mytitle':mytitle,'countby':countby,"plotMin":plt_min,'plotMax':plt_max} subkwargs.update(newkwargs) - #print(newkwargs) - #print(FirstLast[0],FirstLast[1]) - #result,energy_data,flux_data=id_calibration_fit(FirstLast[0],FirstLast[1],det,poly_order,subkwargs) - result,energy_data,flux_data=id_calibration_fit(FirstLast[0],FirstLast[1],det,poly_order,mytitle=mytitle,countby=countby,plotMin=plt_min,plotMax=plt_max,plotType=['fit','dif'],filepath=filepath,prefix=prefix) + + result,energy_data,flux_data=id_calibration_fit(FirstLast[0],FirstLast[1],detNum,poly_order,mytitle=mytitle,countby=countby,plotMin=plt_min,plotMax=plt_max,plotType=['fit','dif'],filepath=filepath,prefix=prefix) + tmp0.append([cutoff[1],result.tolist()]) tmp1.append([cutoff[1],energy_data]) tmp2.append([cutoff[1],flux_data]) @@ -477,37 +410,3 @@ def update_id_dict(first,last,det,update_file,ignore,**kwargs): flux_dict[grt]=id_flux ##### Read & update old dictionary: - - try: - id_function=read_dict('Dict_IDCal.txt') - print(id_function) - id_function.update(new_id_function) - - except KeyError: - print("Unable to read previous dictionary") - - - if update_file == True and foo == 'yes': - filepath = "/home/beams22/29IDUSER/Documents/User_Macros/Macros_29id/IEX_Dictionaries/" - filename = 'Dict_IDCal.txt' - - with open(join(filepath, filename), "a+") as f: - f.write('\n======= '+today()+': \n') - f.write(str(id_function)) - f.write('\n') - print('\nWriting dictionary to:',join(filepath, filename)) - - if update_file == True and foo == 'yes': - filepath = "/home/beams22/29IDUSER/Documents/User_Macros/Macros_29id/IEX_Dictionaries/" - filename = 'Flux_Curves.txt' - - with open(join(filepath, filename), "a+") as f: - f.write('\n======= '+today()+': \n') - f.write('\n----- flux_x:\n') - f.write(str(energy_dict)) - f.write('\n----- flux_y:\n') - f.write(str(flux_dict)) - f.write('\n') - print('\nWriting flux curves to:',join(filepath, filename)) - - return id_function,energy_dict,flux_dict