Plot peaks for each column in a csv then fit a curve.

I am trying to process some data obtained from an oscilloscope. The columns in the csv represent the amplitude of the signal. I expected to see what is in the picture. I want to find the peaks for each signal and then plot the peaks as a function of tau (ms) ranging from 0.1 and up. I'd also like to fit a curve to each line following the exponential function y = a*e(-x/T1) +c, where a is a constant (-2) and T1 is found from the curve. Finding T1 is the ultimate goal of plotting the signals. I have tried a few different methods: here is my most recent attempt and it did not give me anywhere near what I was hoping for:

filename = 'D:\Grad Lab\NMR\Data\T1 Data\compiledT1nolabel.csv';
data = readtable(filename);
y1= table2array(data(:,1));
y2= table2array(data(:,2));
y3= table2array(data(:,3));
y4= table2array(data(:,4));
y5= table2array(data(:,5));
y6= table2array(data(:,6));
y7= table2array(data(:,7));
Fs = 2e-10;
T = 1/Fs;
L = 2500;
t = (0:L-1)*T;
Fn = Fs/2;
Fy1 = fft(y1)/L;
Fv = linspace(0, 1, fix(L/2)+1)*Fn;
Iv = 1:length(Fv);
figure(1)
plot(Fv, abs(Fy1(Iv))*2)
grid
title("Fourier Transfrom of Mineral Oil Original Signal")
xlabel('Frequency (Hz)')
ylabel('Amplitude')
Fy1dcoc = fft(y1-mean(y1))/L;
figure(2)
plot(Fv, abs(Fy1dcoc(Iv))*2)
grid
title('Fourier Transform of D-C Offset Corrected Signal')
xlabel('Frequency (Hz)')
ylabel('Amplitude')
[Fy1n_max,Iv_max] = max(abs(Fy1dcoc(Iv))*2);
Wp = 2*Fv(Iv_max)/Fn;
Ws = Wp*2;
Rp = 10
Rs = 30;
[n,Wn] = buttord(Wp,Ws,Rp,Rs);
[b,a] = butter(n,Wn);
[SOS,G] = tf2sos(b,a);
S = filtfilt(SOS,G,y1);
figure(4)
plot(t, y1)
hold on
plot(t, S, '-r', 'LineWidth',1.5)
hold off
grid
legend('Mineral oil', "'S'", 'Location', 'N')
title('Original Signal of Mineral Oil and Uncorrupted Signal (S)')
xlabel('Tau (ms)')
ylabel('Amplitude')

I appreciate any help with this!

Answer 1

[Star Strider]

There are no visible peaks in the time-domain signal.

Doing the Foureir transform and finding the peaks is easy enough, however I have no idea how to do the plot in the image because I have no idea how to calculate whatever the variables are that it depicts. The independent variable appears to be time.

If your intent is to do a time-frequency analysis, it would be best to use the pspectrum function with the 'spectrogram' option, rather than fft.

data = readtable('https://www.mathworks.com/matlabcentral/answers/uploaded_files/1178628/compiledT1.csv', 'VariableNamingRule','preserve')

data = 2500×7 table
    Mineral Oil    Glycerol    CuSO4 0.625%    CuSO4 1.25%    CuSO4 2.5%    CuSO4 5%    CuSO4 10%
    ___________    ________    ____________    ___________    __________    ________    _________

       -0.002       -0.089        0.0032         -0.0136        -0.0036      -0.0036     -0.0118 
     -0.00196      -0.0889       0.00324        -0.01356       -0.00358     -0.00358    -0.01178 
     -0.00192      -0.0888       0.00328        -0.01352       -0.00356     -0.00356    -0.01176 
     -0.00188      -0.0887       0.00332        -0.01348       -0.00354     -0.00354    -0.01174 
     -0.00184      -0.0886       0.00336        -0.01344       -0.00352     -0.00352    -0.01172 
      -0.0018      -0.0885        0.0034         -0.0134        -0.0035      -0.0035     -0.0117 
     -0.00176      -0.0884       0.00344        -0.01336       -0.00348     -0.00348    -0.01168 
     -0.00172      -0.0883       0.00348        -0.01332       -0.00346     -0.00346    -0.01166 
     -0.00168      -0.0882       0.00352        -0.01328       -0.00344     -0.00344    -0.01164 
     -0.00164      -0.0881       0.00356        -0.01324       -0.00342     -0.00342    -0.01162 
      -0.0016       -0.088        0.0036         -0.0132        -0.0034      -0.0034     -0.0116 
     -0.00156      -0.0879       0.00364        -0.01316       -0.00338     -0.00338    -0.01158 
     -0.00152      -0.0878       0.00368        -0.01312       -0.00336     -0.00336    -0.01156 
     -0.00148      -0.0877       0.00372        -0.01308       -0.00334     -0.00334    -0.01154 
     -0.00144      -0.0876       0.00376        -0.01304       -0.00332     -0.00332    -0.01152 
      -0.0014      -0.0875        0.0038          -0.013        -0.0033      -0.0033     -0.0115 

VN = data.Properties.VariableNames;

Fs = 2e-10;
T = 1/Fs;
L = 2500;
t = (0:L-1)*T;
Fn = Fs/2;

y = table2array(data)

y = 2500×7

-0.0020 -0.0890 0.0032 -0.0136 -0.0036 -0.0036 -0.0118 -0.0020 -0.0889 0.0032 -0.0136 -0.0036 -0.0036 -0.0118 -0.0019 -0.0888 0.0033 -0.0135 -0.0036 -0.0036 -0.0118 -0.0019 -0.0887 0.0033 -0.0135 -0.0035 -0.0035 -0.0117 -0.0018 -0.0886 0.0034 -0.0134 -0.0035 -0.0035 -0.0117 -0.0018 -0.0885 0.0034 -0.0134 -0.0035 -0.0035 -0.0117 -0.0018 -0.0884 0.0034 -0.0134 -0.0035 -0.0035 -0.0117 -0.0017 -0.0883 0.0035 -0.0133 -0.0035 -0.0035 -0.0117 -0.0017 -0.0882 0.0035 -0.0133 -0.0034 -0.0034 -0.0116 -0.0016 -0.0881 0.0036 -0.0132 -0.0034 -0.0034 -0.0116

figure
plot(t, y)
grid
xlabel('Time')
ylabel('Amplitude')
legend(VN, 'Location','best')

ydt1 = detrend(y(:,1),1);

figure
plot(t, ydt1)
grid
xlabel('Time')
ylabel('Amplitude')
title(VN{1})

[p,f,t] = pspectrum(ydt1,Fs,'spectrogram');

figure
waterfall(f,t,p')
xlabel('Frequency (Hz)')
ylabel('Time (seconds)')
wtf = gca;
wtf.XDir = 'reverse';
colormap(turbo)
view([30 45])

NFFT = 2^nextpow2(L);
ym = y - mean(y);
FTy = fft(ym, NFFT)/L;
Fv = linspace(0, 1, NFFT/2+1)*Fn;
Iv = 1:numel(Fv);

figure
plot(Fv, abs(FTy(Iv,:))*2)
grid
xlabel('Frequency')
ylabel('Magnitude')
legend(VN, 'Location','best')
xlim([0 5E-12])

[pks1,locs1] = findpeaks(abs(FTy(Iv,1))*2);

figure
plot(Fv, abs(FTy(Iv,1))*2, 'DisplayName','Data')
hold on
plot(Fv(locs1), pks1, '^r', 'DisplayName','Peaks')
hold off
grid
xlabel('Frequency')
ylabel('Magnitude')
title(VN{1})
legend()
% legend(VN{1}, 'Location','best')
xlim([0 5E-12])

.

[Alex]

data = readtable(amplitudeT1.xlsx)
VN = data.Properties.VariableNames;
Fs = 1e-10;
T = 1/Fs;
L = 2500;
t = (0:L-1)*T;
Fn = Fs/2;
y = table2array(data);
figure(1)
plot(t, y)
grid
xlabel('Time (ms)')
ylabel('Amplitude')
legend(VN, 'Location', 'best')
ydt1 = detrend(y(:,1),1);
figure(2)
plot(t, ydt1)
grid
xlabel('Time')
ylabel('Amplitude')
title(VN{1})
[p,f,t] = pspectrum(ydt1,Fs,'spectrogram');
figure(3)
waterfall(f,t,p')
xlabel('Frequency (Hz)')
ylabel('Time')
wtf = gca;
wtf.xDir = 'reverse';
colormap(turbo)
view([30 45])
NFFT = 2^nextpow2(L);
ym = y - mean(y);
FTy = fft(ym, NFFT)/L;
Fv = linspace(0, 1, NFFT/2+1)*Fn;
Iv = 1:numel(Fv);
figure(4)
plot(Fv, abs(FTy(Iv,:))*2)
grid
xlabel('Frequency')
ylabel('Magnitude')
legend(VN, 'Location', 'best')
xlim([0 5e-12])
[pks1,locs1] = findpeaks(abs(FTy(Iv,1))*2);
figure(5)
plot(Fv, abs(FTy(Iv,1))*2), 'DisplayName','Data')
hold on
plot(Fv(locs1), pks1, '^r', 'DisplayName', 'Peaks')
hold off
grid
xlabel('Frequency')
ylabel('Magnitude')
title(VN{1})
legend(VN{1}, 'Location', 'best')
xlim([0 5E-12])

Running this it looks great (all the plots), I get the same plots you saw. I attached the incorrect file so I attached the correct one. I think the best thing for me to do would be to just plot each one separately though seeing how busy these plots. For the plot attached I want to grab that highest peak and then every other peak after that. Then if I can fit a curve to those peaks I can find T1 (the relaxation time of the atom.) Thank you!

[Star Strider]

It took a while to get this to work. This is the best I can do with these data.

Part of the problem is that the magnitude of the ‘t’ vector puts many of the calculations near the limit of floating-point precision for addition and subtraction (the eps value). That is also the reason that it’s not possible to design a functioning digital filter for these data. (I gave that a shot as well.)

Also, the ‘a’ value of -2 will absolutely not work unless the data themsellves are also negative. I let that float here.

Editing these data further, for example starting at the maximum peak, (the fourth peak for the Mineral Oil data) may produce a better fit. I don’t know what these data are or what they represent, so I didn’t edit them. I also didn’t detrend them, although they have a wandering baseline.

The Fourier transform data are interesting. All the data have peaks at either or , and some have both.

Try this —

data = readtable('https://www.mathworks.com/matlabcentral/answers/uploaded_files/1178803/amplitudeT1.xlsx', 'VariableNamingRule','preserve')

data = 2500×7 table
    Mineral Oil    Glycerol    CuSO4 0.625%    CuSO4 1.25%    CuSO4 2.5%    CuSO4 5%    CuSO4 10%
    ___________    ________    ____________    ___________    __________    ________    _________

           0        -0.08         -0.08               0             0            0           0   
           0        -0.08         -0.16               0             0            0           0   
           0        -0.08         -0.16               0             0            0           0   
           0        -0.08         -0.16               0             0            0           0   
           0        -0.08         -0.16               0             0            0           0   
           0        -0.08         -0.16               0         -0.08            0        0.08   
           0        -0.08         -0.08               0             0            0           0   
           0            0             0               0             0        -0.08           0   
           0        -0.08             0               0         -0.08            0           0   
       -0.08        -0.08             0               0         -0.08            0           0   
           0        -0.08             0           -0.08         -0.08        -0.08           0   
           0            0          0.08               0             0            0           0   
           0        -0.08          0.16               0         -0.08        -0.08           0   
           0            0          0.24               0             0        -0.08           0   
           0        -0.08          0.32           -0.08             0        -0.08           0   
        0.08            0          0.56               0             0            0           0   

VN = data.Properties.VariableNames;

Fs = 2e-10;
T = 1/Fs;
L = 2500;
t = (0:L-1)*T;
Fn = Fs/2;

y = table2array(data);

figure
plot(t, y)
grid
xlabel('Time')
ylabel('Amplitude')
legend(VN, 'Location','best')

% Lminv = islocalmin(y(:,1), 'MinProminence',1);
% p = polyfit(t(Lminv).', y(Lminv,1), 9);
% ydt1 = y(:,1) - polyval(p,t(:));
% ydt1 = detrend(y(:,1),15);

[pkst1, locst1] = findpeaks(y(:,1), 'MinPeakProminence',4.5);

figure
% plot(t, ydt1)
plot(t, y(:,1))
hold on
plot(t(locst1), pkst1, '^r')
% plot(t(Lminv), y(Lminv,1), '+r')
hold off
grid
xlabel('Time')
ylabel('Amplitude')
title(VN{1})

format longE
p = polyfit(t(locst1), log(pkst1), 1)                   % Estimate Initial Parameters For Exponential Fit

p = 1×2

-2.381757003380055e-13 2.247436395545789e+00

ypv = exp(polyval(p, t(locst1)));

figure
plot(t(locst1), pkst1,'.b',   t(locst1),ypv, '-r')
grid
xlabel('Time')
ylabel('Amplitude')
title('Initial Parameter Estimate Result')

% y = a*e(-x/T1) +c
objfcn = @(b,t) b(1).*exp(-t./b(2)) + b(3);
mdl = fitnlm(t(locst1), pkst1, objfcn, [p(2); -1/(p(1)); 1])

Warning: Rank deficient, rank = 2, tol = 1.547956e-14.

Warning: Some columns of the Jacobian are effectively zero at the solution, indicating that the model is insensitive to some of its parameters. That may be because those parameters are not present in the model, or otherwise do not affect the predicted values. It may also be due to numerical underflow in the model function, which can sometimes be avoided by choosing better initial parameter values, or by rescaling or recentering. Parameter estimates may be unreliable.

mdl = 
Nonlinear regression model:
    y ~ b1*exp( - t/b2) + b3

Estimated Coefficients:
                Estimate                    SE                     tStat                   pValue        
          _____________________    ____________________    _____________________    _____________________

    b1     6.42850666045294e+01    1.20005318013054e+01     5.35685148532637e+00     1.30448378489792e-04
    b2     4.39486428989028e+13    3.18060840239682e-11     1.38176843354197e+24    5.64414374960479e-308
    b3    -5.53061960933816e+01    1.15303083591107e+01    -4.79659297660338e+00     3.48966691203665e-04


Number of observations: 15, Error degrees of freedom: 13
Root Mean Squared Error: 1.06
R-Squared: 0.688,  Adjusted R-Squared 0.664
F-statistic vs. constant model: 28.7, p-value = 0.00013

tv = linspace(min(t(locst1)), max(t(locst1)));
[yv,yci] = predict(mdl, tv(:));

figure
hp{1} = plot(t(locst1), pkst1,'.b', 'DisplayName','Data');
hold on
hp{2}  = plot(tv,yv,'-r', 'DisplayName','Exponential Fit');
hp{3} = plot(tv,yci,'--r', 'DisplayName','±95% CI');
hold off
grid
xlabel('Time')
ylabel('Amplitude')
title(VN{1})
legend([hp{1};hp{2};hp{3}(1)],'Location','best')

[p,f,t] = pspectrum(y(:,1),Fs,'spectrogram');

figure
waterfall(f,t,p')
xlabel('Frequency (Hz)')
ylabel('Time (seconds)')
wtf = gca;
wtf.XDir = 'reverse';
colormap(turbo)
xlim([0 3E-11])
view([135 45])

NFFT = 2^nextpow2(L);
ym = y - mean(y);
FTy = fft(ym, NFFT)/L;
Fv = linspace(0, 1, NFFT/2+1)*Fn;
Iv = 1:numel(Fv);

figure
plot(Fv, abs(FTy(Iv,:))*2)
grid
xlabel('Frequency')
ylabel('Magnitude')
legend(VN, 'Location','best')
xlim([0 5E-12])

[pksf1,locsf1] = findpeaks(abs(FTy(Iv,1))*2);

figure
plot(Fv, abs(FTy(Iv,1))*2, 'DisplayName','Data')
hold on
plot(Fv(locsf1), pksf1, '^r', 'DisplayName','Peaks')
hold off
grid
xlabel('Frequency')
ylabel('Magnitude')
title(VN{1})
legend()
% legend(VN{1}, 'Location','best')
xlim([0 5E-12])

.

[Alex]

The time step is quite large, I will change a few parameters to match the time step. It is increasing by a step size of 0.1 (ms) per peak. So it would be possible I think to cut down the size of the data and only go out to 1.0 ms.

I will take a look at this further and play around with these changes this evening! Will update on how it goes, thank you!

[Star Strider]

My pleasure!

[Alex]

I figured out what is up with the time step. I was thinking about the data incorrectly. So for each peak that is a one tau (0.1ms.) The reason being, this is an NMR experiment so there is a pulse A and a pulse B at tau after pulse A, then there is an echo peak at 2 tau, then pulse B at 3 tau, echo at 4 tau and so on.

I'm not sure if that will help with any of the coding. I'm not sure how to break that down!

[Alex]

Sorry, see my comment just added below!

[Star Strider]

The figure in your last Comment makes this a bit clearer.

If the peaks are supposed to be at 1 ms intervals, and you want to fit from the highest peak to the peaks greater than 2, this works —

data = readtable('https://www.mathworks.com/matlabcentral/answers/uploaded_files/1178803/amplitudeT1.xlsx', 'VariableNamingRule','preserve')

data = 2500×7 table
    Mineral Oil    Glycerol    CuSO4 0.625%    CuSO4 1.25%    CuSO4 2.5%    CuSO4 5%    CuSO4 10%
    ___________    ________    ____________    ___________    __________    ________    _________

           0        -0.08         -0.08               0             0            0           0   
           0        -0.08         -0.16               0             0            0           0   
           0        -0.08         -0.16               0             0            0           0   
           0        -0.08         -0.16               0             0            0           0   
           0        -0.08         -0.16               0             0            0           0   
           0        -0.08         -0.16               0         -0.08            0        0.08   
           0        -0.08         -0.08               0             0            0           0   
           0            0             0               0             0        -0.08           0   
           0        -0.08             0               0         -0.08            0           0   
       -0.08        -0.08             0               0         -0.08            0           0   
           0        -0.08             0           -0.08         -0.08        -0.08           0   
           0            0          0.08               0             0            0           0   
           0        -0.08          0.16               0         -0.08        -0.08           0   
           0            0          0.24               0             0        -0.08           0   
           0        -0.08          0.32           -0.08             0        -0.08           0   
        0.08            0          0.56               0             0            0           0   

VN = data.Properties.VariableNames;

Fs = 2e-10;
T = 1/Fs;
L = 2500;
t = (0:L-1)*T;
Fn = Fs/2;

y = table2array(data);

figure
plot(t, y)
grid
xlabel('Time')
ylabel('Amplitude')
legend(VN, 'Location','best')

[pkst1, locst1] = findpeaks(y(:,1));                    % Detect All Peaks

figure
plot(t, y(:,1))
hold on
plot(t(locst1), pkst1, '^r')
hold off
grid
xlabel('Time')
ylabel('Amplitude')
title(VN{1})

[mxpk,mxpkidx] = max(pkst1)

mxpk = 9.6000

mxpkidx = 10

Lvpeaks = (pkst1 <= max(pkst1)) & (pkst1 >= 2.45) & (locst1 >= locst1(mxpkidx));     % Select Peaks
locst1e = locst1(Lvpeaks)                          % Select ‘locs'

locst1e = 28×1

152 176 202 226 228 252 279 302 352 376

pkst1e = pkst1(Lvpeaks);                            % Select Corresponding 'pks'
tauv = (1:numel(pkst1e))*1E-3;                      % Create 'tau' Time Vector

format longE
% p = polyfit(t(locst1), log(pkst1), 1)                   % Estimate Initial Parameters For Exponential Fit
p = polyfit(tauv, log(pkst1e), 1)                   % Estimate Initial Parameters For Exponential Fit

p = 1×2

-3.600759660774221e+01 1.928444839224613e+00

ypv = exp(polyval(p, tauv));

figure
plot(tauv, pkst1e,'.b',   tauv,ypv, '-r')
grid
xlabel('\tau')
ylabel('Amplitude')
title('Initial Parameter Estimate Result')

% y = a*e(-x/T1) +c
objfcn = @(b,t) b(1).*exp(-t./b(2)) + b(3);
mdl = fitnlm(tauv, pkst1e, objfcn, [p(2); -1/(p(1)); 1])

mdl = 
Nonlinear regression model:
    y ~ b1*exp( - t/b2) + b3

Estimated Coefficients:
                Estimate                   SE                    tStat                   pValue       
          ____________________    ____________________    ____________________    ____________________

    b1    6.44645526000870e+00    2.75304770820915e+00    2.34157048596957e+00    2.74736774017030e-02
    b2    1.76121018054069e-02    1.92648742030917e-02    9.14207983905779e-01    3.69341137202206e-01
    b3    1.33724950949722e+00    3.37545923932747e+00    3.96168169923941e-01    6.95342079540181e-01


Number of observations: 28, Error degrees of freedom: 25
Root Mean Squared Error: 1.58
R-Squared: 0.473,  Adjusted R-Squared 0.43
F-statistic vs. constant model: 11.2, p-value = 0.000337

tv = linspace(min(tauv), max(tauv));
[yv,yci] = predict(mdl, tv(:));

figure
hp{1} = plot(tauv, pkst1e,'.b', 'DisplayName','Data');
hold on
hp{2} = plot(tv,yv,'-r', 'DisplayName','Exponential Fit');
hp{3} = plot(tv,yci,'--r', 'DisplayName','±95% CI');
hold off
grid
xlabel('\tau')
ylabel('Amplitude')
title(VN{1})
legend([hp{1};hp{2};hp{3}(1)],'Location','best')

This is the best I can do. I am not certain how to put ‘Peak 2’ specifically into the set of peaks to be regressed, or eliminate the others.

.

[Alex]

Thank you, I will continue to play around with it some more. This should work fine though for my data, i think it will really work well when I look at T2 and the free induction decay signals. I appreciate all the help truly!

[Star Strider]

My pleasure!

I also don’t understand the time values on the x-axis in the figure in your last Comment. In my code, I assigned each peak x-value a serial 0.1 millisecond value. That may not be correct, however that’s how I interpreted ‘So for each peak that is a one tau (0.1ms.)’ although it doesn’t make sense to me.

Answer 2

[Image Analyst]

See attached demo. It will be easy for you to adapt it to your data. Just replace demo data with your actual data.

[Alex]

Thank you, I am still getting just a straight line for my own data. How can I extract the peaks from my data and then plot only the peaks?

[Image Analyst]

The curve does not look like an inverted exponential:

Why do you want to fit the red line to an inverted exponential instead of some sigmoid function, like the rate equation (demo atached)?

[Alex]

I’m following some past work so I may need to tweak my parameters or I may have the equation values incorrect. I will triple check the equations I’m given and see what I’m missing. I believe exponential growth is what I’m after not the decay, decay would be for a different variable T2 (which changes the experiment setup so the equation would be different). I will try this out again this evening when I am home

[Image Analyst]

I just plotted s vs t, and y1 vs t like you did. I'm not sure why the plot I posted above looks so drastically different than @Star Strider's.

Well if s and y1 are your measured values (I doubt it since they're so smooth) and if you say that the theory behind the situation says that it's supposed to be an inverted exponential decay, then use an exponential decay.

I thought that my demo was pretty well commented. Did you have trouble replacing y with s, and x with your t? If you can't do it, I can do it for you. But essentially it's just passing x and y and the theoretical model into fitnlm() - just a few lines of code if you ignore the comments and fancy plotting stuff.

[Alex]

I am running into an error when I'm adding in my data, it plots it but there is so much noise I can't read the chart. Attached is the correct signal as well (I attached the wrong spreadsheet initially!) Also something to note, the time step is tau, and tau is sampled at 0.1 ms per peak after the first peak. So peak 1 is 0 peak 2 is 0.1, peak 3 0.2, and so on...

Thank you for your help!

[Image Analyst]

There are 7 columns in that workbook. What do the different columns represent? Which column is the x value and which are/is the y value(s)? What column are you trying to fit with an exponential?

You've accepted @Star Strider's answer already so are you still having trouble doing the fit?

Answer 3

[Alex]

Sorry for the confusion (if it's any consolation by going through this I am understanding the experiment more.) The signal in the 'attached figure' (not the picture) that is a signal from an oscilloscope. I drew arrows indicating 'peak 1', 'peak 2' and 'peak 3' (the first couple peaks before the labeled peak 1 are noise and can be deleted.) What is done we subject a magnetic field by emitting a pulse (peak 1) then we flip that field with a second pulse at a time (tau = 0.2ms) (peak 2) then there is an echo (peak 3) at 2 tau = 0.4 (the next echo would be at 0.8ms and so on.) I am after the decaying echo peaks that I drew arrows to the majority of them (they follow a decreasing trend. This trend will allow us to find T1 our relaxation time of when the field reverts back to its normal axis.

the equations have been different in each different paper i have read, however, y = A-B *exp(-tau/T1) is used the most. in the indicates the decaying of the echo as tau increases. The field should level out so I only need to really draw it out to tau 100 ms (or shorter if it levels off before then.)

I hope this helps clarify things a bit!

[Star Strider]

Viewing the figure —

F = openfig(websave('peaksmineraloil','https://www.mathworks.com/matlabcentral/answers/uploaded_files/1180283/peaksmineraloil.fig'));

I’ll edit my latest Comment to use only the tallest peak and then all peaks greater than 2 ...

.

[Alex]

Thank you so much, I played around with the time step and reviewed my lab notes as well. I needed to account for different time steps for each sample, but they can still be plotted on the same access as the signal decays after a time!

I appreciate all the help, I am relearning Matlab and this has helped me out quite a bit!

[Star Strider]

As always, my pleasure!