Deconvolution results a single value

33 vues (au cours des 30 derniers jours)
lakshitha Lathpandura
lakshitha Lathpandura le 6 Fév 2023
Commenté : John D'Errico le 7 Fév 2023
Hi There!
I am implementing a MATLAB code for the deconvolution of a spectrum. Initially, I am testing my code on the data provided in a published paper.
According to the text mentioned above, by deconvoluting the Fermi function (represented by the dash line, 'Fermi_func' in the MATLAB code), from the spectrum (represented by the black open circles, 'Spectra' in the code), I should obtain the instrumental function, which is represented by a solid line and is a Gaussian function with a width of 0.55 eV.
I digitized the spectrum in the paper (data.txt file attached) and wrote a MATLAB code for its deconvolution (deconvolution.m file attached). However, I am now getting a single value instead of an array for the instrumental function. I was expecting a Gaussian function. I'm not sure what's going wrong.
I am new to MATLAB and any advice regarding this issue would be highly appreciated.
Thank You.
Lakshitha
===============================
Implemented MATLAB code is given below (the .m file is also attached):
clear all;
close all;
clc
A = readmatrix('Data.txt');
Energy = A(:,1);
Spectrum = A(:,2);
Fermi_func = 1-(1./(1+exp((Energy-8.18)/0.025852)));
Fermi_func(1)= 1e-6; % I got an error from deconv that "First coefficient must be non-zero". So I replaced it with a very small value.
Instrumental_Function = deconv(Spectrum,Fermi_func);
% Instrumental_Function =
% deconv(Counts(find(Counts,1):end),Fermi_func(find(Fermi_func,1):end)) %
% This was also implemeted for removing the error of "First coefficient must be non-zero".
display(Instrumental_Function)
f1 = figure('Name','Inverse Photoemission Spectroscopy Fitting');
h1 = plot(Energy, Spectrum,'ko');
hold on;
xlabel('Energy (eV)')
ylabel('Counts (arb. units)')
h2 = plot(Energy, Fermi_func,'k--');
% h3 = plot(Energy, Instrumental_Function,'k-');
hold off;
savefig(f1,'IPES_Fitting.fig')

Réponse acceptée

John D'Errico
John D'Errico le 6 Fév 2023
Convolution of two vectors yields a new vector that is longer than either of them. So if the vectors have length N and M, then the convolved result will have length N+M-1. For example,
X = rand(1,10);
Y = [1 2 1]/4;
Z = conv(X,Y);
numel(Z)
ans = 12
Deconvolution does the opposite. It produces a shorter vector.
Xhat = deconv(Z,Y);
plot(X,'-rx')
hold on
plot(Xhat,'-bo')
As you can see, the deconvolution yields a result that is shorter than the vector Z. The deconvolved result of deconv(Z,Y) will have length
length(Z) - length(Y) + 1
ans = 10
This is as expected. However, IF you perform a deconvolution, where the two vectors have the same lengths, then the result will have length 1. For example...
deconv(randn(1,10),randn(1,10))
ans = -0.0278
So it will be a scalar result. And that makes complete sense. Did you not get a scalar result in what you did?
  4 commentaires
lakshitha Lathpandura
lakshitha Lathpandura le 7 Fév 2023
Modifié(e) : lakshitha Lathpandura le 7 Fév 2023
Thank you very much, John and Walter, for your feedback regarding my post. I understand now that deconvoluting two functions with the same length results in a single value. John, your explanation was very helpful. Based on your suggestion, I worked for the last 24 hours (without a good sleep :) ), and I came to the conclusion that the 'deconv' function is not going to help me since deconvoluting a function from experimental noise data does not result in a valid function.
So, I started using the 'conv' function to solve it in the other direction, and it works pretty well. However, now I am facing a new issue with the 'lsqcurvefit' function for the convoluted function. I have posted a new question in the link below, and you can see how things went after your suggestions as well.
Thank You.
John D'Errico
John D'Errico le 7 Fév 2023
Deconvolution is a difficult problem. It can be ill-posed, which suggests it will be difficult to solve.

Connectez-vous pour commenter.

Plus de réponses (0)

Catégories

En savoir plus sur Historical Contests dans Help Center et File Exchange

Produits


Version

R2022b

Community Treasure Hunt

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

Start Hunting!

Translated by