我们从Python开源项目中,提取了以下26个代码示例,用于说明如何使用matplotlib.pylab.axis()。
def plot_position(self, pos_true, pos_est): N = pos_est.shape[1] pos_true = pos_true[:, :N] pos_est = pos_est[:, :N] # Figure plt.figure() plt.suptitle("Position") # Ground truth plt.plot(pos_true[0, :], pos_true[1, :], color="red", marker="o", label="Grouth truth") # Estimated plt.plot(pos_est[0, :], pos_est[1, :], color="blue", marker="o", label="Estimated") # Plot labels and legends plt.xlabel("East (m)") plt.ylabel("North (m)") plt.axis("equal") plt.legend(loc=0)
def show_samples(y, ndim, nb=10, cmap=''): if ndim == 4: for i in range(nb**2): plt.subplot(nb, nb, i+1) plt.imshow(y[i], cmap=cmap, interpolation='none') plt.axis('off') else: x = y[0] y = y[1] plt.figure(0) for i in range(10): plt.subplot(2, 5, i+1) plt.imshow(x[i], cmap=cmap, interpolation='none') plt.axis('off') plt.figure(1) for i in range(10): plt.subplot(2, 5, i+1) plt.imshow(y[i], cmap=cmap, interpolation='none') plt.axis('off') plt.show()
def linear_testing(): x_axis = np.linspace(1, 51, 100) x_nice = np.linspace(x_axis[0], x_axis[-1], 100) mod, params = qudi_fitting.make_linear_model() print('Parameters of the model', mod.param_names, ' with the independet variable', mod.independent_vars) params['slope'].value = 2 # + abs(np.random.normal(0,1)) params['offset'].value = 50 #+ abs(np.random.normal(0, 200)) #print('\n', 'beta', params['beta'].value, '\n', 'lifetime', #params['lifetime'].value) data_noisy = (mod.eval(x=x_axis, params=params) + 10 * np.random.normal(size=x_axis.shape)) result = qudi_fitting.make_linear_fit(axis=x_axis, data=data_noisy, add_parameters=None) plt.plot(x_axis, data_noisy, 'ob') plt.plot(x_nice, mod.eval(x=x_nice, params=params), '-g') print(result.fit_report()) plt.plot(x_axis, result.best_fit, '-r', linewidth=2.0) plt.plot(x_axis, result.init_fit, '-y', linewidth=2.0) plt.show()
def fit_data(): data=np.loadtxt('data.dat') print(data) params = dict() params["c"] = {"min" : -np.inf,"max" : np.inf} result = qudi_fitting.make_lorentzian_fit(axis=data[:,0], data=data[:,3], add_parameters=params) print(result.fit_report()) plt.plot(data[:,0],-data[:,3]+2,"b-o",label="data mean") # plt.plot(data[:,0],data[:,1],label="data") # plt.plot(data[:,0],data[:,2],label="data") plt.plot(data[:,0],-result.best_fit+2,"r-",linewidth=2.,label="fit") # plt.plot(data[:,0],result.init_fit,label="init") plt.xlabel("time (ns)") plt.ylabel("polarization transfer (arb. u.)") plt.legend(loc=1) # plt.savefig("pol20_24repetition_pol.pdf") # plt.savefig("pol20_24repetition_pol.png") plt.show() savedata=[[data[ii,0],-data[ii,3]+2,-result.best_fit[ii]+2] for ii in range(len(data[:,0]))] np.savetxt("pol_data_fit.csv",savedata) # print(result.params) print(result.params)
def plotSequenceForRotatingState(degPerStep, Sigma, T=1000): theta = degPerStep * np.pi / 180. # radPerStep A = np.asarray([ [np.cos(theta), -np.sin(theta)], [np.sin(theta), np.cos(theta)], ]) Sigma = np.asarray(Sigma) if Sigma.size == 1: Sigma = Sigma * np.eye(2) elif Sigma.size == 2: Sigma = np.diag(Sigma) X = np.zeros((T, 2)) X[0, :] = [1, 0] for t in xrange(1, T): X[t] = np.random.multivariate_normal(np.dot(A, X[t - 1]), Sigma) pylab.plot(X[:, 0], X[:, 1], '.') pylab.axis([-4, 4, -4, 4]) pylab.axis('equal')
def showExampleDocs(pylab=None, nrows=3, ncols=3): if pylab is None: from matplotlib import pylab Data = get_data(seed=0, nObsPerDoc=200) PRNG = np.random.RandomState(0) chosenDocs = PRNG.choice(Data.nDoc, nrows * ncols, replace=False) for ii, d in enumerate(chosenDocs): start = Data.doc_range[d] stop = Data.doc_range[d + 1] Xd = Data.X[start:stop] pylab.subplot(nrows, ncols, ii + 1) pylab.plot(Xd[:, 0], Xd[:, 1], 'k.') pylab.axis('image') pylab.xlim([-1.5, 1.5]) pylab.ylim([-1.5, 1.5]) pylab.xticks([]) pylab.yticks([]) pylab.tight_layout() # Set Toy Parameters ###########################################################
def plot1D_mat(a, b, M, title=''): """ Plot matrix M with the source and target 1D distribution Creates a subplot with the source distribution a on the left and target distribution b on the tot. The matrix M is shown in between. Parameters ---------- a : np.array, shape (na,) Source distribution b : np.array, shape (nb,) Target distribution M : np.array, shape (na,nb) Matrix to plot """ na, nb = M.shape gs = gridspec.GridSpec(3, 3) xa = np.arange(na) xb = np.arange(nb) ax1 = pl.subplot(gs[0, 1:]) pl.plot(xb, b, 'r', label='Target distribution') pl.yticks(()) pl.title(title) ax2 = pl.subplot(gs[1:, 0]) pl.plot(a, xa, 'b', label='Source distribution') pl.gca().invert_xaxis() pl.gca().invert_yaxis() pl.xticks(()) pl.subplot(gs[1:, 1:], sharex=ax1, sharey=ax2) pl.imshow(M, interpolation='nearest') pl.axis('off') pl.xlim((0, nb)) pl.tight_layout() pl.subplots_adjust(wspace=0., hspace=0.2)
def plot_position(self, pos_true, pos_est, cam_states): N = pos_est.shape[1] pos_true = pos_true[:, :N] pos_est = pos_est[:, :N] # Figure plt.figure() plt.suptitle("Position") # Ground truth plt.plot(pos_true[0, :], pos_true[1, :], color="red", label="Grouth truth") # color="red", marker="x", label="Grouth truth") # Estimated plt.plot(pos_est[0, :], pos_est[1, :], color="blue", label="Estimated") # color="blue", marker="o", label="Estimated") # Sliding window cam_pos = [] for cam_state in cam_states: cam_pos.append(cam_state.p_G) cam_pos = np.array(cam_pos).reshape((len(cam_pos), 3)).T plt.plot(cam_pos[0, :], cam_pos[1, :], color="green", label="Camera Poses") # color="green", marker="o", label="Camera Poses") # Plot labels and legends plt.xlabel("East (m)") plt.ylabel("North (m)") plt.axis("equal") plt.legend(loc=0)
def plot_velocity(self, timestamps, vel_true, vel_est): N = vel_est.shape[1] t = timestamps[:N] vel_true = vel_true[:, :N] vel_est = vel_est[:, :N] # Figure plt.figure() plt.suptitle("Velocity") # X axis plt.subplot(311) plt.plot(t, vel_true[0, :], color="red", label="Ground_truth") plt.plot(t, vel_est[0, :], color="blue", label="Estimate") plt.title("x-axis") plt.xlabel("Date Time") plt.ylabel("ms^-1") plt.legend(loc=0) # Y axis plt.subplot(312) plt.plot(t, vel_true[1, :], color="red", label="Ground_truth") plt.plot(t, vel_est[1, :], color="blue", label="Estimate") plt.title("y-axis") plt.xlabel("Date Time") plt.ylabel("ms^-1") plt.legend(loc=0) # Z axis plt.subplot(313) plt.plot(t, vel_true[2, :], color="red", label="Ground_truth") plt.plot(t, vel_est[2, :], color="blue", label="Estimate") plt.title("z-axis") plt.xlabel("Date Time") plt.ylabel("ms^-1") plt.legend(loc=0)
def plot_attitude(self, timestamps, att_true, att_est): # Setup N = att_est.shape[1] t = timestamps[:N] att_true = att_true[:, :N] att_est = att_est[:, :N] # Figure plt.figure() plt.suptitle("Attitude") # X axis plt.subplot(311) plt.plot(t, att_true[0, :], color="red", label="Ground_truth") plt.plot(t, att_est[0, :], color="blue", label="Estimate") plt.title("x-axis") plt.legend(loc=0) plt.xlabel("Date Time") plt.ylabel("rad s^-1") # Y axis plt.subplot(312) plt.plot(t, att_true[1, :], color="red", label="Ground_truth") plt.plot(t, att_est[1, :], color="blue", label="Estimate") plt.title("y-axis") plt.legend(loc=0) plt.xlabel("Date Time") plt.ylabel("rad s^-1") # Z axis plt.subplot(313) plt.plot(t, att_true[2, :], color="red", label="Ground_truth") plt.plot(t, att_est[2, :], color="blue", label="Estimate") plt.title("z-axis") plt.legend(loc=0) plt.xlabel("Date Time") plt.ylabel("rad s^-1")
def axis(self, imin, imax, jmin, jmax, xlabels=[], ylabels=[]): imid=(imin+imax)/2 # midpoint of X-axis jmid=(jmin+jmax)/2 # midpoint of Y-axis # Create axis lines self.create_line((imin, jmax, imax, jmax)) self.create_line((imin, jmin, imin, jmax)) self.create_line((imin, jmin, imax, jmin)) self.create_line((imax, jmin, imax, jmax)) self.create_line((imid, jmin, imid, jmax)) self.create_line((imin, jmid, imax, jmid)) # Create tick marks and labels tic = imin for label in xlabels: self.create_line((tic, jmax+ 5, tic, jmax)) self.create_text( tic, jmax+10, text=label) if len(xlabels)!=1: tic+=(imax-imin)/(len(xlabels)-1) tic = jmax for label in ylabels: self.create_line((imin , tic, imin-5, tic)) self.create_text( imin-20, tic, text=label) if len(ylabels)!=1: tic-=(jmax-jmin)/(len(ylabels)-1)
def fancy_show(y, cmap=''): x = y[0] y = y[1] plt.figure(0) for i in range(100): plt.subplot(10, 10, i+1) plt.imshow(x[i], cmap=cmap, interpolation='none') plt.axis('off') plt.figure(1) for i in range(100): plt.subplot(10, 10, i+1) plt.imshow(y[i], cmap=cmap, interpolation='none') plt.axis('off') plt.show()
def imshow_(x, **kwargs): if x.ndim == 2: plt.imshow(x, interpolation="nearest", **kwargs) elif x.ndim == 1: plt.imshow(x[:,None].T, interpolation="nearest", **kwargs) plt.yticks([]) plt.axis("tight") # ------------- Data -------------
def draw(name, feature4096): plt.figure(name) feature4096 = feature4096.reshape(64,64) for i, x in enumerate(feature4096): for j, y in enumerate(x): #if y <= 0: # print ' ', #else: # print '%3.1f'%y, plt.scatter([j],[i], s=[y*1000]) #print plt.axis([-1, 65, -1, 65]) plt.show()
def plot_validation_cost(train_error, val_error, class_rate=None, savefilename=None): epochs = range(len(train_error)) fig, ax1 = plt.subplots() ax1.plot(epochs, train_error, label='train error') ax1.plot(epochs, val_error, label='validation error') ax1.set_xlabel('epoch') ax1.set_ylabel('cost') plt.title('Validation Cost') lines = ax1.get_lines() # Shrink current axis's height by 10% on the bottom box = ax1.get_position() ax1.set_position([box.x0, box.y0 + box.height * 0.1, box.width, box.height * 0.9]) if class_rate is not None: ax2 = plt.twinx(ax1) ax2.plot(epochs, class_rate, label='classification rate', color='r') ax2.set_ylabel('classification rate') lines.extend(ax2.get_lines()) ax2.set_position([box.x0, box.y0 + box.height * 0.1, box.width, box.height * 0.9]) labels = [l.get_label() for l in lines] # Put a legend below current axis ax1.legend(lines, labels, loc='upper center', bbox_to_anchor=(0.5, -0.05), fancybox=False, shadow=False, ncol=5) # ax1.legend(lines, labels, loc='lower right') if savefilename: plt.savefig(savefilename) plt.show()
def visualize_sequence(sequence, shape=(30, 40), savefilename=None, title='sequence'): cols = int(math.ceil(len(sequence) / 2.0)) vis = tile_raster_images(sequence, shape, (2, cols), tile_spacing=(1, 1)) plt.title(title) plt.axis('off') plt.show() if savefilename: o = Image.fromarray(vis) o.save('{}.png'.format(savefilename))
def plotSequenceForRotatingState3D(degPerStep, Sigma, stationaryDim=0, T=1000): A = makeA_3DRotationMatrix(degPerStep, stationaryDim) Sigma = Sigma * np.eye(3) assert Sigma.shape == (3, 3) X = np.zeros((T, 3)) X[0, :] = [1, 1, 1] for t in xrange(1, T): X[t] = np.random.multivariate_normal(np.dot(A, X[t - 1]), Sigma) ax = Axes3D(pylab.figure()) pylab.plot(X[:, 0], X[:, 1], X[:, 2], '.') pylab.axis('equal')
def ricker_matrix(width, resolution, n_components): """Dictionary of Ricker (Mexican hat) wavelets""" centers = np.linspace(0, resolution - 1, n_components) D = np.empty((n_components, resolution)) for i, center in enumerate(centers): D[i] = ricker_function(resolution, center, width) D /= np.sqrt(np.sum(D ** 2, axis=1))[:, np.newaxis] return D
def get_2DsersicIeff(value,reff,sersicindex,axisratio,boxiness=0.0,returnFtot=False): """ Get the surface brightness value at the effective radius of a 2D sersic profile (given GALFIT Sersic parameters). Ieff is calculated using ewuations (4) and (5) in Peng et al. (2010), AJ 139:2097. This Ieff is what is referred to as 'amplitude' in astropy.modeling.models.Sersic2D used in tdose_utilities.gen_2Dsersic() --- INPUT --- value If returnFtot=False "value" corresponds to Ftot of the profile (total flux for profile integrated til r=infty) and Ieff will be returned. If instead returnFtot=True "value" should provide Ieff so Ftot can be returned reff Effective radius sersicindex Sersic index of profile axisratio Ratio between the minor and major axis (0<axisratio<1) boxiness The boxiness of the profile returnFtot If Ftot is not known, but Ieff is, set returnFtot=True to return Ftot instead (providing Ieff to "value") --- EXAMPLE OF USE --- Ieff = 1.0 reff = 25.0 sersicindex = 4.0 axisratio = 1.0 Ftot_calc = tu.get_2DsersicIeff(Ieff,reff,sersicindex,axisratio,returnFtot=True) Ieff_calc = tu.get_2DsersicIeff(Ftot_calc,reff,sersicindex,axisratio) size = 1000 x,y = np.meshgrid(np.arange(size), np.arange(size)) mod = Sersic2D(amplitude = Ieff, r_eff = reff, n=sersicindex, x_0=size/2.0, y_0=size/2.0, ellip=1-axisratio, theta=-1) img = mod(x, y) hducube = pyfits.PrimaryHDU(img) hdus = [hducube] hdulist = pyfits.HDUList(hdus) hdulist.writeto('/Volumes/DATABCKUP2/TDOSEextractions/models_cutouts/model_sersic_spherical.fits',clobber=True) """ gam2n = scipy.special.gamma(2.0*sersicindex) kappa = scipy.special.gammaincinv(2.0*sersicindex,0.5) Rfct = np.pi * (boxiness + 2.) / (4. * scipy.special.beta(1./(boxiness+2.),1.+1./(boxiness+2.)) ) factor = 2.0 * np.pi * reff**2.0 * np.exp(kappa) * sersicindex * kappa**(-2*sersicindex) * gam2n * axisratio / Rfct if returnFtot: Ftot = value * factor return Ftot else: Ieff = value / factor return Ieff # = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
def roll_2Dprofile(profile,position,padvalue=0.0,showprofiles=False): """ Move 2D profile to given position in array by rolling it in x and y. Note that the roll does not handle sub-pixel moves. tu.shift_2Dprofile() does this using interpolation --- INPUT --- profile profile to shift position position to move center of image (profile) to: [ypos,xpos] padvalue the values to padd the images with when shifting profile showprofiles Show profile when shifted? --- EXAMPLE OF USE --- tu.roll_2Dprofile(gauss2D,) """ profile_dim = profile.shape yroll = np.int(np.round(position[0]-profile_dim[0]/2.)) xroll = np.int(np.round(position[1]-profile_dim[1]/2.)) profile_shifted = np.roll(np.roll(profile,yroll,axis=0),xroll,axis=1) if showprofiles: vmaxval = np.max(profile_shifted) plt.imshow(profile_shifted,interpolation='none',vmin=-vmaxval, vmax=vmaxval) plt.title('Positioned Source') plt.show() if yroll < 0: profile_shifted[yroll:,:] = padvalue else: profile_shifted[:yroll,:] = padvalue if xroll < 0: profile_shifted[:,xroll:] = padvalue else: profile_shifted[:,:xroll] = padvalue if showprofiles: plt.imshow(profile_shifted,interpolation='none',vmin=-vmaxval, vmax=vmaxval) plt.title('Positioned Source with 0s inserted') plt.show() return profile_shifted # = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
def save_image(nifti, anat, cluster_dict, out_path, f, image_threshold=2, texcol=1, bgcol=0, iscale=2, text=None, **kwargs): '''Saves a single nifti image. Args: nifti (str or nipy.core.api.image.image.Image): nifti file to visualize. anat (nipy.core.api.image.image.Image): anatomical nifti file. cluster_dict (dict): dictionary of clusters. f (int): index. image_threshold (float): treshold for `plot_map`. texcol (float): text color. bgcol (float): background color. iscale (float): image scale. text (Optional[str]): text for figure. **kwargs: extra keyword arguments ''' if isinstance(nifti, str): nifti = load_image(nifti) feature = nifti.get_data() elif isinstance(nifti, nipy.core.image.image.Image): feature = nifti.get_data() font = {'size': 8} rc('font', **font) coords = cluster_dict['top_clust']['coords'] if coords == None: return feature /= feature.std() imax = np.max(np.absolute(feature)) imin = -imax imshow_args = dict( vmax=imax, vmin=imin, alpha=0.7 ) coords = ([-coords[0], -coords[1], coords[2]]) plt.axis('off') plt.text(0.05, 0.8, text, horizontalalignment='center', color=(texcol, texcol, texcol)) try: plot_map(feature, xyz_affine(nifti), anat=anat.get_data(), anat_affine=xyz_affine(anat), threshold=image_threshold, cut_coords=coords, annotate=False, cmap=cmap, draw_cross=False, **imshow_args) except Exception as e: return plt.savefig(out_path, transparent=True, facecolor=(bgcol, bgcol, bgcol))
def plot_colat_slice(self, component, colat, valmin, valmax, iteration=0, verbose=True): #- Some initialisations. ------------------------------------------------------------------ colat = np.pi * colat / 180.0 n_procs = self.setup["procs"]["px"] * self.setup["procs"]["py"] * self.setup["procs"]["pz"] vmax = float("-inf") vmin = float("inf") fig, ax = plt.subplots() #- Loop over processor boxes and check if colat falls within the volume. ------------------ for p in range(n_procs): if (colat >= self.theta[p,:].min()) & (colat <= self.theta[p,:].max()): #- Read this field and make lats & lons. ------------------------------------------ field = self.read_single_box(component,p,iteration) r, lon = np.meshgrid(self.z[p,:], self.phi[p,:]) x = r * np.cos(lon) y = r * np.sin(lon) #- Find the colat index and plot for this one box. -------------------------------- idx=min(np.where(min(np.abs(self.theta[p,:]-colat))==np.abs(self.theta[p,:]-colat))[0]) colat_effective = self.theta[p,idx]*180.0/np.pi #- Find min and max values. ------------------------------------------------------- vmax = max(vmax, field[idx,:,:].max()) vmin = min(vmin, field[idx,:,:].min()) #- Make a nice colourmap and plot. ------------------------------------------------ my_colormap=cm.make_colormap({0.0:[0.1,0.0,0.0], 0.2:[0.8,0.0,0.0], 0.3:[1.0,0.7,0.0],0.48:[0.92,0.92,0.92], 0.5:[0.92,0.92,0.92], 0.52:[0.92,0.92,0.92], 0.7:[0.0,0.6,0.7], 0.8:[0.0,0.0,0.8], 1.0:[0.0,0.0,0.1]}) cax = ax.pcolor(x, y, field[idx,:,:], cmap=my_colormap, vmin=valmin,vmax=valmax) #- Add colobar and title. ------------------------------------------------------------------ cbar = fig.colorbar(cax) if component in UNIT_DICT: cb.set_label(UNIT_DICT[component], fontsize="x-large", rotation=0) plt.suptitle("Vertical slice of %s at %i degree colatitude" % (component, colat_effective), size="large") plt.axis('equal') plt.show() #============================================================================================== #- Plot depth slice. #==============================================================================================
def double_lorentzian_fixedsplitting_testing(): # This method does not work and has to be fixed!!! for ii in range(1): # time.sleep(0.51) start=2800 stop=2950 num_points=int((stop-start)/2) x = np.linspace(start, stop, num_points) mod,params = qudi_fitting.make_multiplelorentzian_model(no_of_lor=2) p=Parameters() #============ Create data ========== p.add('c',value=100) p.add('lorentz0_amplitude',value=-abs(np.random.random(1)*50+100)) p.add('lorentz0_center',value=np.random.random(1)*150.0+2800) p.add('lorentz0_sigma',value=abs(np.random.random(1)*2.+1.)) p.add('lorentz1_center',value=p['lorentz0_center']+20) p.add('lorentz1_sigma',value=abs(np.random.random(1)*2.+1.)) p.add('lorentz1_amplitude',value=-abs(np.random.random(1)*50+100)) data_noisy=(mod.eval(x=x,params=p) + 2*np.random.normal(size=x.shape)) para=Parameters() result=qudi_fitting.make_doublelorentzian_fit(axis=x,data=data_noisy,add_parameters=para) data_smooth, offset = qudi_fitting.find_offset_parameter(x,data_noisy) data_level=data_smooth-offset #search for double lorentzian error, \ sigma0_argleft, dip0_arg, sigma0_argright, \ sigma1_argleft, dip1_arg , sigma1_argright = \ qudi_fitting._search_double_dip(x, data_level,make_prints=False) print(x[sigma0_argleft], x[dip0_arg], x[sigma0_argright], x[sigma1_argleft], x[dip1_arg], x[sigma1_argright]) print(x[dip0_arg], x[dip1_arg]) plt.plot((x[sigma0_argleft], x[sigma0_argleft]), ( data_noisy.min() ,data_noisy.max()), 'b-') plt.plot((x[sigma0_argright], x[sigma0_argright]), (data_noisy.min() ,data_noisy.max()), 'b-') plt.plot((x[sigma1_argleft], x[sigma1_argleft]), ( data_noisy.min() ,data_noisy.max()), 'k-') plt.plot((x[sigma1_argright], x[sigma1_argright]), ( data_noisy.min() ,data_noisy.max()), 'k-') try: plt.plot(x,data_noisy,'o') plt.plot(x,result.init_fit,'-y') plt.plot(x,result.best_fit,'-r',linewidth=2.0,) plt.plot(x,data_smooth,'-g') except: print('exception') plt.show()
def compute_inv_dft(x_val, y_val, zeropad_num=0): """ Compute the inverse Discrete fourier Transform @param numpy.array x_val: 1D array @param numpy.array y_val: 1D array of same size as x_val @param int zeropad_num: zeropadding (adding zeros to the end of the array). zeropad_num >= 0, the size of the array, which is add to the end of the y_val before performing the dft. The resulting array will have the length (len(y_val)/2)*(zeropad_num+1) Note that zeropadding will not change or add more information to the dft, it will solely interpolate between the dft_y values. @return: tuple(dft_x, dft_y): be aware that the return arrays' length depend on the zeropad number like len(dft_x) = len(dft_y) = (len(y_val)/2)*(zeropad_num+1) """ x_val = np.array(x_val) y_val = np.array(y_val) corrected_y = np.abs(y_val - y_val.max()) # The absolute values contain the fourier transformed y values # zeropad_arr = np.zeros(len(corrected_y)*(zeropad_num+1)) # zeropad_arr[:len(corrected_y)] = corrected_y fft_y = np.fft.ifft(corrected_y) # Due to the sampling theorem you can only identify frequencies at half # of the sample rate, therefore the FT contains an almost symmetric # spectrum (the asymmetry results from aliasing effects). Therefore take # the half of the values for the display. middle = int((len(corrected_y)+1)//2) # sample spacing of x_axis, if x is a time axis than it corresponds to a # timestep: x_spacing = np.round(x_val[-1] - x_val[-2], 12) # use the helper function of numpy to calculate the x_values for the # fourier space. That function will handle an occuring devision by 0: fft_x = np.fft.fftfreq(len(corrected_y), d=x_spacing) # return abs(fft_x[:middle]), fft_y[:middle] return fft_x[:middle], fft_y[:middle] # return fft_x, fft_y
def illustrate(Colors=Colors): if hasattr(Colors, 'colors'): Colors = Colors.colors from matplotlib import pylab rcParams = pylab.rcParams rcParams['pdf.fonttype'] = 42 rcParams['ps.fonttype'] = 42 rcParams['text.usetex'] = False rcParams['xtick.labelsize'] = 20 rcParams['ytick.labelsize'] = 20 rcParams['legend.fontsize'] = 25 import bnpy Data = get_data(T=1000, nDocTotal=8) for k in xrange(K): zmask = Data.TrueParams['Z'] == k pylab.plot(Data.X[zmask, 0], Data.X[zmask, 1], '.', color=Colors[k], markeredgecolor=Colors[k], alpha=0.4) sigEdges = np.flatnonzero(transPi[k] > 0.0001) for j in sigEdges: if j == k: continue dx = mus[j, 0] - mus[k, 0] dy = mus[j, 1] - mus[k, 1] pylab.arrow(mus[k, 0], mus[k, 1], 0.8 * dx, 0.8 * dy, head_width=2, head_length=4, facecolor=Colors[k], edgecolor=Colors[k]) tx = 0 - mus[k, 0] ty = 0 - mus[k, 1] xy = (mus[k, 0] - 0.2 * tx, mus[k, 1] - 0.2 * ty) ''' pylab.annotate( u'\u27F2', xy=(mus[k,0], mus[k,1]), color=Colors[k], fontsize=35, ) ''' pylab.gca().yaxis.set_ticks_position('left') pylab.gca().xaxis.set_ticks_position('bottom') pylab.axis('image') pylab.ylim([-38, 38]) pylab.xlim([-38, 38])