From b06276664a1d807a237691be4e31078c917cc925 Mon Sep 17 00:00:00 2001 From: Carlos Gomez Date: Thu, 15 Dec 2016 15:19:30 +0100 Subject: [PATCH] VIP v0.6.2: Changes in code style (starting migration to Py3). Functions annular PCA and contrast curve avoid rounding the FWHM. --- vip/__init__.py | 2 +- vip/conf/mem.py | 2 +- vip/pca/pca_fullfr.py | 96 ++++++++++++++++++------------------- vip/pca/pca_local.py | 55 ++++++++++----------- vip/phot/contrcurve.py | 104 ++++++++++++++++++++-------------------- vip/stats/clip_sigma.py | 5 +- 6 files changed, 129 insertions(+), 135 deletions(-) diff --git a/vip/__init__.py b/vip/__init__.py index dc933a69..b5186331 100644 --- a/vip/__init__.py +++ b/vip/__init__.py @@ -12,7 +12,7 @@ from . import stats from . import var -__version__ = "0.6.1" +__version__ = "0.6.2" print("---------------------------------------------------") print(" oooooo oooo ooooo ooooooooo. ") diff --git a/vip/conf/mem.py b/vip/conf/mem.py index 50097ea4..49023422 100644 --- a/vip/conf/mem.py +++ b/vip/conf/mem.py @@ -30,7 +30,7 @@ def check_enough_memory(input_bytes, factor=1, verbose=True): ---------- input_bytes : float The size in bytes of the inputs of a given function. - factor : int + factor : float Scales how much memory is needed in terms of the size of input_bytes. """ mem = virtual_memory() diff --git a/vip/pca/pca_fullfr.py b/vip/pca/pca_fullfr.py index ccb277f5..227bd3a8 100644 --- a/vip/pca/pca_fullfr.py +++ b/vip/pca/pca_fullfr.py @@ -212,11 +212,11 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, if not cube.ndim>2: raise TypeError('Input array is not a 3d or 4d array') if angle_list is not None: - if (cube.ndim==3 and not (cube.shape[0] == angle_list.shape[0])) or \ + if (cube.ndim==3 and not (cube.shape[0] == angle_list.shape[0])) or \ (cube.ndim==4 and not (cube.shape[1] == angle_list.shape[0])): - msg = "Angle list vector has wrong length. It must equal the number" - msg += " frames in the cube." - raise TypeError(msg) + msg = "Angle list vector has wrong length. It must equal the number" + msg += " frames in the cube." + raise TypeError(msg) if source_xy is not None and delta_rot is None or fwhm is None: msg = 'Delta_rot or fwhm parameters missing. They are needed for the ' msg += 'PA-based rejection of frames from the library' @@ -268,7 +268,7 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, if ncomp > z: ncomp = min(10, z) msg = 'Number of PCs too high (max PCs={}), using instead {:} PCs.' - print msg.format(z, ncomp) + print(msg.format(z, ncomp)) scale_list = check_scal_vector(scale_list) #*********************************************************************** @@ -276,7 +276,7 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, #*********************************************************************** if cube.ndim==3: if verbose: - print '{:} spectral channels in IFS cube'.format(z) + print('{:} spectral channels in IFS cube'.format(z)) # cube has been re-scaled to have the planets moving radially cube, _, y, x, _, _ = scale_cube_for_pca(cube, scale_list) residuals_result = subtract_projection(cube, None, ncomp, scaling, @@ -299,7 +299,7 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, full_output=full_output, inverse=True, y_in=y_in, x_in=x_in) if verbose: - print 'Done re-scaling and combining' + print('Done re-scaling and combining') timing(start_time) #*********************************************************************** @@ -308,13 +308,13 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, #*********************************************************************** elif cube.ndim==4 and angle_list is not None: if verbose: - print '{:} spectral channels in IFS cube'.format(z) - print '{:} ADI frames in all channels'.format(n) + print('{:} spectral channels in IFS cube'.format(z)) + print('{:} ADI frames in all channels'.format(n)) residuals_cube_channels = np.zeros((n, y_in, x_in)) bar = pyprind.ProgBar(n, stream=1, title='Looping through ADI frames') - for i in xrange(n): + for i in range(n): cube_res, _, y, x, _, _ = scale_cube_for_pca(cube[:,i,:,:], scale_list) residuals_result = subtract_projection(cube_res, None, ncomp, @@ -341,15 +341,14 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, if ncomp2 > n: ncomp2 = min(10, n) msg = 'Number of PCs too high (max PCs={}), using instead {:} PCs.' - print msg.format(n, ncomp2) + print(msg.format(n, ncomp2)) elif ncomp2 is None: residuals_cube_channels_ = cube_derotate(residuals_cube_channels, angle_list) frame = cube_collapse(residuals_cube_channels_, mode=collapse) if verbose: - msg = 'De-rotating and combining' - print msg + print('De-rotating and combining') timing(start_time) else: @@ -360,9 +359,9 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, residuals_cube_channels_ = cube_derotate(res_ifs_adi, angle_list) frame = cube_collapse(residuals_cube_channels_, mode=collapse) if verbose: - msg = 'Done PCA per ADI multi-spectral frame, de-rotating and ' - msg += 'combining' - print msg + msg = 'Done PCA per ADI multi-spectral frame, de-rotating ' + msg += 'and combining' + print(msg) timing(start_time) #*************************************************************************** @@ -372,7 +371,7 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, if ncomp > n: ncomp = min(ncomp,n) msg = 'Number of PCs too high (max PCs={}), using instead {:} PCs.' - print msg.format(n, ncomp) + print(msg.format(n, ncomp)) residuals_result = subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, debug, svd_mode, verbose, full_output) @@ -388,7 +387,7 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, frame = cube_collapse(residuals_cube_, mode=collapse) if verbose: - print 'Done de-rotating and combining' + print('Done de-rotating and combining') timing(start_time) #*************************************************************************** @@ -398,7 +397,7 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, if ncomp > n: ncomp = min(ncomp,n) msg = 'Number of PCs too high (max PCs={}), using instead {:} PCs.' - print msg.format(n, ncomp) + print(msg.format(n, ncomp)) if source_xy is None: residuals_result = subtract_projection(cube, None, ncomp, scaling, @@ -425,10 +424,10 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, new_pa_th = float(mid_range - mid_range * 0.1) if verbose: msg = 'PA threshold {:.2f} is too big, will be set to {:.2f}' - print msg.format(pa_thr, new_pa_th) + print(msg.format(pa_thr, new_pa_th)) pa_thr = new_pa_th - for frame in xrange(n): + for frame in range(n): if ann_center > fwhm*3: # TODO: 3 optimal value? new parameter? ind = find_indices(angle_list, frame, pa_thr, True) else: @@ -457,7 +456,7 @@ def subtract_projection(cube, cube_ref, ncomp, scaling, mask_center_px, residuals_cube_ = cube_derotate(residuals_cube, angle_list) frame = cube_collapse(residuals_cube_, mode=collapse) if verbose: - print 'Done de-rotating and combining' + print('Done de-rotating and combining') timing(start_time) if full_output and cube.ndim<4: @@ -628,8 +627,8 @@ def grid(matrix, angle_list, y, x, mode, V, fwhm, fmerit, step, inti, intf, if full_output: frlist = [] counter = 0 if debug: - print 'Step current grid:', step - print 'PCs | SNR' + print('Step current grid:', step) + print('PCs | SNR') for pc in range(inti, intf+1, step): if full_output: snr, flux, frame = get_snr(matrix, angle_list, y, x, mode, V, @@ -645,22 +644,23 @@ def grid(matrix, angle_list, y, x, mode, V, fwhm, fmerit, step, inti, intf, if full_output: frlist.append(frame) nsteps += 1 if truncate and nsteps>2 and snr0: batch = hdulist[n].data[intfin:] msg = 'Processing batch [{},{}] with shape {}' if verbose: - print msg.format(intfin, n_frames, batch.shape) - print 'Batch size in memory = {:.3f} MB'.format(batch.nbytes/1e6) + print(msg.format(intfin, n_frames, batch.shape)) + print('Batch size in memory = {:.3f} MB'.format(batch.nbytes/1e6)) matrix = prepare_matrix(batch, verbose=False) ipca.partial_fit(matrix) @@ -953,10 +952,9 @@ def pca_incremental(cubepath, angle_list=None, n=0, batch_size=None, mean = ipca.mean_.reshape(batch.shape[1], batch.shape[2]) if verbose: - print - print 'Reconstructing and obtaining residuals' + print('\nReconstructing and obtaining residuals') medians = [] - for i in xrange(0, int(n_frames/batch_size)): + for i in range(0, int(n_frames/batch_size)): intini = i*batch_size intfin = (i+1)*batch_size batch = hdulist[n].data[intini:intfin] diff --git a/vip/pca/pca_local.py b/vip/pca/pca_local.py index e08cad53..c35ba7cc 100644 --- a/vip/pca/pca_local.py +++ b/vip/pca/pca_local.py @@ -4,7 +4,7 @@ Module with local smart pca (annulus-wise) serial and parallel implementations. """ -from __future__ import division +from __future__ import division, print_function __author__ = 'C. Gomez @ ULg' __all__ = ['pca_adi_annular', @@ -97,13 +97,13 @@ def define_annuli(angle_list, ann, n_annuli, fwhm, radius_int, annulus_width, if verbose: msg2 = 'Annulus {:}, Inn radius = {:.2f}, Ann center = {:.2f} ' - print msg2.format(int(ann+1),inner_radius, ann_center) + print(msg2.format(int(ann+1),inner_radius, ann_center)) return inner_radius, ann_center def fr_ref_correlation(vector, matrix): """ Getting the correlations """ lista = [] - for i in xrange(matrix.shape[0]): + for i in range(matrix.shape[0]): pears, _ = stats.pearsonr(vector, matrix[i]) lista.append(pears) @@ -134,13 +134,12 @@ def do_pca_annulus(ncomp, matrix, svd_mode, noise_error, data_ref): angle_list = check_PA_vector(angle_list) - fwhm = int(np.round(fwhm)) - annulus_width = int(np.round(asize * fwhm)) # equal size for all annuli - n_annuli = int(np.floor((y/2-radius_int)/annulus_width)) + annulus_width = asize * fwhm # equal size for all annuli + n_annuli = int(np.floor((y/2-radius_int)/annulus_width)) if verbose: - msg = '# annuli = {:}, Ann width = {:}, rounded FWHM = {:.3f}\n' - print msg.format(n_annuli, annulus_width, fwhm) - print 'PCA will be done locally per annulus and per quadrant.\n' + msg = '# annuli = {:}, Ann width = {:}, FWHM = {:.3f}\n' + print(msg.format(n_annuli, annulus_width, fwhm)) + print('PCA will be done locally per annulus and per quadrant.\n') cube_out = np.zeros_like(array) for ann in range(n_annuli): @@ -172,15 +171,15 @@ def do_pca_annulus(ncomp, matrix, svd_mode, noise_error, data_ref): cube_out[:, yy, xx] = residuals if verbose: - print '# frames in LIB = {}'.format(nfrslib) - print '# PCs = {}'.format(ncomps) - print 'Done PCA with {:} for current annulus'.format(svd_mode) + print('# frames in LIB = {}'.format(nfrslib)) + print('# PCs = {}'.format(ncomps)) + print('Done PCA with {:} for current annulus'.format(svd_mode)) timing(start_time) cube_der = cube_derotate(cube_out, angle_list) frame = cube_collapse(cube_der, mode=collapse) if verbose: - print 'Done derotating and combining.' + print('Done derotating and combining.') timing(start_time) if full_output: return cube_out, cube_der, frame @@ -290,13 +289,12 @@ def pca_adi_annular(cube, angle_list, radius_int=0, fwhm=4, asize=3, angle_list = check_PA_vector(angle_list) - fwhm = int(np.round(fwhm)) - annulus_width = int(np.round(asize * fwhm)) # equal size for all annuli - n_annuli = int(np.floor((y/2-radius_int)/annulus_width)) + annulus_width = asize * fwhm # equal size for all annuli + n_annuli = int(np.floor((y/2-radius_int)/annulus_width)) if verbose: - msg = '# annuli = {:}, Ann width = {:}, rounded FWHM = {:.3f}' - print msg.format(n_annuli, annulus_width, fwhm), '\n' - print 'PCA per annulus (and per quadrant if requested)\n' + msg = '# annuli = {:}, Ann width = {:}, FWHM = {:.3f}' + print(msg.format(n_annuli, annulus_width, fwhm), '\n') + print('PCA per annulus (and per quadrant if requested)\n') if nproc is None: # Hyper-threading "duplicates" the cores -> cpu_count/2 nproc = (cpu_count()/2) @@ -314,7 +312,7 @@ def pca_adi_annular(cube, angle_list, radius_int=0, fwhm=4, asize=3, else: msge = 'If ncomp is a list, it must match the number of annuli' msg = 'NCOMP : {}, N ANN : {}, ANN SIZE : {}, ANN WIDTH : {}' - print msg.format(ncomp, n_annuli, annulus_width, asize) + print(msg.format(ncomp, n_annuli, annulus_width, asize)) raise TypeError(msge) else: ncompann = ncomp @@ -378,7 +376,7 @@ def pca_adi_annular(cube, angle_list, radius_int=0, fwhm=4, asize=3, cube_out[frame][yy, xx] = residuals[frame] if verbose: - print 'Done PCA with {:} for current annulus'.format(svd_mode) + print('Done PCA with {:} for current annulus'.format(svd_mode)) timing(start_time) #*************************************************************************** @@ -387,7 +385,7 @@ def pca_adi_annular(cube, angle_list, radius_int=0, fwhm=4, asize=3, cube_der = cube_derotate(cube_out, angle_list) frame = cube_collapse(cube_der, mode=collapse) if verbose: - print 'Done derotating and combining.' + print('Done derotating and combining.') timing(start_time) if full_output: return cube_out, cube_der, frame @@ -426,12 +424,12 @@ def define_annuli(angle_list, ann, n_annuli, fwhm, radius_int, annulus_width, new_pa_th = float(mid_range - mid_range * 0.1) if verbose: msg = 'PA threshold {:.2f} is too big, will be set to {:.2f}' - print msg.format(pa_threshold, new_pa_th) + print(msg.format(pa_threshold, new_pa_th)) pa_threshold = new_pa_th if verbose: msg2 = 'Annulus {:}, PA thresh = {:.2f}, Inn radius = {:.2f}, Ann center = {:.2f} ' - print msg2.format(int(ann+1),pa_threshold,inner_radius, ann_center) + print(msg2.format(int(ann+1),pa_threshold,inner_radius, ann_center)) return pa_threshold, inner_radius, ann_center @@ -443,13 +441,13 @@ def find_indices(angle_list, frame, thr, truncate): n = angle_list.shape[0] index_prev = 0 index_foll = frame - for i in xrange(0, frame): + for i in range(0, frame): if np.abs(angle_list[frame]-angle_list[i]) < thr: index_prev = i break else: index_prev += 1 - for k in xrange(frame, n): + for k in range(frame, n): if np.abs(angle_list[k]-angle_list[frame]) > thr: index_foll = k break @@ -489,7 +487,7 @@ def do_pca_loop(matrix, yy, xx, nproc, angle_list, fwhm, pa_threshold, scaling, nfrslib = [] if nproc==1: residualarr = [] - for frame in xrange(n): + for frame in range(n): res = do_pca_patch(matrix_ann, frame, angle_list, fwhm, pa_threshold, scaling, ann_center, svd_mode, ncomp, min_frames_pca, tol, @@ -560,7 +558,6 @@ def do_pca_patch(matrix, frame, angle_list, fwhm, pa_threshold, scaling, data = data_ref #data = data_ref - data_ref.mean(axis=0) - curr_frame = matrix[frame] # current frame V = get_eigenvectors(ncomp, data, svd_mode, noise_error=tol, debug=False) @@ -604,7 +601,7 @@ def get_eigenvectors(ncomp, data, svd_mode, noise_error=1e-3, max_evs=200, if ncomp>1: px_noise_decay = px_noise[-2] - px_noise[-1] #print 'ncomp {:} {:.4f} {:.4f}'.format(ncomp,px_noise[-1],px_noise_decay) - if debug: print 'ncomp', ncomp + if debug: print('ncomp', ncomp) else: # Performing SVD/PCA according to "svd_mode" flag diff --git a/vip/phot/contrcurve.py b/vip/phot/contrcurve.py index b4bc6b73..ac81ead6 100644 --- a/vip/phot/contrcurve.py +++ b/vip/phot/contrcurve.py @@ -4,8 +4,7 @@ Module with contrast curve generation function. """ -from __future__ import division -from __future__ import print_function +from __future__ import division, print_function __author__ = 'C. Gomez, O. Absil @ ULg' __all__ = ['contrast_curve', @@ -16,14 +15,14 @@ import numpy as np import photutils import inspect -from scipy.interpolate import interp1d, InterpolatedUnivariateSpline +from scipy.interpolate import InterpolatedUnivariateSpline from scipy import stats from scipy.signal import savgol_filter from skimage.draw import circle from matplotlib import pyplot as plt from pandas import DataFrame as DF from .fakecomp import inject_fcs_cube, inject_fc_frame, psf_norm -from ..conf import timeInit, timing +from ..conf import timeInit, timing, sep from ..var import frame_center, dist @@ -113,13 +112,20 @@ def contrast_curve(cube, angle_list, psf_template, fwhm, pxscale, starphot, if not starphot.shape[0] == cube.shape[0]: raise TypeError('Correction vector has bad size') cube = cube.copy() - for i in xrange(cube.shape[0]): + for i in range(cube.shape[0]): cube[i] = cube[i] / starphot[i] - fwhm = int(np.round(fwhm)) + if verbose: + start_time = timeInit() + if isinstance(starphot, float) or isinstance(starphot, int): + msg0 = 'ALGO : {}, FWHM = {}, # BRANCHES = {}, SIGMA = {},' + msg0 += ' STARPHOT = {}' + print(msg0.format(algo.func_name, fwhm, nbranch, sigma, starphot)) + else: + msg0 = 'ALGO : {}, FWHM = {}, # BRANCHES = {}, SIGMA = {}' + print(msg0.format(algo.func_name, fwhm, nbranch, sigma)) + print(sep) - if verbose: start_time = timeInit() - # throughput verbose_thru = False if verbose==2: verbose_thru = True @@ -133,8 +139,7 @@ def contrast_curve(cube, angle_list, psf_template, fwhm, pxscale, starphot, frame_nofc = res_throug[5] if verbose: - msg1 = 'Finished the throughput calculation' - print(msg1) + print('Finished the throughput calculation') timing(start_time) if thruput_mean[-1]==0: @@ -143,7 +148,8 @@ def contrast_curve(cube, angle_list, psf_template, fwhm, pxscale, starphot, # noise measured in the empty PP-frame with better sampling, every px # starting from 1*FWHM - noise_samp, rad_samp = noise_per_annulus(frame_nofc, 1, fwhm, wedge, False) + noise_samp, rad_samp = noise_per_annulus(frame_nofc, separation=1, fwhm=fwhm, + init_rad=fwhm, wedge=wedge) cutin1 = np.where(rad_samp.astype(int)==vector_radd.astype(int).min())[0] noise_samp = noise_samp[cutin1:] rad_samp = rad_samp[cutin1:] @@ -170,10 +176,6 @@ def contrast_curve(cube, angle_list, psf_template, fwhm, pxscale, starphot, window_length=win) if debug: - print('SIGMA = {}'.format(sigma)) - if isinstance(starphot, float) or isinstance(starphot, int): - print('STARPHOT = {}'.format(starphot)) - plt.rc("savefig", dpi=dpi) plt.figure(figsize=(8,4)) plt.plot(vector_radd*pxscale, thruput_mean, '.', label='computed', @@ -365,7 +367,6 @@ def throughput(cube, angle_list, psf_template, fwhm, pxscale, algo, nbranch=1, """ array = cube parangles = angle_list - fwhm = int(np.round(fwhm)) if not array.ndim == 3: raise TypeError('The input array is not a cube') @@ -403,7 +404,7 @@ def throughput(cube, angle_list, psf_template, fwhm, pxscale, algo, nbranch=1, timing(start_time) noise, vector_radd = noise_per_annulus(frame_nofc, separation=fwhm, - fwhm=fwhm, wedge=wedge,verbose=False) + fwhm=fwhm, wedge=wedge) vector_radd = vector_radd[inner_rad-1:] noise = noise[inner_rad-1:] if verbose: @@ -494,8 +495,8 @@ def throughput(cube, angle_list, psf_template, fwhm, pxscale, algo, nbranch=1, -def noise_per_annulus(array, separation, fwhm, wedge=(0,360), verbose=False, - debug=False): +def noise_per_annulus(array, separation, fwhm, init_rad=None, wedge=(0,360), + verbose=False, debug=False): """ Measures the noise as the standard deviation of apertures defined in each annulus with a given separation. @@ -508,6 +509,8 @@ def noise_per_annulus(array, separation, fwhm, wedge=(0,360), verbose=False, center of the frame. fwhm : float FWHM in pixels. + init_rad : float + Initial radial distance to be used. If None then the init_rad = FWHM. wedge : tuple of floats, optional Initial and Final angles for using a wedge. For example (-90,90) only considers the right side of an image. Be careful when using small @@ -530,8 +533,6 @@ def find_coords(rad, sep, init_angle, fin_angle): angular_range = fin_angle-init_angle npoints = (np.deg2rad(angular_range)*rad)/sep #(2*np.pi*rad)/sep ang_step = angular_range/npoints #360/npoints - #print (2*np.pi*rad)/sep, 360/((2*np.pi*rad)/sep) - #print angular_range, npoints, ang_step x = [] y = [] for i in range(int(npoints)): @@ -556,44 +557,41 @@ def find_coords(rad, sep, init_angle, fin_angle): noise = [] vector_radd = [] if verbose: print('{} annuli'.format(n_annuli-1)) - + + if init_rad is None: init_rad = fwhm + if debug: _, ax = plt.subplots(figsize=(6,6)) ax.imshow(array, origin='lower', interpolation='nearest', alpha=0.5, cmap='gray') - for _ in range(n_annuli-1): - y -= separation + for i in range(n_annuli-1): + y = centery + init_rad + separation*(i) rad = dist(centery, centerx, y, x) - if rad>=fwhm: - yy, xx = find_coords(rad, sep=fwhm, init_angle=init_angle, - fin_angle=fin_angle) - yy += centery - xx += centerx - - #print yy, xx - - fluxes = [] - apertures = photutils.CircularAperture((xx, yy), fwhm/2.) - fluxes = photutils.aperture_photometry(array, apertures) - fluxes = np.array(fluxes['aperture_sum']) - - noise_ann = np.std(fluxes) - noise.append(noise_ann) - vector_radd.append(rad) - - if debug: - for i in range(xx.shape[0]): - # Circle takes coordinates as (X,Y) - aper = plt.Circle((xx[i], yy[i]), radius=fwhm/2., color='r', - fill=False, alpha=0.8) - ax.add_patch(aper) - cent = plt.Circle((xx[i], yy[i]), radius=0.8, color='r', - fill=True, alpha=0.5) - ax.add_patch(cent) - - if verbose: - print('Radius(px) = {:}, Noise = {:.3f} '.format(rad, noise_ann)) + yy, xx = find_coords(rad, fwhm, init_angle, fin_angle) + yy += centery + xx += centerx + + apertures = photutils.CircularAperture((xx, yy), fwhm/2.) + fluxes = photutils.aperture_photometry(array, apertures) + fluxes = np.array(fluxes['aperture_sum']) + + noise_ann = np.std(fluxes) + noise.append(noise_ann) + vector_radd.append(rad) + + if debug: + for i in range(xx.shape[0]): + # Circle takes coordinates as (X,Y) + aper = plt.Circle((xx[i], yy[i]), radius=fwhm/2., color='r', + fill=False, alpha=0.8) + ax.add_patch(aper) + cent = plt.Circle((xx[i], yy[i]), radius=0.8, color='r', + fill=True, alpha=0.5) + ax.add_patch(cent) + + if verbose: + print('Radius(px) = {:}, Noise = {:.3f} '.format(rad, noise_ann)) return np.array(noise), np.array(vector_radd) diff --git a/vip/stats/clip_sigma.py b/vip/stats/clip_sigma.py index 1fc1c766..08575e26 100644 --- a/vip/stats/clip_sigma.py +++ b/vip/stats/clip_sigma.py @@ -4,13 +4,14 @@ Module with sigma clipping functions. """ +from __future__ import division, print_function + __author__ = 'C. Gomez @ ULg', 'V. Christiaens' __all__ = ['clip_array', 'sigma_filter'] import numpy as np from scipy.ndimage.filters import generic_filter -from scipy.stats import norm as Gaussian from astropy.stats import median_absolute_deviation @@ -80,7 +81,7 @@ def sigma_filter(frame_tmp, bpix_map, neighbor_box=3, min_neighbors=3, verbose=F bp[wb[0][n],wb[1][n]] = 0 nb = int(np.sum(bp)) if verbose == True: - print 'Required number of iterations in the sigma filter: ', nit + print('Required number of iterations in the sigma filter: ', nit) return im