How to plot an inflection point from a temperature profile data source?

Question

Robin Patel 2020-12-13

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/692485-how-to-plot-an-inflection-point-from-a-temperature-profile-data-source

评论： Mathieu NOE 2020-12-15

final.mat

Hello everyone,

I have been trying to plot a tangent to a curve. My data source is a temperature profile from a TCL (temperature control lab). It has two heaters and the data from the arduino Leonardo is collected every second. I have a set of 2000 data points from the curve. While plotting the inflection point, I am not able to get the tangent to the curve. I tried to smooth the curve but no proper tangent could be found. I have the file attached. Can anyone please help me?

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Mathieu NOE 2020-12-13

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/692485-how-to-plot-an-inflection-point-from-a-temperature-profile-data-source#answer_574340

在 MATLAB Online 中打开

myslidingavg.m

hello

here you are

enjoy it !

% sliding avg method
temp_out = myslidingavg(temp, 50);
figure(1),
semilogx(time,temp,'b',time,temp_out,'r','linewidth',2);
% first derivative
dx = mean(diff(time));
dy = [0; diff(temp_out)./dx];
dy_out = myslidingavg(dy, 25);
figure(2),
plot(time,dy,'b',time,dy_out,'r','linewidth',2);
% point of inflection = fin peak of first derivative (but not the first
% sample transient)
ind = find(time>time(1)+5 & time<time(1)+200);
[peak,loc] = max(dy_out(ind));
time_inflection = time(ind(loc));
y_inflection = temp_out(ind(loc));
dy_inflection = dy_out(ind(loc));
% tangent equation
temp_tang = y_inflection+dy_inflection*(time-time_inflection);
figure(3),
plot(time,temp,'b',time,temp_out,'r',time(ind),temp_tang(ind),'--k','linewidth',2);
legend('raw temp','smoothed','tangent at inlection point');

4 个评论
显示 2更早的评论隐藏 2更早的评论

Mathieu NOE 2020-12-14

在 MATLAB Online 中打开

hello

this is a new code

I used the second part of the curve to compute the time constant

it is based on the slope measured at the inflection point, and divided by the delta of temperature (final - initial)

I am not sure what you want to do with the tangent crossing the y axis (at y = 0)

load('data.mat')
% sliding avg method
temp_out = myslidingavg(temp, 25);
figure(1),
semilogx(time,temp,'b',time,temp_out,'r','linewidth',2);
% first derivative
dx = mean(diff(time));
dy = [0; diff(temp_out)./dx];
dy_out = myslidingavg(dy, 25);
figure(2),
plot(time,dy,'b',time,dy_out,'r','linewidth',2);
% point of inflection = fin peak of first derivative (but not the first
% sample transient)
ind = find(time>time(1)+500 & time<time(1)+700);
[peak,loc] = max(dy_out(ind));
time_inflection = time(ind(loc));
y_inflection = temp_out(ind(loc));
dy_inflection = dy_out(ind(loc));
% tangent equation
temp_tang = y_inflection+dy_inflection*(time-time_inflection);
% find time point where tangent crosses the y axis (y = 0)
% 0 = y_inflection+dy_inflection*(time_zero-time_inflection);
time_zero = round(time_inflection - y_inflection/dy_inflection);
ind2 = find(time>time_zero & time<time(1)+700);
figure(3),
plot(time,temp,'b',time,temp_out,'r',time(ind2),temp_tang(ind2),'--k','linewidth',2);
legend('raw temp','smoothed','tangent at inlection point');
%% time constant computation
% NB : input step amplitude is unknown but we can scale the data 
% second step has a temperature amplitude of 61° - 36.5° = 24.5 °
% so time constant = computed slope (dy) at inflection point divided by temperature amplitude 
Tc = dy_inflection/(61 - 36.5)      % unit : 1/s

Robin Patel 2020-12-15

Thank you for your response once again!

I checked the code and it seems pretty well to what I wanted. I tried with new data file as data.mat and the following code with Code.txt file but I can not figure out how can I let the tangent to the x-axis. I think there is some small mistake which I cannot figure out. Can you please help me?

Thank you for the help in advance.

Mathieu NOE 2020-12-15

在 MATLAB Online 中打开

hello

so this is a modified code

I noticed that your time vector was not evenly sampled (diff(time) fluctuates with values between 0 and 1)

so first thing is resampling

then I add a bit of smoothing

and I prefered not to use the second derivative (data must be super smooth to use second derivative), but simply use the fact that the inflection point is the peak of the first derivative

code below :

clc
clear all
load data.mat
h1 = diff(temp,2);
K = temp(end);
% resample data because of uneven time samples (sometimes diff(time)
% contains zeros)
dt = 1;
t = min(time):dt:max(time);
temp = interp1(time,temp,t);
temp = myslidingavg(temp, 10);% smoothing
time = t;
%%%%%
L_index = find(temp>=0.1*K,1);
L = time(L_index);
T_index = find(time>=(1-exp(-1))*K,1);
T = time(T_index);
D = diff(temp)./dt;
D = myslidingavg(D, 10); % smoothing
figure,plot(D)
% inflex = find((diff(D)/dt)<0.5,1);
[peak,inflex] = max(D); % inflection point defined as max value of first derivative
% return
A = D(inflex)*time(inflex)-temp(inflex);
tangent = D(inflex)*time - A;
plot(time,temp);
hold on;
plot(time,tangent);
hold on;
plot(L,temp(L_index),'*');
plot(T,temp(T_index),'o')
plot(time,temp);
hold off;
xlim([0 max(time)])
ylim([30 max(temp)])