Python signal.find_peaks_cwt方法代碼示例

# 需要導入模塊: from scipy import signal [as 別名]
# 或者: from scipy.signal import find_peaks_cwt [as 別名]
def detect_lines(x, y, kernal_size=3, centroid_kernal=default_kernal, 
                 center=False):
    """Detect lines goes through a 1-D spectra and detect peaks

    Parameters
    ----------
    x: ~numpy.ndarray
        Array describing the x-position 

    y: ~numpy.ndarray
        Array describing the counts in each x-position

    kernal_size: int
        Size for the detection kernal

    centroid_kernal: 
        Kernal to be used for centroiding

    center: boolean
        If True, centroid for detected peaks will be calculated

    Returns
    -------
    xp: ~numpy.ndarray
        Array of x-positions of peaks in the spectrum
    """
    # find all peaks
    xp = signal.find_peaks_cwt(y, np.array([kernal_size]))
    xp = np.array(xp)

    # set the output values
    if center:
        xdiff = int(0.5 * len(centroid_kernal) + 1)
        x_arr = np.arange(len(x))
        xp = xp * 1.0
        for i in range(len(xp)):
            xp[i] = mcentroid(x, y, kern=centroid_kernal, xdiff=xdiff, xc=x[xp[i]])

    return xp 

# 需要導入模塊: from scipy import signal [as 別名]
# 或者: from scipy.signal import find_peaks_cwt [as 別名]
def get_frequencies_from_correlation(correlation_vector,test_frequencies_range):

    frequencies = []
    for branch in range(correlation_vector.shape[1]):
        peakind = signal.find_peaks_cwt(correlation_vector[:,branch].real, np.arange(1,200) )

     #   plt.plot(test_frequencies_range,correlation_vector[:,branch].real)
     #   plt.plot([test_frequencies_range[i] for i in peakind],[correlation_vector[i,branch].real for i in peakind],'ro')
     #   plt.show()

        heights = [correlation_vector[i,branch] for i in peakind]
        max_height_index = heights.index(max(heights))
        frequencies.append(test_frequencies_range[peakind[max_height_index]])

    return np.array(frequencies) 

# 需要導入模塊: from scipy import signal [as 別名]
# 或者: from scipy.signal import find_peaks_cwt [as 別名]
def get_peaks(x, n):
    return find_peaks_cwt(x, widths=np.arange(1, n + 1), wavelet=ricker) 

# 需要導入模塊: from scipy import signal [as 別名]
# 或者: from scipy.signal import find_peaks_cwt [as 別名]
def test_single_continuous_property(self):
        # Tested on a water molecule
        system = ase.build.molecule("H2O")

        # Descriptor setup
        std = 0.1
        elements = ["H", "O"]
        peaks = [0.3, 2.0]
        values = dict(zip(elements, peaks))
        elemdist = ElementalDistribution(
            properties={
                "first_property": {
                    "type": "continuous",
                    "min": 0,
                    "max": 2.5,
                    "std": std,
                    "n": 50,
                    "values": values
                }
            }
        )

        # Features
        y = elemdist.create(system)
        y = y.todense().A1
        x = elemdist.get_axis("first_property")

        # Test that the peak positions match
        from scipy.signal import find_peaks_cwt
        peak_indices = find_peaks_cwt(y, [std])
        peak_loc = x[peak_indices]

        # Test that the peak locations match within some tolerance
        self.assertTrue(np.allclose(peaks, peak_loc, rtol=0, atol=0.05))

        # Plot for visual inspection
        # mpl.plot(x, y)
        # mpl.show() 

# 需要導入模塊: from scipy import signal [as 別名]
# 或者: from scipy.signal import find_peaks_cwt [as 別名]
def fragment_length_qc(data, prefix):
    results = []

    NFR_UPPER_LIMIT = 150
    MONO_NUC_LOWER_LIMIT = 150
    MONO_NUC_UPPER_LIMIT = 300

    # % of NFR vs res
    nfr_reads = data[data[:, 0] < NFR_UPPER_LIMIT][:, 1]
    percent_nfr = nfr_reads.sum() / data[:, 1].sum()
    results.append(
        QCGreaterThanEqualCheck('Fraction of reads in NFR', 0.4)(percent_nfr))

    # % of NFR vs mononucleosome
    mono_nuc_reads = data[
        (data[:, 0] > MONO_NUC_LOWER_LIMIT) &
        (data[:, 0] <= MONO_NUC_UPPER_LIMIT)][:, 1]

    percent_nfr_vs_mono_nuc = (
        nfr_reads.sum() /
        mono_nuc_reads.sum())
    results.append(
        QCGreaterThanEqualCheck('NFR / mono-nuc reads', 2.5)(
            percent_nfr_vs_mono_nuc))

    # peak locations
    pos_start_val = data[0, 0]  # this may be greater than 0
    peaks = find_peaks_cwt(data[:, 1], np.array([25]))
    nuc_range_metrics = [
        ('Presence of NFR peak', 20 - pos_start_val, 90 - pos_start_val),
        ('Presence of Mono-Nuc peak',
         120 - pos_start_val, 250 - pos_start_val),
        ('Presence of Di-Nuc peak',
         300 - pos_start_val, 500 - pos_start_val)]
    for range_metric in nuc_range_metrics:
        results.append(QCHasElementInRange(*range_metric)(peaks))

    out = prefix + '.nucleosomal.qc'
    with open(out, 'w') as fp:
        for elem in results:
            fp.write(
                '\t'.join(
                    [elem.metric, str(elem.qc_pass), elem.message]) + '\n')

    return out