API

NMRBase

NMRTools.NMRBase.AbstractNMRDataType
AbstractNMRData <: DimensionalData.AbstractDimArray

Abstract supertype for objects that wrap an array of NMR data, and metadata about its contents.

AbstractNMRDatas inherit from AbstractDimArray from DimensionalData.jl. They can be indexed as regular Julia arrays or with DimensionalData.jl Dimensions.

source
NMRTools.NMRBase.GaussWindowType
GaussWindow(expHz, gaussHz, center, tmax)

Abstract representation of Lorentz-to-Gauss window functions, applying an inverse exponential of expHz Hz, and a gaussian broadening of gaussHz Hz, with maximum at center (between 0 and 1). Acquisition time is tmax.

Specialises to LorentzToGaussWindow when center is zero, otherwise GeneralGaussWindow.

source
NMRTools.NMRBase.GradientDimensionType
GradientDimension <: NonFrequencyDimension <: NMRDimension

Abstract supertype for gradient-encoded dimensions used in NMRData objects. Concrete types G1Dim, G2Dim, G3Dim and G4Dim are generated for use in creating objects.

source
NMRTools.NMRBase.HasNonFrequencyDimensionType
@traitdef HasNonFrequencyDimension{D}

A trait indicating whether the data object has a non-frequency domain dimension.

Example

@traitfn f(x::X) where {X; HasNonFrequencyDimension{X}} = "This spectrum has a non-frequency domain dimension!"
@traitfn f(x::X) where {X; SimpleTraits.Not{HasNonFrequencyDimension{X}}} = "This is a pure frequency-domain spectrum!"
source
NMRTools.NMRBase.MQType
MQ(coherences, label=="")

Representation of a multiple-quantum coherence. Coherences are specified as a tuple of tuples, of the form (nucleus, coherenceorder)

Examples

julia> MQ(((H1,1), (C13,-1)), "ZQ")
MQ(((H1, 1), (C13, -1)), "ZQ")

julia> MQ(((H1,3), (C13,1)), "QQ")
MQ(((H1, 3), (C13, 1)), "QQ")

See also Nucleus, SQ.

source
NMRTools.NMRBase.NMRDataType
NMRData <: AbstractNMRData
NMRData(A::AbstractArray{T,N}, dims; kw...)
NMRData(A::AbstractNMRData; kw...)

A generic AbstractNMRData for NMR array data. It holds memory-backed arrays.

Keywords

  • dims: Tuple of NMRDimensions for the array.
  • name: Symbol name for the array, which will also retreive named layers if NMRData is used on a multi-layered file like a NetCDF.
  • missingval: value representing missing data, normally detected from the file. Set manually when you know the value is not specified or is incorrect. This will not change any values in the NMRData, it simply assigns which value is treated as missing.

Internal Keywords

In some cases it is possible to set these keywords as well.

  • data: can replace the data in an AbstractNMRData

  • refdims: Tuple of position Dimensions the array was sliced from, defaulting to ().

source
NMRTools.NMRBase.PowerType
Power{T}

A type to represent NMR pulse powers that can be described either in Watts (W) or as dB attenuation (dB). Internally stores the dB value as type T.

Constructors

  • Power(value, :dB) - Create from dB attenuation value
  • Power(value, :W) - Create from Watts value

Getters

  • db(p::Power) - Get dB attenuation value
  • watts(p::Power) - Get Watts value

Examples

# Create power from dB attenuation
p1 = Power(30.0, :dB)
db(p1)    # 30.0
watts(p1) # 0.001

# Create power from Watts
p2 = Power(0.5, :W)
watts(p2) # 0.5
db(p2)    # ≈ 3.01

# Works with integers and other numeric types
p3 = Power(10, :W)
db(p3)    # -10.0

# Zero watts is represented as 120 dB attenuation
p4 = Power(0.0, :W)
db(p4)    # 120.0

Conversion

Converting W to dB: -10 * log10(plW) Converting dB to W: 10^(-dB/10)

source
NMRTools.NMRBase.SineWindowType
SineWindow(offset, endpoint, power, tmax)

Abstract window function representing multiplication by sine/cosine functions. Acquisition time is tmax.

\[\left[\sin\left( \pi\cdot\mathrm{offset} + \frac{\left(\mathrm{end} - \mathrm{offset}\right)\pi t}{\mathrm{tmax}} \right)\right]^\mathrm{power}\]

Specialises to CosWindow, Cos²Window or GeneralSineWindow.

Arguments

  • offset: initial value is $\sin(\mathrm{offset}\cdot\pi)$ (0 to 1)
  • endpoint: initial value is $\sin(\mathrm{endpoint}\cdot\pi)$ (0 to 1)
  • pow: sine exponent
source
NMRTools.NMRBase.TimeDimensionType
TimeDimension <: NonFrequencyDimension <: NMRDimension

Abstract supertype for time dimensions used in NMRData objects. Concrete types T1Dim, T2Dim, T3Dim and T4Dim are generated for time-domains representing frequency evolution, and TrelaxDim and TkinDim are generated for representing relaxation and real-time kinetics.

source
NMRTools.NMRBase.UnknownDimensionType
UnknownDimension <: NonFrequencyDimension <: NMRDimension

Abstract supertype for unknown, non-frequency dimensions used in NMRData objects. Concrete types X1Dim, X2Dim, X3Dim and X4Dim are generated for use in creating objects.

source
NMRTools.NMRBase.WindowFunctionType
WindowFunction

Abstract type to represent apodization functions.

Window functions are represented by subtypes of the abstract type WindowFunction, each of which contain appropriate parameters to specify the particular function applied. In addition, the acquisition time tmax is also stored (calculated at the point the window function is applied, i.e. after linear prediction but before zero filling).

source
Base.stackMethod
stack(expts::Vector{NMRData})

Combine a collection of equally-sized NMRData into one larger array, by arranging them along a new dimension, of type UnknownDimension.

Throws a DimensionMismatch if data are not of compatible shapes.

source
DimensionalData.Dimensions.Lookups.metadataFunction
metadata(nmrdata, key)
metadata(nmrdata, dim, key)
metadata(nmrdimension, key)

Return the metadata for specified key, or nothing if not found. Keys are passed as symbols.

Examples (spectrum metadata)

  • :ns: number of scans
  • :ds: number of dummy scans
  • :rg: receiver gain
  • :ndim: number of dimensions
  • :title: spectrum title (contents of title pdata file)
  • :filename: spectrum filename
  • :pulseprogram: title of pulse program used for acquisition
  • :experimentfolder: path to experiment
  • :noise: RMS noise level

Examples (dimension metadata)

  • :pseudodim: flag indicating non-frequency domain data
  • :npoints: final number of (real) data points in dimension (after extraction)
  • :td: number of complex points acquired
  • :tdzf: number of complex points when FT executed, including LP and ZF
  • :bf: base frequency, in Hz
  • :sf: carrier frequency, in Hz
  • :offsethz: carrier offset from bf, in Hz
  • :offsetppm: carrier offset from bf, in ppm
  • :swhz: spectrum width, in Hz
  • :swppm: spectrum width, in ppm
  • :region: extracted region, expressed as a range in points, otherwise missing
  • :window: WindowFunction indicating applied apodization
  • :referenceoffset: referencing (in ppm) applied to the dimension

See also estimatenoise!.

source
DimensionalData.Dimensions.labelFunction
label(nmrdata)
label(nmrdata, dim)
label(nmrdimension)

Return a short label associated with an NMRData structure or an NMRDimension. By default, for a spectrum this is obtained from the first line of the title file. For a frequency dimension, this is normally something of the form 1H chemical shift (ppm).

See also label!.

source
NMRTools.NMRBase._lineshapeFunction
_lineshape(ω, R2, ωaxis, window, complexity)

Internal function to calculate a resonance lineshape with frequency ω and relaxation rate R2, calculated at frequencies ωaxis and with apodization according to the specified window function.

source
NMRTools.NMRBase.acqusMethod
acqus(nmrdata)
acqus(nmrdata, key)
acqus(nmrdata, key, index)

Return data from a Bruker acqus file, or nothing if it does not exist. Keys can be passed as symbols or strings. If no key is specified, a dictionary is returned representing the entire acqus file.

If present, the contents of auxilliary files such as vclist and vdlist can be accessed using this function.

Examples

julia> acqus(expt, :pulprog)
"zgesgp"
julia> acqus(expt, "TE")
276.9988
julia> acqus(expt, :p, 1)
9.2
julia> acqus(expt, "D", 1)
0.1
julia> acqus(expt, :vclist)
11-element Vector{Int64}:
[...]

See also metadata.

source
NMRTools.NMRBase.annotationsMethod
annotations(nmrdata)
annotations(nmrdata, key)
annotations(nmrdata, key, index)
annotations(nmrdata, key1, key2)
annotations(nmrdata, keys...)

Return annotation data from pulse programme, or nothing if it does not exist. Keys can be passed as strings or symbols. If no key is specified, the entire annotations dictionary is returned.

Dotted names are automatically split and dereferenced, e.g. annotations(spec, "cest.duration") is equivalent to annotations(spec, "cest", "duration").

This function provides nested access to annotations stored in metadata[:annotations]. Multiple keys can be chained to access nested dictionaries or specific array elements.

Examples

julia> annotations(spec, "title")
"19F CEST"

julia> annotations(spec, "dimensions")
2-element Vector{String}:
 "cest.offset"
 "f1"

julia> annotations(spec, "dimensions", 1)
"cest.offset"

julia> annotations(spec, "cest", "duration")
0.5

julia> annotations(spec, "cest.duration")
0.5

julia> annotations(spec, :reference_pulse)
1-element Vector{Dict{String, Any}}:
 Dict("channel" => "f1", "pulse" => 9.2e-6, "power" => -3.0)

julia> annotations(spec, "calibration.duration.start")
0.001

See also metadata, acqus.

source
NMRTools.NMRBase.apodFunction
apod(spec::NMRData, dimension, zerofill=true)

Return the time-domain apodization function for the specified axis, as a vector of values.

source
NMRTools.NMRBase.coherenceorderFunction
coherenceorder(coherence)

Calculate the total coherence order.

Examples

julia> coherenceorder(SQ(H1))
1

julia> coherenceorder(MQ(((H1,1),(C13,1))))
2

julia> coherenceorder(MQ(((H1,1),(C13,-1))))
0

julia> coherenceorder(MQ(((H1,3),(C13,1))))
4

julia> coherenceorder(MQ(((H1,0),)))
0

See also Nucleus, SQ, MQ.

source
NMRTools.NMRBase.dbMethod
db(p::Power)

Get the power level in dB attenuation.

Example

p = Power(20.0, :dB)
db(p)  # 20.0
source
NMRTools.NMRBase.decimateFunction
decimate(data, n dims=1)

Decimate NMR data into n-point averages along the specified dimension. Note that data are averaged and not summed. Noise metadata is not updated.

source
NMRTools.NMRBase.estimatenoise!Method
estimatenoise!(nmrdata)

Estimate the rms noise level in the data and update :noise metadata.

If called on an Array of data, each item will be updated.

Algorithm

Data are sorted into numerical order, and the highest and lowest 12.5% of data are discarded (so that 75% of the data remain). These values are then fitted to a truncated gaussian distribution via maximum likelihood analysis.

The log-likelihood function is:

\[\log L(\mu, \sigma) = \sum_i{\log P(y_i, \mu, \sigma)}\]

where the likelihood of an individual data point is:

\[\log P(y,\mu,\sigma) = \log\frac{ \phi\left(\frac{x-\mu}{\sigma}\right) }{ \sigma \cdot \left[\Phi\left(\frac{b-\mu}{\sigma}\right) - \Phi\left(\frac{a-\mu}{\sigma}\right)\right]}\]

and $\phi(x)$ and $\Phi(x)$ are the standard normal pdf and cdf functions.

source
NMRTools.NMRBase.gyromagneticratioFunction
gyromagneticratio(n::Nucleus)
gyromagneticratio(c::Coherence)

Return the gyromagnetic ratio in Hz/T of a nucleus, or calculate the effective gyromagnetic ratio of a coherence. This is equal to the product of the individual gyromagnetic ratios with their coherence orders.

Returns nothing if not defined.

Examples

julia> gyromagneticratio(H1)
2.6752218744e8

julia> gyromagneticratio(SQ(H1))
2.6752218744e8

julia> gyromagneticratio(MQ(((H1,1),(C13,1))))
3.3480498744e8

julia> gyromagneticratio(MQ(((H1,0),)))
0.0

See also Nucleus, Coherence.

source
NMRTools.NMRBase.hasnonfrequencydimensionMethod
hasnonfrequencydimension(spectrum)

Return true if the spectrum contains a non-frequency domain dimension.

Example

julia> y2=loadnmr("exampledata/2D_HN/test.ft2");
julia> hasnonfrequencydimension(y2)
false
julia> y3=loadnmr("exampledata/pseudo3D_HN_R2/ft/test%03d.ft2");
julia> hasnonfrequencydimension(y3)
true
source
NMRTools.NMRBase.hzMethod
hz(p::Power, ref_p::Power, ref_pulselength, ref_pulseangle_deg)

Convert power to radiofrequency strength in Hz using reference pulse parameters.

Arguments

  • p::Power: Power to convert
  • ref_p::Power: Reference power
  • ref_pulselength: Reference pulse length in microseconds
  • ref_pulseangle_deg: Reference pulse flip angle in degrees

Returns

Radiofrequency strength in Hz for power p

source
NMRTools.NMRBase.hzMethod
hz(p::Power, ref_p::Power, ref_Hz)

Convert power to radiofrequency strength in Hz using a reference power and known Hz value.

Uses the relationship: Hz = ref_Hz * 10^(-ΔdB/20) where ΔdB is the power difference.

Arguments

  • p::Power: Power to convert
  • ref_p::Power: Reference power with known Hz value
  • ref_Hz: Radiofrequency strength in Hz at the reference power

Returns

Radiofrequency strength in Hz for power p

source
NMRTools.NMRBase.hzMethod
hz(δ, axis)

Return the offset (in Hz) for a chemical shift (or list of shifts) on a frequency axis.

source
NMRTools.NMRBase.label!Function
label!(nmrdata, labeltext)
label!(nmrdata, dim, labeltext)
label!(nmrdimension, labeltext)

Set the label associated with an NMRData structure or an NMRDimension.

See also label.

source
NMRTools.NMRBase.lineshapeMethod
lineshape(axis, δ, R2, complexity=RealLineshape())

Return a simulated real- or complex-valued spectrum for a resonance with chemical shift δ and relaxation rate R2, using the parameters and window function associated with the specified axis.

source
NMRTools.NMRBase.nucleusMethod
nucleus(label::AbstractString) -> Nucleus

Parse a nucleus from a string label, returning the corresponding Nucleus enum value.

The function accepts common NMR nucleus notation formats:

  • Mass number followed by element symbol: "1H", "13C", "15N", "19F", "31P"
  • Element symbol followed by mass number: "H1", "C13", "N15", "F19", "P31"

Examples

nucleus("19F")  # returns F19
nucleus("1H")   # returns H1
nucleus("13C")  # returns C13
nucleus("15N")  # returns N15

Throws ArgumentError if the nucleus is not recognised or not defined in the Nucleus enum.

source
NMRTools.NMRBase.pclsMethod
pcls(A, y)

Compute the phase-constrained least squares solution:

\[y = A x e^{i\phi}\]

following the algorithm of Bydder (2010) Lin Alg & Apps.

Returns the tuple (x, ϕ), containing the component amplitudes and global phase.

If passed a matrix Y, the function will return a matrix of component amplitudes and a list of global phases corresponding to each row.

Arguments

  • A: (m,n) complex matrix with component spectra
  • y: (m,) complex vector containing the observed spectrum
source
NMRTools.NMRBase.ppmMethod
ppm(offset, axis)

Return the chemical shifts for a given offset or list of offsets along a frequency axis.

source
NMRTools.NMRBase.referencepulseMethod
referencepulse(spec, nucleus) -> (pulse_length, power)

Get the reference pulse calibration for a given nucleus from pulse programme annotations. Returns a tuple of (pulse_length, power) or nothing if not found.

The reference pulse is the calibrated 90° hard pulse used as the basis for all power calculations on that channel.

Arguments

  • spec: NMRData object with annotations
  • nucleus: Nucleus symbol or string (e.g., :19F, "19F", :1H, "1H")

Returns

  • (pulse_length, power): Tuple where pulse_length is in μs and power is a Power object
  • nothing: If no reference pulse found for the specified nucleus

Examples

julia> referencepulse(spec, "19F")
(13.29, Power(-9.03 dB, 8.0 W))

julia> referencepulse(spec, :1H)
(9.2, Power(-3.0 dB, 32.0 W))

See also annotations.

source
NMRTools.NMRBase.sampleMethod
sample(nmrdata)
sample(nmrdata, keys...)

Return the sample metadata associated with the NMRData object, or nothing if not available. Sample metadata follows the schema defined at https://github.com/waudbygroup/nmr-sample-schema.

When keys are provided, navigates nested dictionaries using the given keys. Keys can be strings or symbols and are automatically converted to lowercase strings for consistent access.

Examples

julia> sample(spec)
Dict{String, Any} with 6 entries:
  "nmr_tube" => Dict{String, Any}("diameter"=>"5 mm", "type"=>"regular", "sample_volume_uL"=>600)
  "buffer"   => Dict{String, Any}("solvent"=>"10% D2O", "reference_unit"=>"%w/v", "chemical_shift_reference"=>"none")
  [...]

julia> sample(spec, "notes")
"here's a note!"

julia> sample(spec, :people, :users)
1-element Vector{Any}:
 "Chris"

julia> sample(spec, "sample", "components")
2-element Vector{Any}:
 Dict{String, Any}("name" => "HEWL", "isotopic_labelling" => "unlabelled", "unit" => "mM", "concentration" => 10)
 Dict{String, Any}("name" => "gadodiamide", "isotopic_labelling" => "unlabelled", "unit" => "mM", "concentration" => 1)

See also metadata, hassample.

source
NMRTools.NMRBase.scaleMethod
scale(d::AbstractNMRData)

Return a scaling factor for the data combining the number of scans, receiver gain, and, if specified, the sample concentration.

\[\mathrm{scale} = \mathrm{ns} \cdot \mathrm{rg} \cdot \mathrm{conc}\]

source
NMRTools.NMRBase.setgradientlistMethod
setgradientlist(A::NMRData, [dimnumber], relativegradientlist, Gmax=nothing)

Return a new NMRData with a gradient axis containing the passed values. If no maximum strength is specified, a default gradient strength of 0.55 T m⁻¹ will be set, but a warning raised for the user.

If a dimension number is specified, that dimension will be replaced. If not, the function will search for a unique non-frequency dimension, and replace that. If there are multiple non-frequency dimensions, the dimension number must be specified explicitly, and the function will throw an error.

source
NMRTools.NMRBase.setkinetictimesMethod
setkinetictimes(A::NMRData, [dimnumber], tvals, units="")

Return a new NMRData with a kinetic time axis containing the passed values (and optionally, units). If a dimension number is specified, that dimension will be replaced with a TkinDim. If not, the function will search for a unique non-frequency dimension, and replace that. If there are multiple non-frequency dimensions, the dimension number must be specified explicitly, and the function will throw an error.

source
NMRTools.NMRBase.setrelaxtimesMethod
setrelaxtimes(A::NMRData, [dimnumber], tvals, units="")

Return a new NMRData with a relaxation time axis containing the passed values (and optionally, units). If a dimension number is specified, that dimension will be replaced with a TrelaxDim. If not, the function will search for a unique non-frequency dimension, and replace that. If there are multiple non-frequency dimensions, the dimension number must be specified explicitly, and the function will throw an error.

source
NMRTools.NMRBase.shiftdimMethod
shiftdim!(data::NMRData, dim_ref, offset)

Add an offset to a frequency dimension in an NMRData object. The dimension can be specified as a numerical index or an object like F1Dim. The metadata is copied using replacedimension, and an entry is added or updated in the dimension metadata to record the offset change.

source

NMRIO

NMRTools.NMRIO.FQListType
FQList(values, unit::Symbol, relative::Bool)

Represents a frequency list. unit can be :Hz or :ppm, and relative indicates whether the frequency is given relative to SFO (true) or BF (false).

Raw values can be extracted using the data function, or (better) as absolute chemical shifts (in ppm) or relative offsets (in Hz) using ppm and hz functions.

See also: ppm, hz.

source
NMRTools.NMRBase.hzMethod
hz(f::FQList, ax::FrequencyDimension)

Return frequency list values as offsets relative to the spectrometer frequency, in Hz.

See also: ppm

source
NMRTools.NMRBase.ppmMethod
ppm(f::FQList, ax::FrequencyDimension)

Return frequency list values in ppm (in absolute terms, i.e. relative to 0 ppm).

See also: hz

source
NMRTools.NMRIO._find_dimension_pathFunction
_find_dimension_path(key, annotations, prefix="")

Find the dotted path to a key in the annotations dictionary. Returns the path that should match entries in the dimensions array.

source
NMRTools.NMRIO._is_programmatic_patternMethod
_is_programmatic_pattern(value) -> Bool

Check if a value is a programmatic list pattern. Must be a Dict with either:

  • "type" key and either "step" or "end" key (linear/log patterns)
  • "counter" and "scale" keys (counter-based pattern)
source
NMRTools.NMRIO._resolve_parameterMethod
_resolve_parameter(param_str::String, spec::NMRData) -> Union{Nothing, Any}

Attempt to resolve a parameter string using the acqus data. Returns the resolved value or nothing if not a recognized parameter.

Handles parameter patterns (case-insensitive):

  • Indexed parameters: p1, pl1, d20, cnst8, gpz6, etc.
  • Nucleus assignments: f1, f2, f3, f4
  • Any other parameter names that exist in acqus (including lists)
source
NMRTools.NMRIO.exampledataMethod
exampledata(name::String)

Load example NMR data provided with NMRTools.jl

Arguments

  • name::String : Name of example data set to load. Call exampledata() to list available data.

Returns

  • NMR data object(s) corresponding to the requested example data set.

Example

data = exampledata("1D_19F_titration")
source
NMRTools.NMRIO.getformatMethod
getformat(filename)

Take an input filename and return either :ucsf, :nmrpipe, :pdata (bruker processed), or :unknown after checking whether the filename matches any known format.

source
NMRTools.NMRIO.loadnmrMethod
loadnmr(filename, experimentfolder=nothing)

Main function for loading NMR data. experimentfolder contains the path to an experiment directory, for identification of metadata, if the filename is not directly within an experiment.

Returns an NMRData structure, or throws an NMRToolsError is there is a problem.

Examples

nmrPipe import:

loadnmr("exampledata/1D_1H/1/test.ft1");
loadnmr("exampledata/1D_19F/1/test.ft1");
loadnmr("exampledata/2D_HN/1/test.ft2");
loadnmr("exampledata/pseudo2D_XSTE/1/test.ft1");
loadnmr("exampledata/pseudo3D_HN_R2/1/ft/test%03d.ft2");

Bruker pdata import:

loadnmr("exampledata/1D_19F/1");
loadnmr("exampledata/1D_19F/1/");
loadnmr("exampledata/1D_19F/1/pdata/1");
loadnmr("exampledata/1D_19F/1/pdata/1/");

ucsf (Sparky) import:

loadnmr("exampledata/2D_HN/1/hmqc.ucsf");
loadnmr("exampledata/1D_19F/1/");
loadnmr("exampledata/1D_19F/1/pdata/1");
loadnmr("exampledata/1D_19F/1/pdata/1/");
source
NMRTools.NMRIO.parse_annotationsMethod
parse_annotations(content::String) -> Dict{String, Any}

Parse pulse programme annotations from content string. Annotations are embedded within pulse programme comments marked with ;@ at the start of lines.

The function:

  1. Filters out lines starting with ;@
  2. Strips the ;@ prefix from each line
  3. Parses the resulting YAML document into a dictionary
source
NMRTools.NMRIO.parsefqlistMethod
parsefqlist(lines)

Return contents of the specified fqlist.

fqlists can have several different formats:

first line | reference | unit

         | sfo         | Hz

sfo hz | sfo | Hz sfo ppm | sfo | ppm bf hz | bf | Hz bf ppm | bf | ppm p | sfo | ppm P | bf | ppm

source
NMRTools.NMRIO.parsenmrpipeheaderMethod
parsenmrpipeheader(header)

Pass a 512 x 4 byte array containing an nmrPipe header file, and returns dictionaries of metadata. The nmrPipe header format is defined in fdatap.h.

Return values

  • md: dictionary of spectrum metadata
  • mdax: array of dictionaries containing axis metadata

Examples

md, mdax = parsenmrpipeheader(header)
source
NMRTools.NMRIO.parsepulseprogramlists!Method

Example lines to parse: define list<gradient> EA=<EA> define list<gradient> EA3 = { 1.0000 0.8750 } define list<pulse> taulist = <$VPLIST> define list<power> powerlist = <$VALIST> define list<gradient> diff=<Difframp> define list<frequency> F19sat = <$FQ1LIST>

NB gradient files have weird syntax!!

source
NMRTools.NMRIO.resolve_parameter_references!Method
resolve_parameter_references!(annotations::Dict{String, Any}, spec::NMRData)

Recursively traverse the annotations dictionary and resolve parameter references. References are bare parameter names like p1, pl1, d20, VALIST, etc.

Parameters are resolved using the acqus() function to get values from the acqus file.

Arguments

  • annotations::Dict{String, Any}: Annotations dictionary to modify in-place
  • spec::NMRData: NMRData object containing acqus metadata

Examples

# Before resolution: {"pulse": "p1", "power": "pl1", "duration": "VDLIST"}
# After resolution: {"pulse": 9.2, "power": Power(...), "duration": [0.01, 0.02, ...]}
source
NMRTools.NMRIO.resolve_programmatic_lists!Method
resolve_programmatic_lists!(annotations::Dict{String, Any}, spec::NMRData)

Resolve programmatic list patterns in annotations to actual vectors.

Programmatic lists are specified as dictionaries with one of two patterns:

Pattern 1: Type-based (linear/log)

  • type: "linear" or "log" (logarithmic)
  • start: starting value
  • step: step size (for linear with start/step)
  • end: ending value (for linear or log with start/end)

Pattern 2: Counter-based

  • counter: array of values (typically integers from a loop counter)
  • scale: scale factor to multiply each counter element by

For type-based patterns, the function:

  1. Identifies fields containing programmatic list patterns ({type, start, step/end})
  2. Finds the dimension number for that field from the dimensions annotation
  3. Gets the number of points in that dimension
  4. Replaces the pattern with a correctly sized vector

For counter-based patterns, the function:

  1. Multiplies each element of the counter array by the scale value
  2. Returns the resulting vector

Examples

# Linear spacing with start/step
# {type: linear, start: p9, step: p9} with 10 points → [p9, 2*p9, 3*p9, ...]

# Linear spacing with start/end
# {type: linear, start: 0.001, end: 0.1} with 10 points → [0.001, 0.012, ..., 0.1]

# Logarithmic spacing with start/end
# {type: log, start: 0.001, end: 0.1} with 10 points → [0.001, 0.0016, ..., 0.1]

# Counter-based pattern
# {counter: [0, 1, 2, 3], scale: 0.01} → [0.0, 0.01, 0.02, 0.03]
source
NMRTools.NMRIO.sumexptsMethod
sumexpts(inputexpts...; out, weights=[])

Sum a collection of Bruker NMR experiments, with optional weighting factors. Data will be truncated to fit the shortest data file.

Arguments

  • inputexpts...: Input experiment names (strings or integers)
  • out: Output experiment name (string or integer)
  • weights: Optional vector of weights for each input experiment

Details

  • Experiment names should be either strings with filenames or integers converted to strings
  • If weights are supplied, the length must match the number of input experiments
  • First input experiment is used as a template for the output experiment
  • For 1D experiments, 'fid' files are added; for nD experiments, 'ser' files are added
  • Processed data files in pdata/1 are removed, other pdata subdirectories are deleted
  • Updates the title file with information about the summed experiments
  • Updates the number of scans in the acqus file to be the sum of individual experiments
  • Prompts before overwriting any existing experiment
source
NMRTools.NMRIO.topspinversionMethod
parsetopspinversion(acqusfilename)

Return the TopSpin version as a VersionNumber (e.g. v"4.2.0").

This is obtained from the end of the first line of the acqus file, e.g.

##TITLE= Parameter file, TopSpin 4.2.0
##TITLE= Parameter file, TOPSPIN		Version 2.1
##TITLE= Parameter file, TOPSPIN		Version 3.2
source