jjhelmus / nmrglue

A module for working with NMR data in Python
BSD 3-Clause "New" or "Revised" License
208 stars 84 forks source link

Peak picking from a 2D contour plot #217

Open coffeedealer opened 1 month ago

coffeedealer commented 1 month ago

Hey there,

I used the Help plotting Bruker 2D spectra #214 in order to be able to plot the 2D contour plot from NMR experimental data from Bruker. I succeed, but next step I was planning to detect all peaks from that 2D plot using the peak peaking as described in here.

The final code used is here below:

# read data
dic, data = ng.bruker.read_pdata(".../test/11/pdata/1/")

# read in axis ppm parameters
udic = ng.bruker.guess_udic(dic, data, strip_fake=False)
uc = {i: ng.fileiobase.uc_from_udic(udic, dim=i) for i in (0, 1)}
axis_limits = [*uc[1].ppm_limits(), *uc[0].ppm_limits()]

####################
#This part is the one implemented for the peak picking detection from the second source
####################

# detect all peaks with a threshold
peaks = ng.peakpick.pick(data, pthres=0.6, algorithm="thres", msep=[0.1, 0.1])

# add markers for peak positions
x = uc[1].ppm(peaks["X_AXIS"])
y = uc[0].ppm(peaks["Y_AXIS"])

####################

# contour levels
levels = [data.max() / 6.0 * 1.4**i for i in range(10)]  # change as needed

# plot
#fig, ax = plt.subplots(figsize=(9,14))
ax.contour(
    data.real, 
    levels=levels, 
    extent=axis_limits, 
    colors="#0092c8", # change as per your preference 
    linewidths=1.2, # change as per your preference
)

####################
ax.scatter(x, y, marker=".", s=10, color="r")
####################

plt.xlabel('1H Chemical Shift [ppm]', fontsize=14)
plt.ylabel('13C Chemical Shift [ppm]', fontsize=14)
ax.set_xlim(10.0, -2.5)  # change as needed
ax.set_ylim(140, 0)  # change as needed
ax.yaxis.tick_right()
ax.yaxis.set_label_position("right")

However, it does not work and the error I got is the following:

IndexError                                Traceback (most recent call last)
Cell In[4], line 27
     24 dic, data = ng.bruker.read_pdata(".../test/11/pdata/1/")
     26 # detect all peaks with a threshold
---> 27 peaks = ng.peakpick.pick(data, pthres=0.6, algorithm="thres", msep=[0.1, 0.1])
     29 # read in axis ppm parameters
     30 udic = ng.bruker.guess_udic(dic, data, strip_fake=False)

File ~/anaconda3/lib/python3.11/site-packages/nmrglue/analysis/peakpick.py:214, in pick(data, pthres, nthres, msep, algorithm, est_params, lineshapes, edge, diag, c_struc, c_ndil, cluster, table, axis_names)
    210 # scales = np.zeros(np.array(locations).shape,dtype=float)
    211 # amps = np.zeros(len(locations),dtype=float)
    213 for i, (l, seg_slice) in enumerate(zip(locations, seg_slices)):
--> 214     null, scales[i], amps[i] = guess_params_slice(data, l, seg_slice,
    215                                                   ls_classes)
    217 ########################################################
    218 # return locations, scales and amplitudes as requested #
    219 ########################################################
    220 if cluster:

File ~/anaconda3/lib/python3.11/site-packages/nmrglue/analysis/peakpick.py:370, in guess_params_slice(data, location, seg_slice, ls_classes)
    345 """
    346 Guess the parameter of a peak in a segment.
    347 
   (...)
    367 
    368 """
    369 # find the rectangular region around the segment
--> 370 region = data[seg_slice]
    371 edge = [s.start for s in seg_slice]
    372 rlocation = [l - s.start for l, s in zip(location, seg_slice)]

IndexError: only integers, slices (`:`), ellipsis (`...`), numpy.newaxis (`None`) and integer or boolean arrays are valid indices

Any suggestion? Thank you in advance.

kaustubhmote commented 1 month ago

Hopefully #216 has fixed this issue. This change is not in the the version of nmrglue published on pypi. You can install the latest branch on github with python -m pip install git+git://github.com/jjhelmus/nmrglue.