Tibeuleu/FOC_Reduction
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

remove border effect from weighted gaussian smoothing

Browse Source
This commit is contained in:
Tibeuleu
2023-01-31 17:28:20 +01:00
parent 8a0bcef84f
commit 5c0d08a44b
1 changed files with 10 additions and 11 deletions
Download Patch File Download Diff File Expand all files Collapse all files

21
src/lib/reduction.py
View File

@@ -849,30 +849,29 @@ def smooth_data(data_array, error_array, data_mask, headers, FWHM=1.,
error[r,c] = np.where(data_mask[r,c], np.sqrt(np.sum(weight*g_rc**2))/np.sum(weight*g_rc), (np.sqrt(np.sum(error_array**2,axis=0)/error_array.shape[0]))[r,c]) error[r,c] = np.where(data_mask[r,c], np.sqrt(np.sum(weight*g_rc**2))/np.sum(weight*g_rc), (np.sqrt(np.sum(error_array**2,axis=0)/error_array.shape[0]))[r,c])
# Nan handling # Nan handling
error[np.isnan(smoothed)] = 0. error[np.logical_or(np.isnan(smoothed*error),1-data_mask)] = 0.
smoothed[np.isnan(smoothed)] = 0. smoothed[np.logical_or(np.isnan(smoothed*error),1-data_mask)] = 0.
error[np.isnan(error)] = 0.
elif smoothing.lower() in ['weight_gauss','weighted_gaussian','gauss','gaussian']: elif smoothing.lower() in ['weight_gauss','weighted_gaussian','gauss','gaussian']:
# Convolution with gaussian function # Convolution with gaussian function
smoothed = np.zeros(data_array.shape) smoothed = np.zeros(data_array.shape)
error = np.zeros(error_array.shape) error = np.zeros(error_array.shape)
for i,image in enumerate(data_array): for i,(image,image_error) in enumerate(zip(data_array, error_array)):
x, y = np.meshgrid(np.arange(-image.shape[1]/2,image.shape[1]/2),np.arange(-image.shape[0]/2,image.shape[0]/2)) x, y = np.meshgrid(np.arange(-image.shape[1]/2,image.shape[1]/2),np.arange(-image.shape[0]/2,image.shape[0]/2))
weights = np.ones(error_array[i].shape) weights = np.ones(image_error.shape)
if smoothing.lower()[:6] in ['weight']: if smoothing.lower()[:6] in ['weight']:
weights = 1./error_array[i]**2 weights = 1./image_error**2
weights[1-data_mask] = 0. weights[(1-np.isfinite(weights)).astype(bool)] = 0.
weights[(1-data_mask).astype(bool)] = 0.
weights /= weights.sum() weights /= weights.sum()
kernel = gaussian2d(x, y, stdev) kernel = gaussian2d(x, y, stdev)
kernel /= kernel.sum() kernel /= kernel.sum()
smoothed[i] = np.where(data_mask, convolve2d(image*weights,kernel,'same')/convolve2d(weights,kernel,'same'), image) smoothed[i] = np.where(data_mask, convolve2d(image*weights,kernel,'same')/convolve2d(weights,kernel,'same'), image)
error[i] = np.where(data_mask, np.sqrt(convolve2d(error_array[i]**2*weights**2,kernel**2,'same'))/convolve2d(weights,kernel,'same'), error_array[i]) error[i] = np.where(data_mask, np.sqrt(convolve2d(image_error**2*weights**2,kernel**2,'same'))/convolve2d(weights,kernel,'same'), image_error)
# Nan handling # Nan handling
error[i][np.isnan(smoothed[i])] = 0. error[i][np.logical_or(np.isnan(smoothed[i]*error[i]),1-data_mask)] = 0.
smoothed[i][np.isnan(smoothed[i])] = 0. smoothed[i][np.logical_or(np.isnan(smoothed[i]*error[i]),1-data_mask)] = 0.
error[i][np.isnan(error[i])] = 0.
else: else:
raise ValueError("{} is not a valid smoothing option".format(smoothing)) raise ValueError("{} is not a valid smoothing option".format(smoothing))