How to do Fminsearch- Curve fitting

21 views (last 30 days)
Anand Ra
Anand Ra on 29 Apr 2022
Commented: Anand Ra on 29 Apr 2022
I am trying to perfrom Fminsearch to fit the image intensity decaying data points to the wave propogation model , to determine constants
I am not sure how to go about the fminsearch. Kindly requesting help.
Below is my code:
% Objective:
%{
1. Import a CT scan image( CT scan of an Cu wire)
2. Extract the intensity from the CT scan image(CT scan of an Cu wire) :
which will be the observed data
3. Define the exponential decay model describing the Wave propagation through a medium
4. Perform fminsearch to determine the constants ( imaginary constants)
in the decay model.
%}
clc; % Clear the command window.
close all; % Close all figures (except those of imtool.)
clear; % Erase all existing variables. Or clearvars if you want.
workspace; % Make sure the workspace panel is showing.
format long g;
format compact;
%% Importing the CT scan/test image
% fullname = get_full_filename('20220422_5anglestep000.tif')
test_image = imread('test_image.tif');
figure(1)
imshow(test_image)
% Getting the dimensions of the image.
[rows, columns] = size(test_image);
%% % Computung and displaying the histogram.
figure(5)
[pixelCount, grayLevels] = imhist(test_image);
bar(pixelCount);
grid on;
title('Histogram of test image');
xlim([0 grayLevels(end)]); % Scale x axis manually.
%% Extracting the intensities from the medium
% Image represents the wave propgation through a medium in terms of
% intensity vs thickness(displacement) of the medium
% Looking to extract the intensity as a function of displacment from image
% k: displacement ( Z axis locations)
% wave_obs = intensity (Y intenisity)
figure(3)
middleRow = round(size(test_image, 1) / 2) % Find middle row.
wave_obs = test_image(middleRow, :) % Get all columns in the middle row.
% Now plot the intensity.
plot(wave_obs, 'b-', 'LineWidth', 2);
grid on;
xlabel('X Column');
ylabel('Gray Level');
%% Fitting wave_obs data points to the below defined exponential decay model
%% Defining the equation wave = exp(i(1-del)k*z * exp(-beta)*z
%{
i - signifying imaginary part
z - array representing the length of the medium : # of wave_obs data points
k - refractive index
del - constant (to be determined)
beta -constant (to be determined)
wave_obs- data points to fit
%}
k=1.5;
%% fminsearch : HELP NEEDED!
wave_pred = fminsearch(wave_prediction);
predicted_wave = wave_pred(z,k)
function wave_prediction = wave_pred(z,k)
del =5;
beta =7;
k=1.5;
z= 0:0.2:20;
wave1 = exp(-1i*(1-del)*k*z);
wave2 = exp(-beta)*z;
wave_prediction = @(del,beta)wave1.*wave2;
end
  3 Comments
Walter Roberson
Walter Roberson on 29 Apr 2022
wave_pred = fminsearch(wave_prediction);
At that point in the code, you do not have any function or variable named wave_prediction so it is not obvious what you are trying to do there.
Remember also that fminsearch() needs an initial search point.
predicted_wave = wave_pred(z,k)
I wonder if you should be calling that first and assigning the result to wave_prediction, since the output of the function is a function handle ?
wave_prediction = @(del,beta)wave1.*wave2;
That function handle accepts two parameters, del and beta, and ignores them both and always returns the product of wave1 and wave2. But if you are planning to make that the function handle used by fminsearch() then you have the difficulty that fminsearch() will only pass in a single parameter (that is a vector)
wave1 and wave2 are both defined in terms of the vector z, so the result of executing the function handle would be a vector. But fminsearch() needs a scalar to be returned.
You are trying to do a curve fitting, which requires that you have existing known data points and known corresponding value. It is not immediately clear to me which variables in your code hold the existing values ?
Anand Ra
Anand Ra on 29 Apr 2022
@Walter Roberson, Thanks, trying to understand your response. Btw the variable wave_obs is holding the data points to be fitted.

Sign in to comment.

Answers (1)

Matt J
Matt J on 29 Apr 2022
Edited: Matt J on 29 Apr 2022
Currently, you are minimizing your prediction function. You need to be minimizing, not the prediction, but a function that compares your prediction to measured data. You need to provide the measured data as well as an initial gues [del0,beta0] if the unknown parameters
fminsearch(@(p,meas) wave_prediction_Error(p, meas) , [del0,beta0]);
function Error = wave_prediction_Error(params, measurements)
del =params(1);
beta =(2);
k=1.5;
z= 0:0.2:20;
waveprediction = exp(-1i*(1-del)*k*z) .* exp(-beta)*z;
Error = norm( waveprediction - measurements);
end
  23 Comments
Matt J
Matt J on 29 Apr 2022
No, the figure shows that the decaying curves follow equations and which are not the model equation curves you have implemented.
Also, the curves decay starting from 1 at z=0. Your wave_obs data are all well above 1. Presumably, there was some sort of data normalization that you were supposed to do.
Also, the data you have attached are not monotonically decaying. They have oscillations:
Anand Ra
Anand Ra on 29 Apr 2022
@Matt J Thanks a lot for the support. Very much appreciated. I will look at the math again to see if I can simplify to work on the fitting.

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!