Using fsolve: how to plot fit curve/function to data?

조회 수: 9 (최근 30일)
Jamie Al
Jamie Al 2023년 2월 27일
이동: Torsten 2023년 2월 27일
I have some data that I am trying to show that it follows a function (I have the functional form of and solving for here using fsolve) and I want to make a plot of the data with the function fit. How can I do that? Thanks
%Load data from excel sheet
Data = readmatrix('Data.xlsx');
h = Data(:,1); %km height
qE = Data(:,2); %1/vol/s energy deposition
Zenith = 0; %deg solar zenith angle
%Find initial guess from data of star values for solver
qE_star_max_hat = 9192.04;
h_star_max_hat = 159.602;
H = 30.45; %scale height
%Evaluate Chapman prodction func at initial guess
for z = 1:length(h)
tao(z) = secd(Zenith) *exp(-(h(z) - h_star_max_hat)/H); %optical depth
ChapmanInitGuess_qE(z,1) = qE_star_max_hat*exp(1- (h(z) - h_star_max_hat)/H - tao(z)); %Energy
end
%Initial Guess vector
StarVars_hat = [qE_star_max_hat; h_star_max_hat; H];
disp('---------Fitted Values --------')
[qE_star_max, h_star_max, H] = GetChapmanFit(Zenith, StarVars_hat, h, qE);
Functions
%This function solves for the qE_star and h_star of a Chapman production
%function to fit a curve of real data
function [qE_star, h_star, H] = GetChapmanFit(Zenith, StarVars_hat, h, qE)
options = optimoptions('fmincon', 'FunctionTolerance', 1E-5, 'OptimalityTolerance', 1E-5, 'MaxIterations', 1000, 'display', 'none', 'Diagnostics', 'off');
ObjectiveFunction = @(StarVars) ChapmanObjectiveFunction(Zenith, StarVars, h, qE);
%
[solution,~] = fsolve(ObjectiveFunction, StarVars_hat, options);
qE_star = solution(1);
h_star = solution(2);
H = solution(3);
and
%This function is the objective of function which minimizes the difference
%b/w the actual qE data and the qE fitted from the Chapman production
%function
function [ErrorNorm] = ChapmanObjectiveFunction(Zenith, StarVars, h, qE)
qE_star_max = StarVars(1);
h_star_max = StarVars(2);
H = StarVars(3);
%Evaluate Chapman Production function at the initial Guess
for z = 1:length(h)
tao(z) = secd(Zenith) *exp(-(h(z) - h_star_max)/H); %optical depth
ChapmanFit_qE(z,1) = qE_star_max *exp(1- (h(z) - h_star_max)/H - tao(z)); %Energy
end
Error = abs(qE - ChapmanFit_qE);
ErrorNorm = norm(Error,2);
I also get some warning that I am not sure if it's a problem
MINCON options to FSOLVE. FSOLVE will use the common options and ignore the FMINCON options that do
not apply.

답변 (1개)

Torsten
Torsten 2023년 2월 27일
이동: Torsten 2023년 2월 27일
Call "fmincon", not "fsolve".
"fsolve" is not appropriate to solve your problem.
A better option than "fmincon" for your problem would be "lsqcurvefit".

카테고리

Help CenterFile Exchange에서 Interpolation에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by