Hi Edward. Well, I could not figure out what an interpolated time point would be for new points. So I collect all time points at the original offset and dispersion, and make a concatenated list. Then I make this unique.
The decision for the "correct" time point will be in the target function. For minimising a model, the corresponding time point can be extracted. But for interpolated points, I simply don't know. Best Troels 2014-08-22 14:57 GMT+02:00 Edward d'Auvergne <edw...@nmr-relax.com>: > Hi Troels, > > This all looks great. I was just wondering what the relaxation time > interpolation is for? Is this interpolation or simply zero-filling > the end of the relaxation time lists? I don't fully understand what > is happening with the concatenate() and unique() section and what is > created. > > Cheers, > > Edward > > > > On 22 August 2014 14:48, <tlin...@nmr-relax.com> wrote: >> Author: tlinnet >> Date: Fri Aug 22 14:48:11 2014 >> New Revision: 25214 >> >> URL: http://svn.gna.org/viewcvs/relax?rev=25214&view=rev >> Log: >> Modified following functions: >> >> Time points are now saved at the [ei][mi][oi][di] index level. >> At this index Ãlevelall time points are saved for the R2eff point. >> >> - interpolate_disp() >> To interpolate time points, all time points through the original dispersion >> points di, are >> collected and then made unique. This time list can potentially be the >> largest of all time lists. >> >> - interpolate_offset() >> To interpolate time points, all time points through the original offset >> points, and then dispersion points di, are >> collected and then made unique. This time list can potentially be the >> largest of all time lists. >> >> - plot_disp_curves_to_file() >> To acquire the original relax_times points. >> >> - return_r2eff_arrays() >> To save all time points on the level of [ei][mi][oi][di]. >> At this index level, it will be a numpy array list with all time values used >> for fitting. >> >> bug #22461(https://gna.org/bugs/?22461): NS R1rho 2-site_fit_r1 has >> extremely high chi2 value in systemtest >> Relax_disp.test_r1rho_kjaergaard_missing_r1. >> >> Modified: >> trunk/specific_analyses/relax_disp/data.py >> >> Modified: trunk/specific_analyses/relax_disp/data.py >> URL: >> http://svn.gna.org/viewcvs/relax/trunk/specific_analyses/relax_disp/data.py?rev=25214&r1=25213&r2=25214&view=diff >> ============================================================================== >> --- trunk/specific_analyses/relax_disp/data.py (original) >> +++ trunk/specific_analyses/relax_disp/data.py Fri Aug 22 14:48:11 2014 >> @@ -52,7 +52,7 @@ >> >> # Python module imports. >> from math import cos, pi, sin, sqrt >> -from numpy import array, float64, int32, ones, zeros >> +from numpy import array, concatenate, float64, int32, ones, unique, zeros >> from os import F_OK, access >> from os.path import expanduser >> from random import gauss >> @@ -775,7 +775,7 @@ >> desel_spin(spin_id) >> >> >> -def interpolate_disp(spin=None, spin_id=None, si=None, num_points=None, >> extend_hz=None): >> +def interpolate_disp(spin=None, spin_id=None, si=None, num_points=None, >> extend_hz=None, relax_times=None): >> """Interpolate function for 2D Grace plotting function for the >> dispersion curves. >> >> @keyword spin: The specific spin data container. >> @@ -788,6 +788,8 @@ >> @type num_points: int >> @keyword extend_hz: How far to extend the interpolated fitted >> curves to (in Hz). >> @type extend_hz: float >> + @keyword relax_times: The experiment specific fixed time period for >> relaxation (in seconds). The dimensions are {Ei, Mi, Oi, Di, Ti}. >> + @type relax_times: rank-4 list of floats >> @return: The interpolated_flag, list of back calculated >> R2eff/R1rho values in rad/s {Ei, Si, Mi, Oi, Di}, list of interpolated >> frequencies for cpmg_frqs in Hz {Ei, Si, Mi, Oi, Di}, interpolated spin-lock >> offsets in rad/s {Ei, Si, Mi, Oi}, list of interpolated spin-lock field >> strength frequencies for spin_lock_nu1_new in Hz {Ei, Si, Mi, Oi, Di}, >> chemical shifts in rad/s {Ei, Si, Mi}, interpolated rotating frame tilt >> angles theta {Ei, Si, Mi, Oi, Di}, interpolated average resonance offset in >> the rotating frame Omega in rad/s {Ei, Si, Mi, Oi, Di} and the interpolated >> effective field in rotating frame w_eff in rad/s {Ei, Si, Mi, Oi, Di}. >> @rtype: boolean, rank-4 list of numpy rank-1 float >> arrays, rank-4 list of numpy rank-1 float arrays, rank-3 list of numpy >> rank-1 float arrays, rank-4 list of numpy rank-1 float arrays, rank-2 list >> of numpy rank-1 float arrays, rank-4 list of numpy rank-1 float arrays, >> rank-4 list of numpy rank-1 float arrays, rank-4 list of numpy rank-1 float >> arrays >> """ >> @@ -797,6 +799,7 @@ >> # Initialise some structures. >> cpmg_frqs_new = None >> spin_lock_nu1_new = None >> + relax_times_new = None >> >> # Interpolate the CPMG frequencies (numeric models). >> if spin.model in MODEL_LIST_NUMERIC_CPMG or spin.model in [MODEL_B14, >> MODEL_B14_FULL]: >> @@ -871,28 +874,42 @@ >> >> if spin_lock_nu1 != None and len(spin_lock_nu1[0][0][0]): >> spin_lock_nu1_new = [] >> + relax_times_new = [] >> for ei in range(len(spin_lock_nu1)): >> # Add a new dimension. >> spin_lock_nu1_new.append([]) >> + relax_times_new.append([]) >> >> # Then loop over the spectrometer frequencies. >> for mi in range(len(spin_lock_nu1[ei])): >> # Add a new dimension. >> spin_lock_nu1_new[ei].append([]) >> + relax_times_new[ei].append([]) >> >> # Finally the offsets. >> for oi in range(len(spin_lock_nu1[ei][mi])): >> # Add a new dimension. >> spin_lock_nu1_new[ei][mi].append([]) >> + relax_times_new[ei][mi].append([]) >> >> # No data. >> if not len(spin_lock_nu1[ei][mi][oi]): >> continue >> + >> + # There is no way to interpolate the time points >> correct. >> + # The best suggestion is to concatenate all values at >> original offset, and then make a unique list. >> + relax_time_temp = array([]) >> + for di_o, times in enumerate(relax_times[ei][mi][oi]): >> + relax_time_temp = concatenate( (relax_time_temp, >> times) ) >> + >> + # Make a unique list. >> + relax_time_temp = unique(relax_time_temp) >> >> # Interpolate (adding the extended amount to the end). >> for di in range(num_points): >> point = (di + 1) * >> (max(spin_lock_nu1[ei][mi][oi])+extend_hz) / num_points >> spin_lock_nu1_new[ei][mi][oi].append(point) >> + relax_times_new[ei][mi][oi].append(relax_time_temp) >> >> # Convert to a numpy array. >> spin_lock_nu1_new[ei][mi][oi] = >> array(spin_lock_nu1_new[ei][mi][oi], float64) >> @@ -914,12 +931,12 @@ >> back_calc = None >> else: >> # Back calculate R2eff data for the second sets of plots. >> - back_calc = >> specific_analyses.relax_disp.optimisation.back_calc_r2eff(spin=spin, >> spin_id=spin_id, cpmg_frqs=cpmg_frqs_new, spin_lock_nu1=spin_lock_nu1_new) >> + back_calc = >> specific_analyses.relax_disp.optimisation.back_calc_r2eff(spin=spin, >> spin_id=spin_id, cpmg_frqs=cpmg_frqs_new, spin_lock_nu1=spin_lock_nu1_new, >> relax_times_new=relax_times_new) >> >> return interpolated_flag, back_calc, cpmg_frqs_new, offsets, >> spin_lock_fields_inter, chemical_shifts, tilt_angles, Delta_omega, w_eff >> >> >> -def interpolate_offset(spin=None, spin_id=None, si=None, num_points=None, >> extend_ppm=None): >> +def interpolate_offset(spin=None, spin_id=None, si=None, num_points=None, >> extend_ppm=None, relax_times=None): >> """Interpolate function for 2D Grace plotting function for the >> dispersion curves, interpolating through spin-lock offset in rad/s. >> >> @keyword spin: The specific spin data container. >> @@ -932,6 +949,8 @@ >> @type num_points: int >> @keyword extend_ppm: How far to extend the interpolated fitted >> curves to in offset ppm. >> @type extend_ppm: float >> + @keyword relax_times: The experiment specific fixed time period for >> relaxation (in seconds). The dimensions are {Ei, Mi, Oi, Di, Ti}. >> + @type relax_times: rank-4 list of floats >> @return: The interpolated_flag, list of back calculated >> R2eff/R1rho values in rad/s {Ei, Si, Mi, Oi, Di}, list of interpolated >> frequencies for cpmg_frqs in Hz {Ei, Si, Mi, Oi, Di}, interpolated spin-lock >> offsets in rad/s {Ei, Si, Mi, Oi}, list of interpolated spin-lock field >> strength frequencies for spin_lock_nu1_new in Hz {Ei, Si, Mi, Oi, Di}, >> chemical shifts in rad/s {Ei, Si, Mi}, interpolated rotating frame tilt >> angles theta {Ei, Si, Mi, Oi, Di}, interpolated average resonance offset in >> the rotating frame Omega in rad/s {Ei, Si, Mi, Oi, Di} and the interpolated >> effective field in rotating frame w_eff in rad/s {Ei, Si, Mi, Oi, Di}. >> @rtype: boolean, rank-4 list of numpy rank-1 float >> arrays, rank-4 list of numpy rank-1 float arrays, rank-3 list of numpy >> rank-1 float arrays, rank-4 list of numpy rank-1 float arrays, rank-2 list >> of numpy rank-1 float arrays, rank-4 list of numpy rank-1 float arrays, >> rank-4 list of numpy rank-1 float arrays, rank-4 list of numpy rank-1 float >> arrays >> """ >> @@ -941,6 +960,7 @@ >> >> # Initialise some structures. >> spin_lock_offset_new = [] >> + relax_times_new = None >> >> # Get the spin-lock field strengths. >> spin_lock_nu1 = return_spin_lock_nu1(ref_flag=False) >> @@ -985,11 +1005,42 @@ >> # The offset data. >> offsets, spin_lock_fields_inter, chemical_shifts, tilt_angles, >> Delta_omega, w_eff = return_offset_data(spins=[spin], spin_ids=[spin_id], >> field_count=field_count, spin_lock_offset=spin_lock_offset_new, >> fields=spin_lock_nu1) >> >> + # Interpolated relaxation time. >> + if tilt_angles != None and len(tilt_angles[0][0][0]): >> + relax_times_new = [] >> + for ei in range(len(tilt_angles)): >> + # Add a new dimension. >> + relax_times_new.append([]) >> + >> + # Then loop over the spectrometer frequencies. >> + for mi in range(len(tilt_angles[ei])): >> + # Add a new dimension. >> + relax_times_new[ei].append([]) >> + >> + # There is no way to interpolate the time points correct. >> + # The best suggestion is to concatenate all values at >> original offset and dispersion point, and then make a unique list. >> + relax_time_temp = array([]) >> + for oi_o, relax_times_oi in enumerate(relax_times[ei][mi]): >> + for di_o, times in enumerate(relax_times_oi): >> + relax_time_temp = concatenate( (relax_time_temp, >> times) ) >> + >> + # Make a unique list. >> + relax_time_temp = unique(relax_time_temp) >> + >> + # Finally the offsets. >> + for oi in range(len(tilt_angles[ei][0][mi])): >> + # Add a new dimension. >> + relax_times_new[ei][mi].append([]) >> + >> + # Interpolate (adding the extended amount to the end). >> + for di in range(len(tilt_angles[ei][0][mi][oi])): >> + relax_times_new[ei][mi][oi].append(relax_time_temp) >> + >> if spin.model == MODEL_R2EFF: >> back_calc = None >> else: >> # Back calculate R2eff data for the second sets of plots. >> - back_calc = >> specific_analyses.relax_disp.optimisation.back_calc_r2eff(spin=spin, >> spin_id=spin_id, spin_lock_offset=spin_lock_offset_new, >> spin_lock_nu1=spin_lock_fields_inter) >> + back_calc = >> specific_analyses.relax_disp.optimisation.back_calc_r2eff(spin=spin, >> spin_id=spin_id, spin_lock_offset=spin_lock_offset_new, >> spin_lock_nu1=spin_lock_fields_inter, relax_times_new=relax_times_new) >> >> # cpmg_frqs are not interpolated. >> cpmg_frqs_new = None >> @@ -1912,6 +1963,16 @@ >> linetype = [] >> linestyle = [] >> >> + # Number of spectrometer fields. >> + fields = [None] >> + field_count = 1 >> + if hasattr(cdp, 'spectrometer_frq_count'): >> + fields = cdp.spectrometer_frq_list >> + field_count = cdp.spectrometer_frq_count >> + >> + # Get the relax_times. >> + values, errors, missing, frqs, frqs_H, exp_types, relax_times = >> return_r2eff_arrays(spins=[spin], spin_ids=[spin_id], fields=fields, >> field_count=field_count) >> + >> # Set up the interpolated curve data structures. >> interpolated_flag = False >> >> @@ -1920,11 +1981,11 @@ >> >> if interpolate == INTERPOLATE_DISP: >> # Interpolate through disp points. >> - interpolated_flag, back_calc, cpmg_frqs_new, offsets_inter, >> spin_lock_nu1_new, chemical_shifts, tilt_angles_inter, Delta_omega_inter, >> w_eff_inter = interpolate_disp(spin=spin, spin_id=spin_id, si=si, >> num_points=num_points, extend_hz=extend_hz) >> + interpolated_flag, back_calc, cpmg_frqs_new, offsets_inter, >> spin_lock_nu1_new, chemical_shifts, tilt_angles_inter, Delta_omega_inter, >> w_eff_inter = interpolate_disp(spin=spin, spin_id=spin_id, si=si, >> num_points=num_points, extend_hz=extend_hz, relax_times=relax_times) >> >> elif interpolate == INTERPOLATE_OFFSET: >> # Interpolate through disp points. >> - interpolated_flag, back_calc, cpmg_frqs_new, offsets_inter, >> spin_lock_nu1_new, chemical_shifts, tilt_angles_inter, Delta_omega_inter, >> w_eff_inter = interpolate_offset(spin=spin, spin_id=spin_id, si=si, >> num_points=num_points, extend_ppm=extend_ppm) >> + interpolated_flag, back_calc, cpmg_frqs_new, offsets_inter, >> spin_lock_nu1_new, chemical_shifts, tilt_angles_inter, Delta_omega_inter, >> w_eff_inter = interpolate_offset(spin=spin, spin_id=spin_id, si=si, >> num_points=num_points, extend_ppm=extend_ppm, relax_times=relax_times) >> >> # Do not interpolate, if model is R2eff. >> if spin.model == MODEL_R2EFF: >> @@ -4318,12 +4379,14 @@ >> missing[ei][si].append([]) >> frqs[ei][si].append(0.0) >> frqs_H[ei][si].append(0.0) >> + relax_times[ei].append([]) >> for offset, oi in loop_offset(exp_type=exp_type, frq=frq, >> return_indices=True): >> values[ei][si][mi].append([]) >> errors[ei][si][mi].append([]) >> missing[ei][si][mi].append([]) >> - for mi in range(field_count): >> - relax_times[ei].append(None) >> + relax_times[ei][mi].append([]) >> + for point, di in loop_point(exp_type=exp_type, frq=frq, >> offset=offset, return_indices=True): >> + relax_times[ei][mi][oi].append([]) >> >> # Pack the R2eff/R1rho data. >> data_flag = False >> @@ -4405,34 +4468,8 @@ >> errors[ei][si][mi][oi].append(current_spin.r2eff_err[key]) >> >> # The relaxation times. >> - relax_time = [] >> - for id in cdp.spectrum_ids: >> - # Non-matching data. >> - if cdp.spectrometer_frq[id] != frq: >> - continue >> - if cdp.exp_type[id] != exp_type: >> - continue >> - if exp_type in EXP_TYPE_LIST_CPMG: >> - if id not in cdp.cpmg_frqs.keys() or cdp.cpmg_frqs[id] >> != point: >> - continue >> - else: >> - if id not in cdp.spin_lock_nu1.keys() or >> cdp.spin_lock_nu1[id] != point: >> - continue >> - >> - # Found. >> - relax_time.append(cdp.relax_times[id]) >> - >> - # Use the maximum time value found. >> - relax_time = max(relax_time) >> - >> - # Check the value if already set. >> - if relax_times[ei][mi] != None: >> - if relax_times[ei][mi] != relax_time: >> - raise RelaxError("The relaxation times do not match for >> all experiments.") >> - continue >> - >> - # Store the time. >> - relax_times[ei][mi] = relax_time >> + for time, ti in loop_time(exp_type=exp_type, frq=frq, >> offset=offset, point=point, return_indices=True): >> + relax_times[ei][mi][oi][di].append(time) >> >> # Increment the spin index. >> si += 1 >> @@ -4442,7 +4479,6 @@ >> raise RelaxError("No R2eff/R1rho data could be found for the spin >> cluster %s." % spin_ids) >> >> # Convert to numpy arrays. >> - relax_times = array(relax_times, float64) >> for exp_type, ei in loop_exp(return_indices=True): >> for si in range(spin_num): >> for frq, mi in loop_frq(return_indices=True): >> @@ -4450,6 +4486,8 @@ >> values[ei][si][mi][oi] = array(values[ei][si][mi][oi], >> float64) >> errors[ei][si][mi][oi] = array(errors[ei][si][mi][oi], >> float64) >> missing[ei][si][mi][oi] = >> array(missing[ei][si][mi][oi], int32) >> + for point, di in loop_point(exp_type=exp_type, frq=frq, >> offset=offset, return_indices=True): >> + relax_times[ei][mi][oi][di] = >> array(relax_times[ei][mi][oi][di], float64) >> >> # Return the structures. >> return values, errors, missing, frqs, frqs_H, exp_types, relax_times >> >> >> _______________________________________________ >> relax (http://www.nmr-relax.com) >> >> This is the relax-commits mailing list >> relax-comm...@gna.org >> >> To unsubscribe from this list, get a password >> reminder, or change your subscription options, >> visit the list information page at >> https://mail.gna.org/listinfo/relax-commits > > _______________________________________________ > relax (http://www.nmr-relax.com) > > This is the relax-devel mailing list > relax-devel@gna.org > > To unsubscribe from this list, get a password > reminder, or change your subscription options, > visit the list information page at > https://mail.gna.org/listinfo/relax-devel _______________________________________________ relax (http://www.nmr-relax.com) This is the relax-devel mailing list relax-devel@gna.org To unsubscribe from this list, get a password reminder, or change your subscription options, visit the list information page at https://mail.gna.org/listinfo/relax-devel
