Hi, If that works for the back-calculation, then that's ok. I see now that the unique() function call is required as the ID string matching has been removed, so you will end up with 2 times for duplicated measurements.
Cheers, Edward On 22 August 2014 15:13, Troels Emtekær Linnet <tlin...@nmr-relax.com> wrote: > 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