Inversion of trunk dielectric constant from dihedral scattering mechanism ratio (alpha) using Symbolic math toolbox

Question

I'm trying to perform inversion of soil and trunk dielectric constant (es and et respectively) from dihedral scattering mechanism ratio (alpha) using symbolic math toolbox. Initially I've done polarimetric decomposition of a quad pol SAR dataset and alpha and local incidence angle of soil (theta_s) and vegetation trunk (theta_t) were derived from there. I've developed a code to perform the inversion but the result turns out to be wrong. The optimal inversion of dielectric constant is done by taking the value of dielectric constant at min(αm - αD) where αm corresponds to the dihedral scattering mechanism ratio from the model and αD is obtained from the decomposition of the SAR image. There are non-linear equations for this calculations which are given in the code.

`% Read GeoTIFF information
[im, R] = geotiffread('path to the tif file');
data = im2double(im);
% Extract relevant bands
alpha = data(:, :, x);
theta_s = data(:, :, y);
% Calculate theta_t for each pixel
theta_t = 90 - theta_s;
% Define symbolic variables
syms es et theta_s_ij theta_t_ij alpha_ij real

% Define symbolic expressions for Fresnel scattering coefficients
Rsh = ((cosd(theta_s_ij)-sqrt(es-(sind(theta_s_ij)^2)))/(cosd(theta_s_ij)+sqrt((es-sind(theta_s_ij)^2))));
Rth = ((cosd(theta_t_ij)-sqrt(et-(sind(theta_t_ij)^2)))/(cosd(theta_t_ij)+sqrt((et-sind(theta_t_ij)^2))));
Rsv = ((es*cosd(theta_s_ij)-sqrt(es-(sind(theta_s_ij)^2)))/(es*cosd(theta_s_ij)+sqrt((es-sind(theta_s_ij)^2))));
Rtv = ((et*cosd(theta_t_ij)-sqrt(et-(sind(theta_t_ij)^2)))/(et*cosd(theta_t_ij)+sqrt((et-sind(theta_t_ij)^2))));
alpha_m = real((Rsh * Rth - Rsv * Rtv) / (Rsh * Rth + Rsv * Rtv));

% Define objective function
alpha_diff = alpha_m - alpha_ij;

% Convert the objective function to a MATLAB function handle
alpha_diff_func = matlabFunction(alpha_diff, 'vars', {es, et, theta_s_ij, theta_t_ij, alpha_ij});

% Initialize empty arrays for results
es_out = zeros(size(alpha_ij));
et_out = zeros(size(alpha_ij));

% Loop through each pixel of the dummy matrix
for i = 1:size(alpha, 1)
    for j = 1:size(alpha, 2)
        % Extract current pixel value of alpha
        alpha_ij = alpha(i, j);
        theta_s_ij = theta_s(i, j);
        theta_t_ij = 90 - theta_s_ij; % Calculate theta_t for each pixel

        % Find optimal es and et using fminsearch within range 0-40
        es_opt = fminsearch(@(x) alpha_diff_func(x, 0, theta_s_ij, theta_t_ij, alpha_ij), 20); % Start    with initial guess of 10
        et_opt = fminsearch(@(x) alpha_diff_func(0, x, theta_s_ij, theta_t_ij, alpha_ij), 20);

        % Ensure es and et are within 0-40 range
        es_out(i, j) = max(min(es_opt, 40), 0);
        et_out(i, j) = max(min(et_opt, 40), 0);
    end
end

% Create new GeoTIFF with the same metadata as input
geotiffwrite('path to output', es_out, R);
geotiffwrite('path to output', et_out, R);`

I tried to obtain dielectric constants of soil and vegetation trunk using the above code. But the values are either the lowest value or the highest one. There is some problem with the logic of the code. The decomposed outputs consists of various bands out of which one is dihedral scattering mechanism ratio (alpha) and one is local incidence angle of soil (theta_s). From theta_s, theta_t can be calculated and these values goes as input to the equations. The optimal dielectric constant would be the one corresponding to minimum difference between decomposed alpha and modelled alpha (alpha_m). what is the error in this code and how can I modify the code to get proper result.