我们从Python开源项目中，提取了以下50个代码示例，用于说明如何使用numpy.angle()。
def getAngleInDegrees(_fft): _angle = np.angle(_fft, True); _rows, _columns = _fft.shape; for i in range(0, _rows): for j in range(0, _columns): if (_angle[i][j] < 0 and abs(_angle[i][j]) <= 180): #print(_angle[i][j]); #print("\n"); _angle[i][j] = 180 + _angle[i][j]; #print(_angle[i][j]); #print("\n"); elif (_angle[i][j] < 0 and abs(_angle[i][j]) > 180): #print(_angle[i][j]); #print("\n"); _angle[i][j] = 360 + _angle[i][j]; #print(_angle[i][j]); #print("\n"); return _angle;
def appres(F, H, sig, chg, taux, c, mu, eps, n): Res = np.zeros_like(F) Phase = np.zeros_like(F) App_ImpZ= np.zeros_like(F, dtype='complex_') for i in range(0, len(F)): UD, EH, Z , K = Propagate(F[i], H, sig, chg, taux, c, mu, eps, n) App_ImpZ[i] = EH[0, 1]/EH[1, 1] Res[i] = np.abs(App_ImpZ[i])**2./(mu_0*omega(F[i])) Phase[i] = np.angle(App_ImpZ[i], deg = True) return Res, Phase # Evaluate Up, Down components, E and H field, for a frequency range, # a discretized depth range and a time range (use to calculate envelope)
def phase_plot(data, toffset, log_scale, title): """Plot the phase of the data in linear or log scale.""" print("phase") phase = numpy.angle(data) / numpy.pi fig = plt.figure() ax = fig.add_subplot(1, 1, 1) ax.plot(phase) axmx = numpy.max(phase) #ax.axis([-axmx, axmx, -axmx, axmx]) ax.grid(True) ax.set_xlabel('time') ax.set_ylabel('phase') ax.set_title(title) return fig
def unknown_feature_extractor(x, sr, win_len, shift_len, barks, inner_win, inner_shift, win_type, method_version): x_spectrum = stft_extractor(x, win_len, shift_len, win_type) coef = get_fft_bark_mat(sr, win_len, barks, 20, sr//2) bark_spect = np.matmul(coef, x_spectrum) ams = np.zeros((barks, inner_win//2+1, (bark_spect.shape[1] - inner_win)//inner_shift)) for i in range(barks): channel_stft = stft_extractor(bark_spect[i, :], inner_win, inner_shift, 'hanning') if method_version == 'v1': ams[i, :, :] = 20 * np.log(np.abs(channel_stft[:inner_win//2+1, :(bark_spect.shape[1] - inner_win)//inner_shift])) elif method_version == 'v2': channel_amplitude = np.abs(channel_stft[:inner_win//2+1, :(bark_spect.shape[1] - inner_win)//inner_shift]) channel_angle = np.angle(channel_stft[:inner_win//2+1, :(bark_spect.shape[1] - inner_win)//inner_shift]) channel_angle = channel_angle - (np.floor(channel_angle / (2.*np.pi)) * (2.*np.pi)) ams[i, :, :] = np.power(channel_amplitude, 1./3.) * channel_angle else: ams[i, :, :] = np.abs(channel_stft) return ams
def ams_extractor(x, sr, win_len, shift_len, barks, inner_win, inner_shift, win_type, method_version): x_spectrum = stft_extractor(x, win_len, shift_len, win_type) coef = get_fft_bark_mat(sr, win_len, barks, 20, sr//2) bark_spect = np.matmul(coef, x_spectrum) ams = np.zeros((barks, inner_win//2+1, (bark_spect.shape[1] - inner_win)//inner_shift)) for i in range(barks): channel_stft = stft_extractor(bark_spect[i, :], inner_win, inner_shift, 'hanning') if method_version == 'v1': ams[i, :, :] = 20 * np.log(np.abs(channel_stft[:inner_win//2+1, :(bark_spect.shape[1] - inner_win)//inner_shift])) elif method_version == 'v2': channel_amplitude = np.abs(channel_stft[:inner_win//2+1, :(bark_spect.shape[1] - inner_win)//inner_shift]) channel_angle = np.angle(channel_stft[:inner_win//2+1, :(bark_spect.shape[1] - inner_win)//inner_shift]) channel_angle = channel_angle - (np.floor(channel_angle / (2.*np.pi)) * (2.*np.pi)) ams[i, :, :] = np.power(channel_amplitude, 1./3.) * channel_angle else: ams[i, :, :] = np.abs(channel_stft) return ams
def griffin_lim(mag, phase_angle, n_fft, hop, num_iters): """Iterative algorithm for phase retrival from a magnitude spectrogram. Args: mag: Magnitude spectrogram. phase_angle: Initial condition for phase. n_fft: Size of the FFT. hop: Stride of FFT. Defaults to n_fft/2. num_iters: Griffin-Lim iterations to perform. Returns: audio: 1-D array of float32 sound samples. """ fft_config = dict(n_fft=n_fft, win_length=n_fft, hop_length=hop, center=True) ifft_config = dict(win_length=n_fft, hop_length=hop, center=True) complex_specgram = inv_magphase(mag, phase_angle) for i in range(num_iters): audio = librosa.istft(complex_specgram, **ifft_config) if i != num_iters - 1: complex_specgram = librosa.stft(audio, **fft_config) _, phase = librosa.magphase(complex_specgram) phase_angle = np.angle(phase) complex_specgram = inv_magphase(mag, phase_angle) return audio
def __init__(self, qubit_names, quad="real"): super(PulseCalibration, self).__init__() self.qubit_names = qubit_names if isinstance(qubit_names, list) else [qubit_names] self.qubit = [QubitFactory(qubit_name) for qubit_name in qubit_names] if isinstance(qubit_names, list) else QubitFactory(qubit_names) self.filename = 'None' self.exp = None self.axis_descriptor = None self.cw_mode = False self.saved_settings = config.load_meas_file(config.meas_file) self.settings = deepcopy(self.saved_settings) #make a copy for used during calibration self.quad = quad if quad == "real": self.quad_fun = np.real elif quad == "imag": self.quad_fun = np.imag elif quad == "amp": self.quad_fun = np.abs elif quad == "phase": self.quad_fun = np.angle else: raise ValueError('Quadrature to calibrate must be one of ("real", "imag", "amp", "phase").') self.plot = self.init_plot()
def displayResult(self): fig = plt.figure(101) plt.subplot(221) plt.imshow(np.abs(self.reconstruction),origin='lower') plt.draw() plt.title('Reconstruction Magnitude') plt.subplot(222) plt.imshow(np.angle(self.reconstruction),origin='lower') plt.draw() plt.title('Reconstruction Phase') plt.subplot(223) plt.imshow(np.abs(self.aperture),origin='lower') plt.title("Aperture Magnitude") plt.draw() plt.subplot(224) plt.imshow(np.angle(self.aperture),origin='lower') plt.title("Aperture Phase") plt.draw() fig.canvas.draw() #fig.tight_layout() # display.display(fig) # display.clear_output(wait=True) # time.sleep(.00001)
def displayResult2(self): fig = plt.figure() ax = fig.add_subplot(2,2,1 ) ax.imshow(np.abs(self.reconstruction)) ax.set_title('Reconstruction Magnitude') ax = fig.add_subplot(2,2,2 ) ax.imshow(np.angle(self.reconstruction)) ax.set_title('Reconstruction Phase') ax = fig.add_subplot(2,2,3 ) ax.imshow(np.abs(self.aperture)) ax.set_title("Aperture Magnitude") ax = fig.add_subplot(2,2,4 ) ax.imshow(np.angle(self.aperture)) ax.set_title("Aperture Phase") fig.tight_layout()
def phase_to_color_wheel(complex_number): """Map a phase of a complexnumber to a color in (r,g,b). complex_number is phase is first mapped to angle in the range [0, 2pi] and then to a color wheel with blue at zero phase. """ angles = np.angle(complex_number) angle_round = int(((angles + 2 * np.pi) % (2 * np.pi))/np.pi*6) # print(angleround) color_map = {0: (0, 0, 1), # blue, 1: (0.5, 0, 1), # blue-violet 2: (1, 0, 1), # violet 3: (1, 0, 0.5), # red-violet, 4: (1, 0, 0), # red 5: (1, 0.5, 0), # red-oranage, 6: (1, 1, 0), # orange 7: (0.5, 1, 0), # orange-yellow 8: (0, 1, 0), # yellow, 9: (0, 1, 0.5), # yellow-green, 10: (0, 1, 1), # green, 11: (0, 0.5, 1) # green-blue, } return color_map[angle_round]
def time_sync(self, signal, start_est, stop_est): self.symbol_found = 0 self.symbol_start = 0 self.log.debug('time_sync()') corr_res = np.zeros( self.sym_len[self.mode], dtype=complex) for s in range(start_est, self.sym_len[self.mode]-1): corr_res[s] = self.gi_correlation(self.gi_len[self.mode], self.sym_len[self.mode], signal[s: s+self.sym_len[self.mode]] ) corr_max = np.argmax( np.abs(corr_res) ) self.symbol_found = 1 self.symbol_start = corr_max + 256 - 20 self.generated_samples = self.sym_len[self.mode] - self.gi_len[self.mode] # Tracking self.fine_freq_off = np.angle(corr_res[corr_max])/(2*np.pi*(self.sym_len[self.mode] - self.gi_len[self.mode])*(1.0/self.fs)) return
def compute_by_noise_pow(self, signal, n_pow): s_spec = np.fft.fftpack.fft(signal * self._window) s_amp = np.absolute(s_spec) s_phase = np.angle(s_spec) gamma = self._calc_aposteriori_snr(s_amp, n_pow) xi = self._calc_apriori_snr(gamma) self._prevGamma = gamma nu = gamma * xi / (1.0 + xi) self._G = (self._gamma15 * np.sqrt(nu) / gamma) * np.exp(-nu / 2.0) *\ ((1.0 + nu) * spc.i0(nu / 2.0) + nu * spc.i1(nu / 2.0)) idx = np.less(s_amp ** 2.0, n_pow) self._G[idx] = self._constant idx = np.isnan(self._G) + np.isinf(self._G) self._G[idx] = xi[idx] / (xi[idx] + 1.0) idx = np.isnan(self._G) + np.isinf(self._G) self._G[idx] = self._constant self._G = np.maximum(self._G, 0.0) amp = self._G * s_amp amp = np.maximum(amp, 0.0) amp2 = self._ratio * amp + (1.0 - self._ratio) * s_amp self._prevAmp = amp spec = amp2 * np.exp(s_phase * 1j) return np.real(np.fft.fftpack.ifft(spec))
def compute_by_noise_pow(self, signal, n_pow): s_spec = np.fft.fftpack.fft(signal * self._window) s_amp = np.absolute(s_spec) s_phase = np.angle(s_spec) gamma = self._calc_aposteriori_snr(s_amp, n_pow) xi = self._calc_apriori_snr(gamma) # xi = self._calc_apriori_snr2(gamma,n_pow) self._prevGamma = gamma nu = gamma * xi / (1.0 + xi) self._G = xi / (1.0 + xi) * np.exp(0.5 * spc.exp1(nu)) idx = np.less(s_amp ** 2.0, n_pow) self._G[idx] = self._constant idx = np.isnan(self._G) + np.isinf(self._G) self._G[idx] = xi[idx] / (xi[idx] + 1.0) idx = np.isnan(self._G) + np.isinf(self._G) self._G[idx] = self._constant self._G = np.maximum(self._G, 0.0) amp = self._G * s_amp amp = np.maximum(amp, 0.0) amp2 = self._ratio * amp + (1.0 - self._ratio) * s_amp self._prevAmp = amp spec = amp2 * np.exp(s_phase * 1j) return np.real(np.fft.fftpack.ifft(spec))
def compute_by_noise_pow(self, signal, n_pow): s_spec = np.fft.fftpack.fft(signal * self._window) s_amp = np.absolute(s_spec) s_phase = np.angle(s_spec) gamma = self._calc_aposteriori_snr(s_amp, n_pow) # xi = self._calc_apriori_snr2(gamma,n_pow) xi = self._calc_apriori_snr(gamma) self._prevGamma = gamma u = 0.5 - self._mu / (4.0 * np.sqrt(gamma * xi)) self._G = u + np.sqrt(u ** 2.0 + self._tau / (gamma * 2.0)) idx = np.less(s_amp ** 2.0, n_pow) self._G[idx] = self._constant idx = np.isnan(self._G) + np.isinf(self._G) self._G[idx] = xi[idx] / (xi[idx] + 1.0) idx = np.isnan(self._G) + np.isinf(self._G) self._G[idx] = self._constant self._G = np.maximum(self._G, 0.0) amp = self._G * s_amp amp = np.maximum(amp, 0.0) amp2 = self._ratio * amp + (1.0 - self._ratio) * s_amp self._prevAmp = amp spec = amp2 * np.exp(s_phase * 1j) return np.real(np.fft.fftpack.ifft(spec))
def plot_eigen_function(ax, se, n, psi1l, psi1r): # plt.figure(figsize=(8, 8)) for x in range(se.info['L']): for y in range(se.info['W']): i = x + y * se.info['L'] w = np.sqrt(10) * np.abs(psi1r[i, n]) arg = np.angle(psi1r[i, n]) circle = plt.Circle((x, y), w, color='black', zorder=10) ax.add_artist(circle) ax.plot([x, x + w * np.cos(arg)], [y, y + w * np.sin(arg)], color='white', lw=.8, zorder=12) ax.set_xlim([-.5, se.info['L'] - .5]) ax.set_ylim([-.5, se.info['W'] - .5]) # plt.show()
def single_spectrogram(inseq,fs,wlen,h,imag=False): """ imag: Return Imaginary Data of the STFT on True """ NFFT = int(2**(np.ceil(np.log2(wlen)))) K = np.sum(hamming(wlen, False))/wlen raw_data = inseq.astype('float32') raw_data = raw_data/np.amax(np.absolute(raw_data)) stft_data,_,_ = STFT(raw_data,wlen,h,NFFT,fs) s = np.absolute(stft_data)/wlen/K; if np.fmod(NFFT,2): s[1:,:] *=2 else: s[1:-2] *=2 real_data = np.transpose(20*np.log10(s + 10**-6)).astype(np.float32) if imag: imag_data = np.angle(stft_data).astype(np.float32) return real_data,imag_data return real_data
def complex2rgbalin(s): """ Displays complex image with intensity corresponding to the MODULUS and color (hsv) correponging to PHASE. From: pyVincent/ptycho.py """ ph=np.angle(s) t=np.pi/3 nx,ny=s.shape rgba=np.zeros((nx,ny,4)) rgba[:,:,0]=(ph<t)*(ph>-t) + (ph>t)*(ph<2*t)*(2*t-ph)/t + (ph>-2*t)*(ph<-t)*(ph+2*t)/t rgba[:,:,1]=(ph>t) + (ph<-2*t) *(-2*t-ph)/t+ (ph>0)*(ph<t) *ph/t rgba[:,:,2]=(ph<-t) + (ph>-t)*(ph<0) *(-ph)/t + (ph>2*t) *(ph-2*t)/t a=np.abs(s) a/=a.max() rgba[:,:,3]=a return rgba
def autocorr(trace): """This function takes an obspy trace object and performs a phase autocorrelation of the trace with itself. First, the hilbert transform is taken to obtain the analytic signal and hence the instantaneous phase. This is then passed to a fortran script where phase correlation is performed after Schimmel et al., 2011. This function relies on the shared object file phasecorr.so, which is the file containing the fortran subroutine for phase correlation. """ import numpy as np from scipy.signal import hilbert from phasecorr import phasecorr # Hilbert transform to obtain the analytic signal htrans = hilbert(trace) # Perform phase autocorrelation with instantaneous phase pac = phasecorr(np.angle(htrans),np.angle(htrans),len(htrans)) return pac
def crosscorr(tr1,tr2): """This function takes 2 obspy traces and performs a phase crosscorrelation of the traces. First, the hilbert transform is taken to obtain the analytic signal and hence the instantaneous phase. This is then passed to a fortran script where phase correlation is performed after Schimmel et al., 2011. This function relies on the shared object file phasecorr.so, which is the file containing the fortran subroutine for phase correlation. """ import numpy as np from scipy.signal import hilbert from phasecorr import phasecorr # Hilbert transform to obtain the analytic signal htrans1 = hilbert(tr1) htrans2 = hilbert(tr2) # Perform phase autocorrelation with instantaneous phase pcc = phasecorr(np.angle(htrans1),np.angle(htrans2),len(htrans1)) return pcc
def __init__(self, gain): self.gain = gain * np.pi # pi term puts angle output to pi range # components self.conjugate = Conjugate() self.complex_mult = ComplexMultiply() self.angle = Angle() self.out = Sfix() # specify component delay self._delay = self.conjugate._delay + \ self.complex_mult._delay + \ self.angle._delay + 1 # constants self.gain_sfix = Const(Sfix(self.gain, 3, -14))
def NMF(stft, n_sources): """ Sound source separation using NMF :param stft: the short-time Fourier Transform of the signal :param n_sources: the number of sources :returns: - Ys: sources """ print "Computing approximations" W, H = librosa.decompose.decompose(np.abs(stft), n_components=n_sources, sort=True) print "Reconstructing signals" Ys = list(scipy.outer(W[:,i], H[i])*np.exp(1j*np.angle(stft)) for i in xrange(n_sources)) print "Saving sound arrays" ys = [librosa.istft(Y) for Y in Ys] del Ys return ys
def undo_melspect(spect, sample_rate=SAMPLE_RATE, fps=FPS, frame_len=FRAME_LEN, min_freq=MIN_FREQ, max_freq=MAX_FREQ, invert_melbank_method='transpose', phases='random', normalize=False): """ Resynthesizes a mel spectrogram into a numpy array of samples. """ # undo logarithmic scaling spect = undo_logarithmize(spect) # undo Mel filterbank spect = undo_melfilter(spect, sample_rate, frame_len, min_freq, max_freq, invert_melbank_method) # randomize or reuse phases if phases == 'random': spect = spect * np.exp(np.pi*2.j*np.random.random(spect.shape)) elif phases is not None: spect = spect * np.exp(1.j * np.angle(phases)) # undo STFT hop_size = sample_rate / fps samples = undo_stft(spect, hop_size, frame_len) # normalize if needed if normalize: samples -= samples.mean() samples /= np.abs(samples).max() return samples.astype(np.float32)
def __new__(cls, realpart, imagpart=None): """Create a new EMArray.""" # Create complex obj if np.any(imagpart): obj = np.real(realpart) + 1j*np.real(imagpart) else: obj = np.asarray(realpart, dtype=complex) # Ensure its at least a 1D-Array, view cls obj = np.atleast_1d(obj).view(cls) # Store amplitude obj.amp = np.abs(obj) # Calculate phase, unwrap it, transform to degrees obj.pha = np.rad2deg(np.unwrap(np.angle(obj.real + 1j*obj.imag))) return obj # 2. Input parameter checks for modelling # 2.a <Check>s (alphabetically)
def ph_dec(m_phs, m_phc, mode='angle'): if mode == 'sign': m_bs = np.arcsin(m_phs) m_bc = np.arccos(m_phc) m_ph = np.sign(m_bs) * np.abs(m_bc) elif mode == 'angle': m_ph = np.angle(m_phc + m_phs * 1j) return m_ph #============================================================================== # From 'analysis_with_del_comp_and_ph_encoding_from_files' # f0_type: 'f0', 'lf0'
def phase_string(sig): """Take the angle and create a string as \pi if close to some \pi values""" if isinstance(sig, np.ndarray): sig = sig.ravel()[0] if isinstance(sig, np.complex): angle = np.angle(sig) else: angle = sig fractions = [Fraction(i, 12) for i in np.arange(-12, 12 + 1)] pi_multiples = np.array([np.pi * frac_to_float(f) for f in fractions]) strings = [frac_to_str(f) for f in fractions] if np.min(np.abs(angle - pi_multiples)) < 1e-3: return strings[np.argmin(np.abs(angle - pi_multiples))] else: return '%.2f' % angle
def cimshow(f_impulse): plt.figure(figsize=(7,5)) plt.subplot(2,2,1) plt.imshow(np.real(f_impulse)) plt.colorbar() plt.title('Re{}') plt.subplot(2,2,2) plt.imshow(np.imag(f_impulse)) plt.colorbar() plt.title('Img{}') plt.subplot(2,2,2+1) plt.imshow(np.abs(f_impulse)) plt.colorbar() plt.title('Magnitude') plt.subplot(2,2,2+2) plt.imshow(np.angle(f_impulse)) plt.colorbar() plt.title('Phase')
def update_recon_py(Recon1, support): err1 = 1.0 Constraint = np.ones(Recon1.shape) # cdef int R1y, R1x Recon1_abs = np.abs(Recon1) Recon1_pwr2 = np.power(Recon1_abs, 2) R1y, R1x = Recon1.shape for p in range(R1y): for q in range(R1x): if support[p, q] == 1: Constraint[p, q] = Recon1_abs[p, q] err1 += Recon1_pwr2[p, q] if Recon1_abs[p, q] > 1: Constraint[p, q] = 1 Recon1_update = Constraint * np.exp(1j * np.angle(Recon1)) return Recon1_update, err1
def imshow_GfpGbp(Gfp, Gbp): plt.subplot(2, 2, 1) plt.imshow(np.abs(Gfp), cmap='Greys_r') plt.title('abs(Gfp)') plt.subplot(2, 2, 2) plt.imshow(np.angle(Gfp), cmap='Greys_r') plt.title('angle(Gfp)') plt.subplot(2, 2, 1 + 2) plt.imshow(np.abs(Gbp), cmap='Greys_r') plt.title('abs(Gbp)') plt.subplot(2, 2, 2 + 2) plt.imshow(np.angle(Gbp), cmap='Greys_r') plt.title('angle(Gbp)')
def imshow_h(self): Gbp, Gfp = self.Gbp, self.Gfp plt.figure(figsize=[10, 10]) plt.subplot(2, 2, 1) plt.imshow(np.abs(Gbp), cmap='Greys_r', interpolation='none') plt.title('Gbp - Maganitute') plt.subplot(2, 2, 2) plt.imshow(np.angle(Gbp), cmap='Greys_r', interpolation='none') plt.title('Gbp - Angle') plt.subplot(2, 2, 2 + 1) plt.imshow(np.abs(Gfp), cmap='Greys_r', interpolation='none') plt.title('Gbf - Maganitute') plt.subplot(2, 2, 2 + 2) plt.imshow(np.angle(Gfp), cmap='Greys_r', interpolation='none') plt.title('Gbf - Angle')
def run_complex(self, Original): """ diffract original image and recon. ALso, the results will be shown """ Input = self.diffract_full(Original) Recon = self.reconstruct(Input) figure(figsize=(3 * 3 + 2, 3)) subplot(1, 4, 1) imshow(Original, 'Greys_r') title('Original') subplot(1, 4, 2) imshow(np.abs(Input), 'Greys_r') title('|Diffraction|') subplot(1, 4, 3) imshow(np.angle(Input), 'Greys_r') title('Angle(Diffraction)') subplot(1, 4, 4) imshow(np.abs(Recon), 'Greys_r') title('Modulus: |Recon|')
def fourierExtrapolation(x, n_predict): n = len(x) n_harm = 10 # number of harmonics in model t = np.arange(0, n) p = np.polyfit(t, x, 1) # find linear trend in x x_notrend = x - p[0] * t # detrended x x_freqdom = fft.fft(x_notrend) # detrended x in frequency domain f = fft.fftfreq(n) # frequencies indexes = list(range(n)) # sort indexes by frequency, lower -> higher indexes.sort(key = lambda i: np.absolute(f[i])) t = np.arange(0, n + n_predict) restored_sig = np.zeros(t.size) for i in indexes[:1 + n_harm * 2]: ampli = np.absolute(x_freqdom[i]) / n # amplitude phase = np.angle(x_freqdom[i]) # phase restored_sig += ampli * np.cos(2 * np.pi * f[i] * t + phase) return restored_sig + p[0] * t
def griffinlim(spectrogram, n_iter=50, window='hann', n_fft=2048, win_length=2048, hop_length=-1, verbose=False): if hop_length == -1: hop_length = n_fft // 4 angles = np.exp(2j * np.pi * np.random.rand(*spectrogram.shape)) t = tqdm(range(n_iter), ncols=100, mininterval=2.0, disable=not verbose) for i in t: full = np.abs(spectrogram).astype(np.complex) * angles inverse = librosa.istft(full, hop_length = hop_length, win_length = win_length, window = window) rebuilt = librosa.stft(inverse, n_fft = n_fft, hop_length = hop_length, win_length = win_length, window = window) angles = np.exp(1j * np.angle(rebuilt)) if verbose: diff = np.abs(spectrogram) - np.abs(rebuilt) t.set_postfix(loss=np.linalg.norm(diff, 'fro')) full = np.abs(spectrogram).astype(np.complex) * angles inverse = librosa.istft(full, hop_length = hop_length, win_length = win_length, window = window) return inverse
def fourierEx(x, n_predict, harmonics): n = len(x) # number of input samples n_harmonics = harmonics # f, 2*f, 3*f, .... n_harmonics ( 1,2, ) t = np.arange(0, n) # place to store data p = np.polyfit(t, x, 1) # find trend x_no_trend = x - p[0] * t x_frequency_domains = fft.fft(x_no_trend) f = np.fft.fftfreq(n) # frequencies indexes = list(range(n)) indexes.sort(key=lambda i: np.absolute(f[i])) t = np.arange(0, n + n_predict) restored_signal = np.zeros(t.size) for i in indexes[:1 + n_harmonics * 2]: amplitude = np.absolute(x_frequency_domains[i] / n) phase = np.angle(x_frequency_domains[i]) restored_signal += amplitude * np.cos(2 * np.pi * f[i] * t + phase) return restored_signal + p[0] * t
def processData(self, data): times = data.xvals('Time') dt = times[1]-times[0] data1 = data.asarray() ft = np.fft.fft(data1) ## determine frequencies in fft data df = 1.0 / (len(data1) * dt) freqs = np.linspace(0.0, (len(ft)-1) * df, len(ft)) ## flatten spikes at f0 and harmonics f0 = self.ctrls['f0'].value() for i in xrange(1, self.ctrls['harmonics'].value()+2): f = f0 * i # target frequency ## determine index range to check for this frequency ind1 = int(np.floor(f / df)) ind2 = int(np.ceil(f / df)) + (self.ctrls['samples'].value()-1) if ind1 > len(ft)/2.: break mag = (abs(ft[ind1-1]) + abs(ft[ind2+1])) * 0.5 for j in range(ind1, ind2+1): phase = np.angle(ft[j]) ## Must preserve the phase of each point, otherwise any transients in the trace might lead to large artifacts. re = mag * np.cos(phase) im = mag * np.sin(phase) ft[j] = re + im*1j ft[len(ft)-j] = re - im*1j data2 = np.fft.ifft(ft).real ma = metaarray.MetaArray(data2, info=data.infoCopy()) return ma
def removePeriodic(data, f0=60.0, dt=None, harmonics=10, samples=4): if (hasattr(data, 'implements') and data.implements('MetaArray')): data1 = data.asarray() if dt is None: times = data.xvals('Time') dt = times[1]-times[0] else: data1 = data if dt is None: raise Exception('Must specify dt for this data') ft = np.fft.fft(data1) ## determine frequencies in fft data df = 1.0 / (len(data1) * dt) freqs = np.linspace(0.0, (len(ft)-1) * df, len(ft)) ## flatten spikes at f0 and harmonics for i in xrange(1, harmonics + 2): f = f0 * i # target frequency ## determine index range to check for this frequency ind1 = int(np.floor(f / df)) ind2 = int(np.ceil(f / df)) + (samples-1) if ind1 > len(ft)/2.: break mag = (abs(ft[ind1-1]) + abs(ft[ind2+1])) * 0.5 for j in range(ind1, ind2+1): phase = np.angle(ft[j]) ## Must preserve the phase of each point, otherwise any transients in the trace might lead to large artifacts. re = mag * np.cos(phase) im = mag * np.sin(phase) ft[j] = re + im*1j ft[len(ft)-j] = re - im*1j data2 = np.fft.ifft(ft).real if (hasattr(data, 'implements') and data.implements('MetaArray')): return metaarray.MetaArray(data2, info=data.infoCopy()) else: return data2
def getMag(_fft): _mag = abs(_fft); norm_mag = _mag/_mag.max(); #angle = angle/angle.max(); return norm_mag;
def resample_2d_complex(array, sample_pts, query_pts): ''' Resamples a 2D complex array by interpolating the magnitude and phase independently and merging the results into a complex value. Args: array (numpy.ndarray): complex 2D array. sample_pts (tuple): pair of numpy.ndarray objects that contain the x and y sample locations, each array should be 1D. query_pts (tuple): points to interpolate onto, also 1D for each array. Returns: numpy.ndarray array resampled onto query_pts via bivariate spline ''' xq, yq = np.meshgrid(*query_pts) mag = abs(array) phase = np.angle(array) magfunc = interpolate.RegularGridInterpolator(sample_pts, mag) phasefunc = interpolate.RegularGridInterpolator(sample_pts, phase) interp_mag = magfunc((yq, xq)) interp_phase = phasefunc((yq, xq)) return interp_mag * exp(1j * interp_phase)
def phase(self): """ Return the phase spectrum. """ return np.angle(self)
def phase(z): val = np.angle(z) # val = np.rad2deg(np.unwrap(np.angle((z)))) return val
def DFT(x, w, N): """ Discrete Fourier Transformation(Analysis) of a given real input signal via an FFT implementation from scipy. Single channel is being supported. Args: x : (array) Real time domain input signal w : (array) Desired windowing function N : (int) FFT size Returns: magX : (2D ndarray) Magnitude Spectrum phsX : (2D ndarray) Phase Spectrum """ # Half spectrum size containing DC component hlfN = (N/2)+1 # Half window size. Two parameters to perform zero-phase windowing technique hw1 = int(math.floor((w.size+1)/2)) hw2 = int(math.floor(w.size/2)) # Window the input signal winx = x*w # Initialize FFT buffer with zeros and perform zero-phase windowing fftbuffer = np.zeros(N) fftbuffer[:hw1] = winx[hw2:] fftbuffer[-hw2:] = winx[:hw2] # Compute DFT via scipy's FFT implementation X = fft(fftbuffer) # Acquire magnitude and phase spectrum magX = (np.abs(X[:hlfN])) phsX = (np.angle(X[:hlfN])) return magX, phsX
def test_simple(self): x = np.array([1+0j, 1+2j]) y_r = np.log(np.abs(x)) + 1j * np.angle(x) y = np.log(x) for i in range(len(x)): assert_almost_equal(y[i], y_r[i])
def test_basic_ufuncs(self): # Test various functions such as sin, cos. (x, y, a10, m1, m2, xm, ym, z, zm, xf) = self.d assert_equal(np.cos(x), cos(xm)) assert_equal(np.cosh(x), cosh(xm)) assert_equal(np.sin(x), sin(xm)) assert_equal(np.sinh(x), sinh(xm)) assert_equal(np.tan(x), tan(xm)) assert_equal(np.tanh(x), tanh(xm)) assert_equal(np.sqrt(abs(x)), sqrt(xm)) assert_equal(np.log(abs(x)), log(xm)) assert_equal(np.log10(abs(x)), log10(xm)) assert_equal(np.exp(x), exp(xm)) assert_equal(np.arcsin(z), arcsin(zm)) assert_equal(np.arccos(z), arccos(zm)) assert_equal(np.arctan(z), arctan(zm)) assert_equal(np.arctan2(x, y), arctan2(xm, ym)) assert_equal(np.absolute(x), absolute(xm)) assert_equal(np.angle(x + 1j*y), angle(xm + 1j*ym)) assert_equal(np.angle(x + 1j*y, deg=True), angle(xm + 1j*ym, deg=True)) assert_equal(np.equal(x, y), equal(xm, ym)) assert_equal(np.not_equal(x, y), not_equal(xm, ym)) assert_equal(np.less(x, y), less(xm, ym)) assert_equal(np.greater(x, y), greater(xm, ym)) assert_equal(np.less_equal(x, y), less_equal(xm, ym)) assert_equal(np.greater_equal(x, y), greater_equal(xm, ym)) assert_equal(np.conjugate(x), conjugate(xm))
def synthesis_speech(noisy_speech, ideal_mask, win_type, win_len, shift_len, syn_method='A&R'): samples = noisy_speech.shape[0] frames = (samples - win_len) // shift_len if win_type == 'hanning': window = np.hanning(win_len) elif win_type == 'hamming': window = np.hamming(win_len) elif win_type == 'rectangle': window = np.ones(win_len) to_ifft = np.zeros(win_len, dtype=np.complex64) clean_speech = np.zeros((frames-1)*shift_len+win_len, dtype=np.float32) window_sum = np.zeros((frames-1)*shift_len+win_len, dtype=np.float32) for i in range(frames): one_frame = noisy_speech[i * shift_len: i * shift_len + win_len] windowed_frame = np.multiply(one_frame, window) stft = np.fft.fft(windowed_frame, win_len) masked_abs = np.abs(stft[:win_len//2+1]) * ideal_mask[:, i] to_ifft[:win_len//2+1] = masked_abs * np.exp(1j * np.angle(stft[:win_len//2+1])) to_ifft[win_len//2+1:] = np.conj(to_ifft[win_len//2-1:0:-1]) speech_seg = np.real(np.fft.ifft(to_ifft, win_len)) if syn_method == 'A&R' or syn_method == 'ALLEN & RABINER': clean_speech[i*shift_len:i*shift_len+win_len] += speech_seg window_sum[i*shift_len:i*shift_len+win_len] += window elif syn_method == 'G&L' or syn_method == 'GRIFFIN & LIM': speech_seg = np.multiply(speech_seg, window) clean_speech[i * shift_len:i * shift_len + win_len] += speech_seg window_sum[i * shift_len:i * shift_len + win_len] += np.power(window, 2.) # if i > 0: # clean_speech[i*shift_len: (i-1)*shift_len+win_len] *= 0.5 window_sum = np.where(window_sum < 1e-2, 1e-2, window_sum) return clean_speech / window_sum
def synthesis_speech(ns, mk, win_type, win_len, shift_len, syn_method='A&R'): samples = ns.shape[0] frames = (samples - win_len) // shift_len if win_type == 'hanning': window = np.hanning(win_len) elif win_type == 'hamming': window = np.hamming(win_len) elif win_type == 'rectangle': window = np.ones(win_len) to_ifft = np.zeros(win_len, dtype=np.complex64) clean_speech = np.zeros((frames-1)*shift_len+win_len, dtype=np.float32) window_sum = np.zeros((frames-1)*shift_len+win_len, dtype=np.float32) for i in range(frames): one_frame = ns[i * shift_len: i * shift_len + win_len] windowed_frame = np.multiply(one_frame, window) stft = np.fft.fft(windowed_frame, win_len) masked_abs = np.abs(stft[:win_len//2+1]) * mk[:, i] to_ifft[:win_len//2+1] = masked_abs * np.exp(1j * np.angle(stft[:win_len//2+1])) to_ifft[win_len//2+1:] = np.conj(to_ifft[win_len//2-1:0:-1]) speech_seg = np.real(np.fft.ifft(to_ifft, 320)) if syn_method == 'A&R' or syn_method == 'ALLEN & RABINER': clean_speech[i*shift_len:i*shift_len+win_len] += speech_seg window_sum[i*shift_len:i*shift_len+win_len] += window elif syn_method == 'G&L' or syn_method == 'GRIFFIN & LIM': speech_seg = np.multiply(speech_seg, window) clean_speech[i * shift_len:i * shift_len + win_len] += speech_seg window_sum[i * shift_len:i * shift_len + win_len] += np.power(window, 2.) # if i > 0: # clean_speech[i*shift_len: (i-1)*shift_len+win_len] *= 0.5 window_sum = np.where(window_sum < 1e-2, 1e-2, window_sum) return clean_speech / window_sum
