Update polyfit to interpolation

This commit is contained in:
linarphy 2023-07-20 17:10:16 +02:00
parent 00670cf47a
commit 920374b144
No known key found for this signature in database
GPG key ID: 3D4AAAC3AD16E79C

122
ETA.py
View file

@ -8,6 +8,7 @@ import shelve
from scipy.signal import convolve as sp_convolve from scipy.signal import convolve as sp_convolve
from scipy.signal import find_peaks from scipy.signal import find_peaks
from scipy.ndimage import rotate from scipy.ndimage import rotate
import scipy.interpolate as interpolate
cache , filename , output , calibration = None , None , None , None cache , filename , output , calibration = None , None , None , None
verbose , no_cache = False , False verbose , no_cache = False , False
@ -466,7 +467,7 @@ else:
if verbose: if verbose:
print( 'cache saved' ) print( 'cache saved' )
size_x = border[ 'x' ][ 'max' ] - border[ 'x' ][ 'min' ] size_x = border[ 'x' ][ 'max' ] - border[ 'x' ][ 'min' ]
size_y = calibrations[ 'down' ][ 'y' ][ 'min' ] - calibrations[ 'top' ][ 'y' ][ 'max' ] size_y = calibrations[ 'down' ][ 'y' ][ 'min' ] - calibrations[ 'top' ][ 'y' ][ 'max' ]
# Calibration # Calibration
@ -478,12 +479,16 @@ if calibration != None:
peaks_calib = np.loadtxt( calibration ) # sorted list peaks_calib = np.loadtxt( calibration ) # sorted list
peaks_calib = np.sort( peaks_calib ) peaks_calib = np.sort( peaks_calib )
mean_up = np.mean( data[ mean_up = np.mean( data[
calibrations[ 'top' ][ 'y' ][ 'min' ] : calibrations[ 'top' ][ 'y' ][ 'max' ], calibrations[ 'top' ][ 'y' ][ 'min' ]:
calibrations[ 'top' ][ 'x' ][ 'min' ] : calibrations[ 'top' ][ 'x' ][ 'max' ] calibrations[ 'top' ][ 'y' ][ 'max' ],
calibrations[ 'top' ][ 'x' ][ 'min' ]:
calibrations[ 'top' ][ 'x' ][ 'max' ]
] , axis = 0 ) ] , axis = 0 )
mean_down = np.mean( data[ mean_down = np.mean( data[
calibrations[ 'down' ][ 'y' ][ 'min' ] : calibrations[ 'down' ][ 'y' ][ 'max' ], calibrations[ 'down' ][ 'y' ][ 'min' ]:
calibrations[ 'down' ][ 'x' ][ 'min' ] : calibrations[ 'down' ][ 'x' ][ 'max' ] calibrations[ 'down' ][ 'y' ][ 'max' ],
calibrations[ 'down' ][ 'x' ][ 'min' ]:
calibrations[ 'down' ][ 'x' ][ 'max' ]
] , axis = 0 ) ] , axis = 0 )
peaks_up = np.array( peaks_up = np.array(
utils.retrieve_peaks( utils.retrieve_peaks(
@ -559,7 +564,9 @@ if calibration != None:
cons_diff , i = utils.same_value( diff ) , 0 cons_diff , i = utils.same_value( diff ) , 0
for list_same in cons_diff: for list_same in cons_diff:
new_data[ i : i + len( list_same ) , x ] = data[ new_data[ i : i + len( list_same ) , x ] = data[
calibrations[ 'top' ][ 'y' ][ 'max' ] + i : calibrations[ 'top' ][ 'y' ][ 'max' ] + i + len( list_same ), calibrations[ 'top' ][ 'y' ][ 'max' ] + i :
calibrations[ 'top' ][ 'y' ][ 'max' ] + i
+ len( list_same ),
calibrations[ 'top' ][ 'x' ][ 'min' ] + list_same[0] calibrations[ 'top' ][ 'x' ][ 'min' ] + list_same[0]
] ]
i += len( list_same ) i += len( list_same )
@ -572,7 +579,7 @@ if calibration != None:
print( 'y to x pixel calibrated' ) print( 'y to x pixel calibrated' )
polyval_wavelength = np.polyfit( # x = 0 begins at the start of polyval_wavelength = np.polyfit( # x = 0 begins at the start of
peaks_up , # calibrations[ 'top' ][ 'x' ][ 'min' ] peaks_up , # calibrations[ 'top' ][ 'x' ][ 'min' ]
peaks_calib, peaks_calib,
1 , 1 ,
) )
@ -600,38 +607,88 @@ if verbose:
""" """
remove artefact remove artefact
""" """
list_ = np.convolve( first_stripe = np.argmax(
np.gradient( np.convolve(
np.mean( np.gradient(
data[ np.mean(
calibrations[ 'top' ][ 'y' ][ 'max' ] : calibrations[ 'down' ][ 'y' ][ 'min' ], data[
border[ 'x' ][ 'min' ] : border[ 'x' ][ 'max' ] calibrations[ 'top' ][ 'y' ][ 'max' ] :
] , calibrations[ 'down' ][ 'y' ][ 'min' ],
axis = 1, border[ 'x' ][ 'min' ] : border[ 'x' ][ 'max' ]
), ] ,
) , axis = 1,
np.ones( 50 ), ),
'same' , ) ,
np.ones( 50 ),
'same' ,
),
) )
first_stripe = np.argmax( list_ ) range_ = np.arange(
range_ = range(
calibrations[ 'top' ][ 'y' ][ 'max' ] + first_stripe, calibrations[ 'top' ][ 'y' ][ 'max' ] + first_stripe,
calibrations[ 'down' ][ 'y' ][ 'max' ] calibrations[ 'down' ][ 'y' ][ 'max' ] ,
100 ,
) )
for line in range_: for line in range_:
abciss = range( convolve_point = 500
abciss = np.arange(
border[ 'x' ][ 'min' ], border[ 'x' ][ 'min' ],
border[ 'x' ][ 'max' ], border[ 'x' ][ 'max' ],
) )
list_ = data[ list_ = np.convolve(
line , data[
border[ 'x' ][ 'min' ] : border[ 'x' ][ 'max' ] line ,
border[ 'x' ][ 'min' ] : border[ 'x' ][ 'max' ]
] ,
np.ones(
convolve_point,
) ,
'valid', # Need to append value at extremities
) / convolve_point # keep correct value
list_ = np.concatenate(
(
data[
line ,
border[ 'x' ][ 'min' ] :
border[ 'x' ][ 'min' ] + convolve_point //
2 - 1
] ,
list_,
data[
line ,
border[ 'x' ][ 'max' ]
- convolve_point // 2 + convolve_point % 2 :
border[ 'x' ][ 'max' ]
] ,
),
)
number_point = 10
points_taken = abciss[
convolve_point // 2 + convolve_point % 2 :
- convolve_point // 2 + convolve_point % 2:
( abciss.shape[0] - convolve_point + 1 ) // ( number_point - 1 )
] ]
polyfit_ = np.polyfit( values_taken = list_[
abciss, convolve_point // 2 + convolve_point % 2 :
list_ , - convolve_point // 2 + convolve_point % 2:
8 , ( abciss.shape[0] - convolve_point + 1 ) // ( number_point - 1 )
]
sigma = np.std( list_ )
interpolation = interpolate.PchipInterpolator( # Seems the best
# (Have to check to be sure) -> TODO find reference
# value to determine which interpolation is better
points_taken,
values_taken,
)
wrong_points = np.where(
np.abs(
interpolation(
abciss,
) - list_ ,
) > sigma
)[0]
list_[ wrong_points ] = interpolation( # no border point
abciss[ wrong_points ],
) )
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
@ -641,9 +698,8 @@ for line in range_:
) )
plt.plot( plt.plot(
abciss, abciss,
np.polyval( interpolation(
polyfit_, abciss,
abciss ,
) , ) ,
) )
plt.show() plt.show()