mailRe: numerical cpmg fit


Others Months | Index by Date | Thread Index
>>   [Date Prev] [Date Next] [Thread Prev] [Thread Next]

Header


Content

Posted by Paul Schanda on July 16, 2013 - 08:02:
Yes, that's right. As one of the states is assumed to be at zero frequency, the frequency of the other state is the same as the chemical shift difference.

paul


Am 15.07.2013 um 10:56 schrieb "Edward d'Auvergne" <edward@xxxxxxxxxxxxx>:

Hi,

Could it be that the 'fg' parameter in the original code at
http://thread.gmane.org/gmane.science.nmr.relax.devel/4132 is actually
the chemical shift difference?  Is 'fe' assumed to be zero, hence 'fg
= fe + df', where 'df' is the chemical shift difference in rad/second?

Regards,

Edward


On 15 July 2013 10:01, Edward d'Auvergne <edward@xxxxxxxxxxxxx> wrote:
Hi,

fg is the frequency of the ground state (call it "A" or "G"), while fe is the excited state frequency.

Is this the chemical shift in Hz?  Or the chemical shift difference in
Hz?  The code does not use the 'fe' parameter, so how is the chemical
shift of the excited state handled?  Does this mean that chemical
shifts need to be explicitly given?  How does this work with
intermediate to fast exchange where the shift of the peak is a mixture
of the two chemical shifts?

I'm still not sure how to handle this in relax.  I may need to add a
new dispersion parameter(s) to the definitions at the start of the
Relax_disp class in specific_analyses/relax_disp/api.py.  Though very
easy to expand, this currently only supports the following parameters
for chemical shifts:

phi_ex = pA.pB.dw**2 (in ppm^2)
padw2 = pA.dw**2 (in ppm^2)
dw = (in ppm)

where dw is the chemical shift difference.  In all of the analytic
models, only dw is used.


tcpmg is the total duration of the CPMG element. At some instances this is called Tc instead of tcpmg.
Not to be confused with the delay before and after each pulse, which is tcp. I.e., between two 180deg pulses there are two delays tcp.

Thank you for the additional details.  I have now updated the
lib/dispersion/ns_2site_star.py file with the details
(http://article.gmane.org/gmane.science.nmr.relax.scm/18046).


I saw that you keep the convention of "G" and "E" in the functions. Is this also used in the existing functions, or do you use "A" and "B" there? It might be cleaner to have a coherent naming convention.

I am in the process of converting.  I think I'll use the arbitrary mapping:

G -> A
E -> B

This is for consistency within relax between this code and the
analytic models, as well as to minimise user confusion.  It will also
allow relax to more logically support a C state for the 3-site models.


Here are some comments to the function:


70 # Initialise some structures.
71 Rr = -1.0 * matrix([[r20b, 0.0],[0.0, r20a]])
72 Rex = -1.0 * matrix([[kge, -keg],[-kge, keg]])
73 RCS = complex(0.0, -1.0) * matrix([[0.0, 0.0],[0.0, fg]])
74 R = Rr + Rex + RCS
75 cR2 = conj(R) * 2.0
76
77 kex = kge + keg # conversion of kinetic rates
78 IGeq = keg / kex # calculate relative populations - will be used for generating the equilibrium magnetizations
79 IEeq = kge / kex # like preceding line, but for excited state
80 M0 = matrix([[IGeq], [IEeq]]) # This is a vector that contains the initial magnetizations, corresponding to ground (G) and excited (E) state transverse magnetizations
81
82 # Replicated calculations for faster operation.
83 inv_tcpmg = 1.0 / tcpmg
84
85 # Loop over the time points, back calculating the R2eff values.
86 for i in range(num_points):
87 # Replicated calculations for faster operation.
88 tcp = 0.25 / cpmg_frqs[i]
89 expm_R_tcp = expm(R*tcp) # This matrix is a propagator that will evolve the magnetization with the matrix R for a delay tcp
90
91 prop_2 = dot(dot(expm_R_tcp, expm(cR2*tcp)), expm_R_tcp) # This is the propagator for an element of delay tcp - 180degpulse -2 times delay tcp - 180degpulse - delay tau, i.e. for 2 times tau-180-tau
92
93 prop_total = square_matrix_power(prop_2, cpmg_frqs[i]*tcpmg) # Now create the total propagator, that will evolve the magnetization under the CPMG train, ie. it applies the above tau-180-tau-tau-180-tau so many times as required for the CPMG frequency under consideration.
94
95 Moft = prop_total * M0 # now we apply the above propagator to the initial magnetization vector - resulting in the magnetization that remains after the full CPMG pulse train. It is called M of t (t is the time after the CPMG train
96
97 Mgx = Moft[0].real / M0[0] # This and the next line calculate the R2eff, using a two-point approximation, i.e. assuming that the decay is mono-exponential.
98 back_calc[i]= -inv_tcpmg * log(Mgx)

These comments are awesome!!!  It's far better than what I could have
added.  The code in lib/dispersion/ns_2site_star.py is now a perfect
teaching aid that can be shown directly to students.  It will also
allow NMR spectroscopists who are not experts to quickly learn from
the code.  I have added all of these comments and updated your
copyright notice appropriately.


Thanks for the implementation!

Thank you for contributing the code.  I hope this will up and running
soon and, together with the tools built into relax, that the
calculations will even be faster.

Cheers,

Edward



***************************

Paul Schanda





Related Messages


Powered by MHonArc, Updated Wed Jul 17 11:20:10 2013