jjhelmus / nmrglue

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

Problem number of points in 2D experiments #168

Closed carlosbornes closed 2 years ago

carlosbornes commented 2 years ago

Hi I'm trying to get a plot some NMR data using nmrglue and matplotlib and I found that when trying to plot the 2D spectra and the respective F1 and F2 projections I get the following error

x and y must have same first dimension, but have shapes (2048,) and (1024,)

I was able to troubleshoot this in the following way 1st - Import the data and calculate the projections using

raw_2d = os.path.normpath('path to directory')
dic_2d, data_2d = ng.bruker.read_pdata(raw_2d)
udic_2d = ng.bruker.guess_udic(dic_2d, data_2d)

C = ng.convert.converter()
C.from_bruker(dic_2d, data_2d, udic_2d)
ng.pipe.write("data_2d.fid", *C.to_pipe(), overwrite=True)

dic_pipe, data_pipe = ng.pipe.read('data_2d.fid')
uc_13C = ng.pipe.make_uc(dic_pipe, data_pipe, dim=1)
ppm_13C = uc_13C.ppm_scale()
ppm_13C_0, ppm_13C_1 = uc_13C.ppm_limits()
uc_1H = ng.pipe.make_uc(dic_pipe, data_pipe, dim=0)
ppm_1H = uc_1H.ppm_scale()
ppm_1H_0, ppm_1H_1 = uc_1H.ppm_limits()

# Calculate projections
f1proj = np.amax(data_pipe, axis=1)
f2proj = np.amax(data_pipe, axis=0)

2nd - Go to topspin reprocess spectra with half SI (Size of real spectrum) of F1 dimension and run

raw_2d = os.path.normpath('path to directory')
dic_2d, data_2d = ng.bruker.read_pdata(raw_2d)
udic_2d = ng.bruker.guess_udic(dic_2d, data_2d)

C = ng.convert.converter()
C.from_bruker(dic_2d, data_2d, udic_2d)
ng.pipe.write("data_2d.fid", *C.to_pipe(), overwrite=True)

dic_pipe, data_pipe = ng.pipe.read('data_2d.fid')
uc_13C = ng.pipe.make_uc(dic_pipe, data_pipe, dim=1)
ppm_13C = uc_13C.ppm_scale()
ppm_13C_0, ppm_13C_1 = uc_13C.ppm_limits()
uc_1H = ng.pipe.make_uc(dic_pipe, data_pipe, dim=0)
ppm_1H = uc_1H.ppm_scale()
ppm_1H_0, ppm_1H_1 = uc_1H.ppm_limits()

3rd - Plot everything

cmap = cm.Blues_r 
contour_start = 100000     
contour_num = 30 
contour_factor = 1.4
cl = contour_start * contour_factor ** np.arange(contour_num)

ax = plt.figure(constrained_layout= False).subplot_mosaic(
    """
    .a
    bA
    """,
    gridspec_kw={"height_ratios": [1, 3.5], "width_ratios": [1, 3.5], 'wspace' : 0.02, 'hspace' : 0.03},
)

ax['A'].contour(data_pipe, cl, colors='k', extent=(ppm_13C_0, ppm_13C_1, ppm_1H_0, ppm_1H_1))
ax['A'].yaxis.tick_right()
ax['a'].plot(ppm_13C, f2proj)
ax['b'].plot(-f1proj, ppm_1H)

If I do this 3 steps process with reprocessing the data with half SI in the second step the dimensions seem to see correct otherwise I get the error mentioned above

kaustubhmote commented 2 years ago

I am assuming that the error is with the plotting function, where one of these two lines is failing:

ax['a'].plot(ppm_13C, f2proj)
ax['b'].plot(-f1proj, ppm_1H)

My best guess is that the udic that is being generated is not consistent with the data shape. Are the shapes of ppm_13C and f2_proj (or f1proj and ppm_1H) the same? If not, then there is some inconsistency while reading in the parameters and generating the universal dictionary. one place this could occur is if you use non-zero STSR and STSI values for processing in Topspin. If so, you should indicate that when generating a universal dictionary.

udic_2d = ng.bruker.guess_udic(dic_2d, data_2d, strip_fake=True)

Also, as far as I can see, this should also work without the intermediate conversion to nmrpipe formats (unless you need that for something else). You can generate the unit conversion object from a universal dictionary and data using:

uc_13C = ng.fileiobase.uc_from_udic(udic_2d, dim=1)
uc_1H = ng.fileiobase.uc_from_udic(udic_2d, dim=0)
carlosbornes commented 2 years ago

Everything seems fine with the F2 dimension since the shape is the same

f2proj.shape
(2048,)
ppm_13C.shape
(2048,)

The F1, the calculated projection has double the number of points

f1proj.shape
(256,)
ppm_1H.shape
(128,)

I've tried to use

udic_2d = ng.bruker.guess_udic(dic_2d, data_2d, strip_fake=True)

but the problem persists. But anyway both STSR and STSI are zero.

kaustubhmote commented 2 years ago

It is unlikely that the amax function is giving the incorrect array size. So my guess is that the problem should be with the unit conversion. Can you confirm that the shape of data_pipe is incompatible with the shape of ppm_1H.

carlosbornes commented 2 years ago

I believe it is. I was checking the dimension and

data_pipe.shape
(256, 2048)

But then when I run

uc_1H = ng.pipe.make_uc(dic_pipe, data_pipe, dim=0)
ppm_1H = uc_1H.ppm_scale()

The size is cut in half

ppm_1H.shape
(128,)
kaustubhmote commented 2 years ago

Does skipping the conversion back and forth from nmrpipe, and generating the unit conversion object from bruker data (as shown in this comment) give the same results? If yes, will it be possible to share the data for debugging this?

carlosbornes commented 2 years ago

Not converting to NMRpipe solves the problem. For now its ok but maybe it would be worth it to solve this.