Python scipy.ndimage.filters.gaussian_filter1d() Examples

The following are 30 code examples of scipy.ndimage.filters.gaussian_filter1d(). You can vote up the ones you like or vote down the ones you don't like, and go to the original project or source file by following the links above each example. You may also want to check out all available functions/classes of the module scipy.ndimage.filters , or try the search function .
Example #1
Source File: visualization.py    From audio-reactive-led-strip with MIT License 6 votes vote down vote up
def visualize_scroll(y):
    """Effect that originates in the center and scrolls outwards"""
    global p
    y = y**2.0
    gain.update(y)
    y /= gain.value
    y *= 255.0
    r = int(np.max(y[:len(y) // 3]))
    g = int(np.max(y[len(y) // 3: 2 * len(y) // 3]))
    b = int(np.max(y[2 * len(y) // 3:]))
    # Scrolling effect window
    p[:, 1:] = p[:, :-1]
    p *= 0.98
    p = gaussian_filter1d(p, sigma=0.2)
    # Create new color originating at the center
    p[0, 0] = r
    p[1, 0] = g
    p[2, 0] = b
    # Update the LED strip
    return np.concatenate((p[:, ::-1], p), axis=1) 
Example #2
Source File: smooth_bbox.py    From human_dynamics with BSD 2-Clause "Simplified" License 6 votes vote down vote up
def smooth_bbox_params(bbox_params, kernel_size=11, sigma=8):
    """
    Applies median filtering and then gaussian filtering to bounding box
    parameters.

    Args:
        bbox_params (Nx3): [cx, cy, scale].
        kernel_size (int): Kernel size for median filtering (must be odd).
        sigma (float): Sigma for gaussian smoothing.

    Returns:
        Smoothed bounding box parameters (Nx3).
    """
    smoothed = np.array([signal.medfilt(param, kernel_size)
                         for param in bbox_params.T]).T
    return np.array([gaussian_filter1d(traj, sigma) for traj in smoothed.T]).T 
Example #3
Source File: pid_analysis.py    From flight_review with BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
def stackspectrum(self, time, throttle, trace, window):
        ### calculates spectrogram from stack of windows against throttle.
        # slicing off last 2s to get rid of landing
        gyro = trace[:-int(Trace.noise_superpos*2./Trace.noise_framelen),:] * window
        thr = throttle[:-int(Trace.noise_superpos*2./Trace.noise_framelen),:] * window
        time = time[:-int(Trace.noise_superpos*2./Trace.noise_framelen),:]

        freq, spec = self.spectrum(time[0], gyro)

        weights = abs(spec.real)
        avr_thr = np.abs(thr).max(axis=1)

        hist2d=self.hist2d(avr_thr, freq,weights,[101,len(freq)/4])

        filt_width = 3  # width of gaussian smoothing for hist data
        hist2d_sm = gaussian_filter1d(hist2d['hist2d_norm'], filt_width, axis=1, mode='constant')

        # get max value in histogram >100hz
        thresh = 100.
        mask = self.to_mask(freq[:-1:4].clip(thresh-1e-9,thresh))
        maxval = np.max(hist2d_sm.transpose()*mask)

        return {'throt_hist_avr':hist2d['throt_hist'],'throt_axis':hist2d['throt_scale'],'freq_axis':freq[::4],
                'hist2d_norm':hist2d['hist2d_norm'], 'hist2d_sm':hist2d_sm, 'hist2d':hist2d['hist2d'], 'max':maxval} 
Example #4
Source File: segmenter.py    From msaf with MIT License 6 votes vote down vote up
def pick_peaks(nc, L=16):
    """Obtain peaks from a novelty curve using an adaptive threshold."""
    offset = nc.mean() / 20.

    nc = filters.gaussian_filter1d(nc, sigma=4)  # Smooth out nc

    th = filters.median_filter(nc, size=L) + offset
    #th = filters.gaussian_filter(nc, sigma=L/2., mode="nearest") + offset

    peaks = []
    for i in range(1, nc.shape[0] - 1):
        # is it a peak?
        if nc[i - 1] < nc[i] and nc[i] > nc[i + 1]:
            # is it above the threshold?
            if nc[i] > th[i]:
                peaks.append(i)
    #plt.plot(nc)
    #plt.plot(th)
    #for peak in peaks:
        #plt.axvline(peak)
    #plt.show()

    return peaks 
Example #5
Source File: transitionModels.py    From bayesloop with MIT License 6 votes vote down vote up
def computeForwardPrior(self, posterior, t):
        """
        Compute new prior from old posterior (moving forwards in time).

        Args:
            posterior(ndarray): Parameter distribution from current time step
            t(int): integer time step

        Returns:
            ndarray: Prior parameter distribution for subsequent time step
        """
        axisToTransform = self.study.observationModel.parameterNames.index(self.selectedParameter)
        normedSigma = self.hyperParameterValues[0]/self.latticeConstant[axisToTransform]
        
        if normedSigma > 0.:
            newPrior = gaussian_filter1d(posterior, normedSigma, axis=axisToTransform)
        else:
            newPrior = posterior.copy()
        
        return newPrior 
Example #6
Source File: instrument_model.py    From isofit with Apache License 2.0 5 votes vote down vote up
def _high_frequency_vert(X, sigma=4.0):
    """."""

    nl, nb, nr = X.shape
    Xvert = X.copy()
    for r in range(nr):
        for b in range(nb):
            filt = gaussian_filter1d(Xvert[:, b, r], sigma, mode='nearest')
            Xvert[:, b, r] = X[:, b, r] - filt
    return Xvert 
Example #7
Source File: instrument_model.py    From isofit with Apache License 2.0 5 votes vote down vote up
def _low_frequency_horiz(X, sigma=4.0):
    """."""

    nl, nb, nr = X.shape
    Xhoriz = X.copy()
    for l in range(nl):
        for b in range(nb):
            Xhoriz[l, b, :] = gaussian_filter1d(
                Xhoriz[l, b, :], sigma, mode='nearest')
    return Xhoriz 
Example #8
Source File: cnn.py    From anticipating-activities with MIT License 5 votes vote down vote up
def __post_process(self, result, sigma):
        new_res =  gaussian_filter1d(result, sigma=sigma, axis=0)
        return new_res 
Example #9
Source File: util.py    From spm1d with GNU General Public License v3.0 5 votes vote down vote up
def smooth(Y, fwhm=5.0):
	'''
	Smooth a set of 1D continua.
	This method uses **scipy.ndimage.filters.gaussian_filter1d** but uses the *fwhm*
	instead of the standard deviation.
	
	:Parameters:
	
	- *Y* --- a (J x Q) numpy array
	- *fwhm* ---  Full-width at half-maximum of a Gaussian kernel used for smoothing.
	
	:Returns:
	
	- (J x Q) numpy array
	
	:Example:
	
	>>> Y0  = np.random.rand(5, 101)
	>>> Y   = spm1d.util.smooth(Y0, fwhm=10.0)
	
	.. note:: A Gaussian kernel's *fwhm* is related to its standard deviation (*sd*) as follows:
	
	>>> fwhm = sd * sqrt(8*log(2))
	'''
	sd    = fwhm / sqrt(8*log(2))
	return gaussian_filter1d(Y, sd, mode='wrap') 
Example #10
Source File: rdf.py    From vasppy with MIT License 5 votes vote down vote up
def smeared_rdf(self,sigma=0.1):
        """
        Smear the RDF with a Gaussian kernel.

        Args:
            sigma (:obj:`float`, optional): Standard deviation for Gaussian kernel. 
                Optional, default is 0.1.

        Returns:
            (np.array): Smeared RDF data.

        """
        sigma_n_bins = sigma / self.dr
        return gaussian_filter1d(self.rdf, sigma=sigma_n_bins) 
Example #11
Source File: rdf.py    From vasppy with MIT License 5 votes vote down vote up
def smeared_gsrt(self,sigma=0.1):
        """
        Smear the self part of the Van Hove correlation function with a Gaussian kernel.

        Args:
            sigma (:obj:`float`, optional): Standard deviation for Gaussian kernel. Optional, default is 0.1.

        Returns:
            (np.array): Smeared data.

        """
        sigma_n_bins = sigma / self.dr
        return gaussian_filter1d(self.gsrt, sigma=sigma_n_bins) 
Example #12
Source File: rdf.py    From vasppy with MIT License 5 votes vote down vote up
def smeared_gdrt(self,sigma=0.1):
        """
        Smear the distinct part of the Van Hove correlation function with a Gaussian kernel.

        Args:
            sigma (:obj:`float`, optional): Standard deviation for Gaussian kernel. Optional, default is 0.1.
        
        Returns:
            (np.array): Smeared data.

        """
        sigma_n_bins = sigma / self.dr
        return gaussian_filter1d(self.gdrt, sigma=sigma_n_bins) 
Example #13
Source File: main.py    From Systematic-LEDs with MIT License 5 votes vote down vote up
def visualize_scroll(self, y):
        # Effect that scrolls colours corresponding to frequencies across the strip 
        #y = y**1.5
        n_pixels = config.settings["devices"][self.board]["configuration"]["N_PIXELS"]
        y = np.copy(interpolate(y, n_pixels // 2))
        board_manager.signal_processers[self.board].common_mode.update(y)
        diff = y - self.prev_spectrum
        self.prev_spectrum = np.copy(y)

        y = np.clip(y, 0, 1)
        lows = y[:len(y) // 6]
        mids = y[len(y) // 6: 2 * len(y) // 5]
        high = y[2 * len(y) // 5:]
        # max values
        lows_max = np.max(lows)#*config.settings["devices"][self.board]["effect_opts"]["Scroll"]["lows_multiplier"])
        mids_max = float(np.max(mids))#*config.settings["devices"][self.board]["effect_opts"]["Scroll"]["mids_multiplier"])
        high_max = float(np.max(high))#*config.settings["devices"][self.board]["effect_opts"]["Scroll"]["high_multiplier"])
        # indexes of max values
        # map to colour gradient
        lows_val = (np.array(colour_manager.colour(config.settings["devices"][self.board]["effect_opts"]["Scroll"]["lows_color"])) * lows_max).astype(int)
        mids_val = (np.array(colour_manager.colour(config.settings["devices"][self.board]["effect_opts"]["Scroll"]["mids_color"])) * mids_max).astype(int)
        high_val = (np.array(colour_manager.colour(config.settings["devices"][self.board]["effect_opts"]["Scroll"]["high_color"])) * high_max).astype(int)
        # Scrolling effect window
        speed = config.settings["devices"][self.board]["effect_opts"]["Scroll"]["speed"]
        self.output[:, speed:] = self.output[:, :-speed]
        self.output = (self.output * config.settings["devices"][self.board]["effect_opts"]["Scroll"]["decay"]).astype(int)
        self.output = gaussian_filter1d(self.output, sigma=config.settings["devices"][self.board]["effect_opts"]["Scroll"]["blur"])
        # Create new color originating at the center
        self.output[0, :speed] = lows_val[0] + mids_val[0] + high_val[0]
        self.output[1, :speed] = lows_val[1] + mids_val[1] + high_val[1]
        self.output[2, :speed] = lows_val[2] + mids_val[2] + high_val[2]
        # Update the LED strip
        #return np.concatenate((self.prev_spectrum[:, ::-speed], self.prev_spectrum), axis=1)
        if config.settings["devices"][self.board]["effect_opts"]["Scroll"]["mirror"]:
            p = np.concatenate((self.output[:, ::-2], self.output[:, ::2]), axis=1)
        else:
            p = self.output
        return p 
Example #14
Source File: main.py    From Systematic-LEDs with MIT License 5 votes vote down vote up
def visualize_energy(self, y):
        """Effect that expands from the center with increasing sound energy"""
        y = np.copy(y)
        board_manager.signal_processers[self.board].gain.update(y)
        y /= board_manager.signal_processers[self.board].gain.value
        scale = config.settings["devices"][self.board]["effect_opts"]["Energy"]["scale"]
        # Scale by the width of the LED strip
        y *= float((config.settings["devices"][self.board]["configuration"]["N_PIXELS"] * scale) - 1)
        y = np.copy(interpolate(y, config.settings["devices"][self.board]["configuration"]["N_PIXELS"] // 2))
        # Map color channels according to energy in the different freq bands
        #y = np.copy(interpolate(y, config.settings["devices"][self.board]["configuration"]["N_PIXELS"] // 2))
        diff = y - self.prev_spectrum
        self.prev_spectrum = np.copy(y)
        spectrum = np.copy(self.prev_spectrum)
        spectrum = np.array([j for i in zip(spectrum,spectrum) for j in i])
        # Color channel mappings
        r = int(np.mean(spectrum[:len(spectrum) // 3]**scale)*config.settings["devices"][self.board]["effect_opts"]["Energy"]["r_multiplier"])
        g = int(np.mean(spectrum[len(spectrum) // 3: 2 * len(spectrum) // 3]**scale)*config.settings["devices"][self.board]["effect_opts"]["Energy"]["g_multiplier"])
        b = int(np.mean(spectrum[2 * len(spectrum) // 3:]**scale)*config.settings["devices"][self.board]["effect_opts"]["Energy"]["b_multiplier"])
        # Assign color to different frequency regions
        self.output[0, :r] = 255
        self.output[0, r:] = 0
        self.output[1, :g] = 255
        self.output[1, g:] = 0
        self.output[2, :b] = 255
        self.output[2, b:] = 0
        # Apply blur to smooth the edges
        self.output[0, :] = gaussian_filter1d(self.output[0, :], sigma=config.settings["devices"][self.board]["effect_opts"]["Energy"]["blur"])
        self.output[1, :] = gaussian_filter1d(self.output[1, :], sigma=config.settings["devices"][self.board]["effect_opts"]["Energy"]["blur"])
        self.output[2, :] = gaussian_filter1d(self.output[2, :], sigma=config.settings["devices"][self.board]["effect_opts"]["Energy"]["blur"])
        if config.settings["devices"][self.board]["effect_opts"]["Energy"]["mirror"]:
            p = np.concatenate((self.output[:, ::-2], self.output[:, ::2]), axis=1)
        else:
            p = self.output
        return p 
Example #15
Source File: main.py    From Systematic-LEDs with MIT License 5 votes vote down vote up
def visualize_wavelength(self, y):
        y = np.copy(interpolate(y, config.settings["devices"][self.board]["configuration"]["N_PIXELS"] // 2))
        board_manager.signal_processers[self.board].common_mode.update(y)
        diff = y - self.prev_spectrum
        self.prev_spectrum = np.copy(y)
        # Color channel mappings
        r = board_manager.signal_processers[self.board].r_filt.update(y - board_manager.signal_processers[self.board].common_mode.value)
        g = np.abs(diff)
        b = board_manager.signal_processers[self.board].b_filt.update(np.copy(y))
        r = np.array([j for i in zip(r,r) for j in i])
        output = np.array([colour_manager.full_gradients[self.board][config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["color_mode"]][0][
                                    (config.settings["devices"][self.board]["configuration"]["N_PIXELS"] if config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["reverse_grad"] else 0):
                                    (None if config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["reverse_grad"] else config.settings["devices"][self.board]["configuration"]["N_PIXELS"]):]*r,
                           colour_manager.full_gradients[self.board][config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["color_mode"]][1][
                                    (config.settings["devices"][self.board]["configuration"]["N_PIXELS"] if config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["reverse_grad"] else 0):
                                    (None if config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["reverse_grad"] else config.settings["devices"][self.board]["configuration"]["N_PIXELS"]):]*r,
                           colour_manager.full_gradients[self.board][config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["color_mode"]][2][
                                    (config.settings["devices"][self.board]["configuration"]["N_PIXELS"] if config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["reverse_grad"] else 0):
                                    (None if config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["reverse_grad"] else config.settings["devices"][self.board]["configuration"]["N_PIXELS"]):]*r])
        #self.prev_spectrum = y
        colour_manager.full_gradients[self.board][config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["color_mode"]] = np.roll(
                    colour_manager.full_gradients[self.board][config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["color_mode"]],
                    config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["roll_speed"]*(-1 if config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["reverse_roll"] else 1),
                    axis=1)
        output[0] = gaussian_filter1d(output[0], sigma=config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["blur"])
        output[1] = gaussian_filter1d(output[1], sigma=config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["blur"])
        output[2] = gaussian_filter1d(output[2], sigma=config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["blur"])
        if config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["flip_lr"]:
            output = np.fliplr(output)
        if config.settings["devices"][self.board]["effect_opts"]["Wavelength"]["mirror"]:
            output = np.concatenate((output[:, ::-2], output[:, ::2]), axis=1)
        return output 
Example #16
Source File: main.py    From Systematic-LEDs with MIT License 5 votes vote down vote up
def update(self, audio_samples):
        """ Return processed audio data
        Returns mel curve, x/y data
        This is called every time there is a microphone update
        Returns
        -------
        audio_data : dict
            Dict containinng "mel", "vol", "x", and "y"
        """
        audio_data = {}
        # Normalize samples between 0 and 1
        y = audio_samples / 2.0**15
        # Construct a rolling window of audio samples
        self.y_roll[:-1] = self.y_roll[1:]
        self.y_roll[-1, :] = np.copy(y)
        y_data = np.concatenate(self.y_roll, axis=0).astype(np.float32)
        vol = np.max(np.abs(y_data))
        # Transform audio input into the frequency domain
        N = len(y_data)
        N_zeros = 2**int(np.ceil(np.log2(N))) - N
        # Pad with zeros until the next power of two
        y_data *= self.fft_window
        y_padded = np.pad(y_data, (0, N_zeros), mode='constant')
        YS = np.abs(np.fft.rfft(y_padded)[:N // 2])
        # Construct a Mel filterbank from the FFT data
        mel = np.atleast_2d(YS).T * self.mel_y.T
        # Scale data to values more suitable for visualization
        mel = np.sum(mel, axis=0)
        mel = mel**2
        # Gain normalization
        self.mel_gain.update(np.max(gaussian_filter1d(mel, sigma=0.1)))
        mel /= self.mel_gain.value
        mel = self.mel_smoothing.update(mel)
        x = np.linspace(config.settings["devices"][self.board]["configuration"]["MIN_FREQUENCY"], config.settings["devices"][self.board]["configuration"]["MAX_FREQUENCY"], len(mel))
        y = self.fft_plot_filter.update(mel)

        audio_data["mel"] = mel
        audio_data["vol"] = vol
        audio_data["x"]   = x
        audio_data["y"]   = y
        return audio_data 
Example #17
Source File: visualization.py    From audio-reactive-led-strip with MIT License 5 votes vote down vote up
def visualize_energy(y):
    """Effect that expands from the center with increasing sound energy"""
    global p
    y = np.copy(y)
    gain.update(y)
    y /= gain.value
    # Scale by the width of the LED strip
    y *= float((config.N_PIXELS // 2) - 1)
    # Map color channels according to energy in the different freq bands
    scale = 0.9
    r = int(np.mean(y[:len(y) // 3]**scale))
    g = int(np.mean(y[len(y) // 3: 2 * len(y) // 3]**scale))
    b = int(np.mean(y[2 * len(y) // 3:]**scale))
    # Assign color to different frequency regions
    p[0, :r] = 255.0
    p[0, r:] = 0.0
    p[1, :g] = 255.0
    p[1, g:] = 0.0
    p[2, :b] = 255.0
    p[2, b:] = 0.0
    p_filt.update(p)
    p = np.round(p_filt.value)
    # Apply substantial blur to smooth the edges
    p[0, :] = gaussian_filter1d(p[0, :], sigma=4.0)
    p[1, :] = gaussian_filter1d(p[1, :], sigma=4.0)
    p[2, :] = gaussian_filter1d(p[2, :], sigma=4.0)
    # Set the new pixel value
    return np.concatenate((p[:, ::-1], p), axis=1) 
Example #18
Source File: optics.py    From sumo with MIT License 5 votes vote down vote up
def broaden_eps(dielectric, sigma):
    """Apply gaussian broadening to the dielectric response.

    Args:
        dielectric_data (tuple): The high-frequency dielectric data, following
            the same format as
            :attr:`pymatgen.io.vasp.outputs.Vasprun.dielectric`.
            This is a :obj:`tuple` containing the energy, the real part of the
            dielectric tensor, and the imaginary part of the tensor, as a
            :obj:`list` of :obj:`floats`. E.g.::

                (
                    [energies],
                    [[real_xx, real_yy, real_zz, real_xy, real_yz, real_xz]],
                    [[imag_xx, imag_yy, imag_zz, imag_xy, imag_yz, imag_xz]]
                )

        sigma (float): Standard deviation for gaussian broadening.

    Returns:
        :obj:`tuple` of :obj:`list` of :obj:`list` of :obj:`float`: The
        broadened dielectric response. Returned as a tuple containing the
        energy, the real part of the dielectric tensor, and the imaginary
        part of the tensor. E.g.::

            (
                [energies],
                [[real_xx, real_yy, real_zz, real_xy, real_yz, real_xz]],
                [[imag_xx, imag_yy, imag_zz, imag_xy, imag_yz, imag_xz]]
            )
    """
    e = dielectric[0]
    diff = [e[i + 1] - e[i] for i in range(len(e) - 1)]
    diff_avg = sum(diff) / len(diff)
    real = [gaussian_filter1d(np.array(dielectric[1])[:, x], sigma / diff_avg)
            for x in range(6)]
    imag = [gaussian_filter1d(np.array(dielectric[2])[:, x], sigma / diff_avg)
            for x in range(6)]

    return (e, np.array(real).T, np.array(imag).T) 
Example #19
Source File: transformations.py    From ronin with GNU General Public License v3.0 5 votes vote down vote up
def __call__(self, feat, targ, **kwargs):
        sigma = np.random.random() * self.max_sigma
        return gaussian_filter1d(feat, sigma=sigma, axis=0), targ 
Example #20
Source File: bam_cov.py    From basenji with Apache License 2.0 5 votes vote down vote up
def gc_normalize(self, chrom, coverage):
    """ Apply a model to normalize for GC content.

        In
         chrom (str): Chromosome
         coverage ([float]): Coverage array
         pseudocount (int): Coverage pseudocount.

        Out
         model (sklearn object): To control for GC%.
        """

    # trim chromosome strand
    if self.stranded and chrom[-1] in '+-':
      chrom = chrom[:-1]

    # get sequence
    seq = self.fasta.fetch(chrom)
    assert (len(seq) == len(coverage))

    # compute GC boolean vector
    seq_gc = np.array([nt in 'CG' for nt in seq], dtype='float32')

    # gaussian filter1d
    seq_gc_gauss = gaussian_filter1d(seq_gc, sigma=self.fragment_sd, truncate=3)

    # compute norm quantity
    seq_gc_norm = self.gc_model.predict(seq_gc_gauss[:, np.newaxis])

    # apply it
    return coverage * np.exp(-seq_gc_norm + self.gc_base) 
Example #21
Source File: plot_E.py    From fedavgpy with MIT License 5 votes vote down vote up
def smooth(y):
    return gaussian_filter1d(y, sigma=0.6) 
Example #22
Source File: deskew.py    From lpr with Apache License 2.0 5 votes vote down vote up
def skew_detection(image_gray):
    h, w = image_gray.shape[:2]
    eigen = cv2.cornerEigenValsAndVecs(image_gray,12, 5)
    angle_sur = np.zeros(180,np.uint)
    eigen = eigen.reshape(h, w, 3, 2)
    flow = eigen[:,:,2]
    vis = image_gray.copy()
    vis[:] = (192 + np.uint32(vis)) / 2
    d = 12
    points =  np.dstack( np.mgrid[d/2:w:d, d/2:h:d] ).reshape(-1, 2)
    for x, y in points:
        vx, vy = np.int32(flow[int(y), int(x)]*d)
        # cv2.line(rgb, (x-vx, y-vy), (x+vx, y+vy), (0, 355, 0), 1, cv2.LINE_AA)
        ang = angle(vx,vy)
        angle_sur[(ang+180)%180] +=1

    # torr_bin = 30
    angle_sur = angle_sur.astype(np.float)
    angle_sur = (angle_sur-angle_sur.min())/(angle_sur.max()-angle_sur.min())
    angle_sur = filters.gaussian_filter1d(angle_sur,5)
    skew_v_val =  angle_sur[20:180-20].max()
    skew_v = angle_sur[30:180-30].argmax() + 30
    skew_h_A = angle_sur[0:30].max()
    skew_h_B = angle_sur[150:180].max()
    skew_h = 0
    if (skew_h_A > skew_v_val*0.3 or skew_h_B > skew_v_val*0.3):
        if skew_h_A>=skew_h_B:
            skew_h = angle_sur[0:20].argmax()
        else:
            skew_h = - angle_sur[160:180].argmax()
    return skew_h,skew_v 
Example #23
Source File: smooth_bbox.py    From phd with BSD 2-Clause "Simplified" License 5 votes vote down vote up
def smooth_bbox_params(bbox_params, kernel_size=11, sigma=8):
    """
    Applies median filtering and then gaussian filtering to bounding box
    parameters.
    Args:
        bbox_params (Nx3): [cx, cy, scale].
        kernel_size (int): Kernel size for median filtering (must be odd).
        sigma (float): Sigma for gaussian smoothing.
    Returns:
        Smoothed bounding box parameters (Nx3).
    """
    smoothed = np.array([signal.medfilt(param, kernel_size)
                         for param in bbox_params.T]).T
    return np.array([gaussian_filter1d(traj, sigma) for traj in smoothed.T]).T 
Example #24
Source File: contrib.py    From norbert with MIT License 5 votes vote down vote up
def smooth(v, width=1, temporal=False):
    """
    smoothes a ndarray with a Gaussian blur.

    Parameters
    ----------
    v: np.ndarray [shape=(nb_frames, ...)]
        input array

    sigma: int [scalar]
        lengthscale of the gaussian blur

    temporal: boolean
        if True, will smooth only along time through 1d blur. Will use a
        multidimensional Gaussian blur otherwise.

    Returns
    -------
    result: np.ndarray [shape=(nb_frames, ...)]
        filtered array

    """
    if temporal:
        return gaussian_filter1d(v, sigma=width, axis=0)
    else:
        return gaussian_filter(v, sigma=width, truncate=width) 
Example #25
Source File: spectools.py    From specidentify with BSD 3-Clause "New" or "Revised" License 5 votes vote down vote up
def smooth_spectra(xarr, farr, sigma=3, nkern=20):
    """Given a xarr and flux, smooth the spectrum"""
    xkern = np.arange(nkern)
    kern = np.exp(-(xkern - 0.5 * nkern) ** 2 / (sigma) ** 2)

    return gaussian_filter1d(farr, sigma) 
Example #26
Source File: util.py    From aitom with GNU General Public License v3.0 5 votes vote down vote up
def cluster_formation_alignment_fsc__by_global_maximum(self, dj, op=None):
    if ('debug' not in op):
        op['debug'] = False
    dj = copy.deepcopy(dj)
    djm = defaultdict(list)
    for r in dj:
        if ('template' not in r):
            continue
        djm[str(r['template']['subtomogram'])].append(r)
    djm_org = copy.deepcopy(djm)
    for k in djm:
        djmt = djm[k]
        djmt = sorted(djmt, key=(lambda _: float(_['score'])), reverse=True)
        if (('max_expansion_size' in op) and (len(djmt) > op['max_expansion_size'])):
            djmt = djmt[:op['max_expansion_size']]
        djm[k] = djmt
    ssnr_sequential_op = copy.deepcopy(op['ssnr_sequential'])
    ssnr_sequential_op['n_chunk'] = op['n_chunk']
    ssnr_s = SS.ssnr_sequential_parallel(self=self, data_json_dict=djm, op=ssnr_sequential_op)
    fsc_sum = {}
    for k in ssnr_s:
        fsc_sum[k] = N.array([N.sum(_) for _ in ssnr_s[k]['fsc']])
    import scipy.ndimage.filters as SDF
    if ('gaussian_smooth_sigma' in op):
        for k in fsc_sum:
            fsc_sum[k] = SDF.gaussian_filter1d(fsc_sum[k], op['gaussian_smooth_sigma'])
    if ('min_expansion_size' in op):
        for k in copy.deepcopy(list(fsc_sum.keys())):
            if (len(fsc_sum[k]) < op['min_expansion_size']):
                del fsc_sum[k]
                continue
            fsc_sum[k][:op['min_expansion_size']] = (N.min(fsc_sum[k]) - 1)
    dj_gm = {}
    for k in fsc_sum:
        i = N.argmax(fsc_sum[k])
        if op['debug']:
            print('template', k, 'original subtomogram num', len(djm_org[k]), 'global maximum', i)
        dj_gm[k] = {'k': k, 'i': i, 'data_json': copy.deepcopy(djm[k][:(i + 1)]), 'fsc': ssnr_s[k]['fsc'][i], 'fsc_sum': fsc_sum[k][i], }
    return {'dj_gm': dj_gm, 'djm': djm, 'ssnr_s': ssnr_s, } 
Example #27
Source File: pid_analysis.py    From flight_review with BSD 3-Clause "New" or "Revised" License 5 votes vote down vote up
def wiener_deconvolution(self, input, output, cutfreq):      # input/output are two-dimensional
        pad = 1024 - (len(input[0]) % 1024)                     # padding to power of 2, increases transform speed
        input = np.pad(input, [[0,0],[0,pad]], mode='constant')
        output = np.pad(output, [[0, 0], [0, pad]], mode='constant')
        H = np.fft.fft(input, axis=-1)
        G = np.fft.fft(output,axis=-1)
        freq = np.abs(np.fft.fftfreq(len(input[0]), self.dt))
        sn = self.to_mask(np.clip(np.abs(freq), cutfreq-1e-9, cutfreq))
        len_lpf=np.sum(np.ones_like(sn)-sn)
        sn=self.to_mask(gaussian_filter1d(sn,len_lpf/6.))
        sn= 10.*(-sn+1.+1e-9)       # +1e-9 to prohibit 0/0 situations
        Hcon = np.conj(H)
        deconvolved_sm = np.real(np.fft.ifft(G * Hcon / (H * Hcon + 1./sn),axis=-1))
        return deconvolved_sm 
Example #28
Source File: pid_analysis.py    From flight_review with BSD 3-Clause "New" or "Revised" License 5 votes vote down vote up
def weighted_mode_avr(self, values, weights, vertrange, vertbins):
        ### finds the most common trace and std
        threshold = 0.5  # threshold for std calculation
        filt_width = 7  # width of gaussian smoothing for hist data

        resp_y = np.linspace(vertrange[0], vertrange[-1], vertbins, dtype=np.float64)
        times = np.repeat(np.array([self.time_resp],dtype=np.float64), len(values), axis=0)
        weights = np.repeat(weights, len(values[0]))

        hist2d = np.histogram2d(times.flatten(), values.flatten(),
                                range=[[self.time_resp[0], self.time_resp[-1]], vertrange],
                                bins=[len(times[0]), vertbins], weights=weights.flatten())[0].transpose()
        ### shift outer edges by +-1e-5 (10us) bacause of dtype32. Otherwise different precisions lead to artefacting.
        ### solution to this --> somethings strage here. In outer most edges some bins are doubled, some are empty.
        ### Hence sometimes produces "divide by 0 error" in "/=" operation.

        if hist2d.sum():
            hist2d_sm = gaussian_filter1d(hist2d, filt_width, axis=0, mode='constant')
            hist2d_sm /= np.max(hist2d_sm, 0)


            pixelpos = np.repeat(resp_y.reshape(len(resp_y), 1), len(times[0]), axis=1)
            avr = np.average(pixelpos, 0, weights=hist2d_sm * hist2d_sm)
        else:
            hist2d_sm = hist2d
            avr = np.zeros_like(self.time_resp)
        # only used for monochrome error width
        hist2d[hist2d <= threshold] = 0.
        hist2d[hist2d > threshold] = 0.5 / (vertbins / (vertrange[-1] - vertrange[0]))

        std = np.sum(hist2d, 0)

        return avr, std, [self.time_resp, resp_y, hist2d_sm]

    ### calculates weighted avverage and resulting errors 
Example #29
Source File: data_stabilized_local_speed.py    From ronin with GNU General Public License v3.0 5 votes vote down vote up
def __init__(self, seq_type, root_dir, data_list, step_size=10, window_size=200,
                 random_shift=0, transform=None, **kwargs):
        super(StabilizedLocalSpeedDataset, self).__init__()
        self.shift = 0
        assert seq_type.aux_dim == 11
        self.feature_dim = seq_type.feature_dim
        self.aux_dim = seq_type.aux_dim
        self.target_dim = 2
        self.window_size = window_size
        self.step_size = step_size
        self.random_shift = random_shift
        self.transform = transform

        self.data_path = [osp.join(root_dir, data) for data in data_list]

        self.ts, self.features, self.targets, self.gravity, self.orientations, self.gt_pos = [], [], [], [], [], []
        self.index_map = []

        feature_sigma = kwargs.get('feature_sigma', None)
        target_sigma = kwargs.get('target_sigma', None)

        for i in range(len(data_list)):
            seq = seq_type(osp.join(root_dir, data_list[i]), interval=1, **kwargs)
            self.features.append(seq.get_feature())
            self.targets.append(seq.get_target())
            aux = seq.get_aux()
            print(seq.get_meta())
            self.ts.append(aux[:, 0])
            self.orientations.append(aux[:, 1:5])
            self.gravity.append(aux[:, 5:8])
            self.gt_pos.append(aux[:, -3:])

            if feature_sigma is not None:
                self.features[i] = gaussian_filter1d(self.features[i], sigma=feature_sigma, axis=0)
            if target_sigma is not None:
                self.targets[i] = gaussian_filter1d(self.targets[i], sigma=target_sigma, axis=0)
            self.index_map += [[i, j] for j in range(window_size, self.targets[i].shape[0], step_size)]

        if kwargs.get('shuffle', True):
            random.shuffle(self.index_map) 
Example #30
Source File: im_text_rnn_model.py    From tumblr-emotions with Apache License 2.0 5 votes vote down vote up
def blur_image(np_image, sigma=1):
    np_image = gaussian_filter1d(np_image, sigma, axis=1)
    np_image = gaussian_filter1d(np_image, sigma, axis=2)
    return np_image