Hi,
I am trying to implement the 1D DMUSIC algorithm described in:
[1] Y. Li, J. Razavilar, and K. Liu, "A Super-Resolution Parameter
Estimation Algorithm for Multi- Dimensional NMR Spectroscopy," University
of Maryland, College Park, MD 1995.
I am not sure what is wrong with my MATLAB code. Does anyone know if there
is an implementation available for download? My MATLAB code is:
function [f s] = dmusic(y,K,J)
% function dmusic(y,K,J)
%
% 1D Damped MUSIC (DMUSIC) Algorithm
%
% y = signal
% K = # of damped sinusoids
% J = Subsample of prediction matrix
%
% Ex:
% dmusic(y,2,100);
%
% Reference:
% [1] Y. Li, J. Razavilar, and K. Liu, "A Super-Resolution Parameter
% Estimation Algorithm for Multi- Dimensional NMR Spectroscopy,"
% University of Maryland, College Park, MD 1995.
%
N = length(y)-1;
if((N-J) < K || J < K)
error('J has to be between K and N-K');
end
% Prediction Matrix
for ind1=1:J
i = ind1-1;
endPoint = N-J-1+i;
A(1:J,i+1) = y(i+1:endPoint+1);
end
% Singular Value Decomposition
[U D V] = svd(A);
% MUSIC Algorithm
Vn = 0;
for k=K+1:J
Vn = Vn + conj(V(:,k))*(V(:,k).');
end
indB = 1;
indW = 1;
for B=0:0.05:0.5
for w=0:0.01:2*pi
s(indB,indW) = B + j*w;
r = 0;
for i=0:J-1
n=i+1;
r(n,1) = exp(i*s(indB,indW));
end
rt = r/norm(r);
f(indB,indW) = 1/((rt')*Vn*rt);
indW = indW + 1;
end
indW = 1;
indB = indB + 1;
end
contour(imag(s),real(s),abs(f));
xlabel('Frequency');
ylabel('Amplitude');
Damped MUSIC (DMUSIC) Algorithm
Started by ●January 7, 2008
Reply by ●January 7, 20082008-01-07
On 7 Jan, 12:41, "jfaller30" <jfalle...@hotmail.com> wrote:> Hi, > > I am trying to implement the 1D DMUSIC algorithm described in: > > [1] � � Y. Li, J. Razavilar, and K. Liu, "A Super-Resolution Parameter > Estimation Algorithm for Multi- Dimensional NMR Spectroscopy," University > of Maryland, College Park, MD 1995.Not the easiest reference to get hold of... do you know if a copy is available on-line?> I am not sure what is wrong with my MATLAB code.You need to test your code with ideal signals. Start with one damped sinusoidal with known parameters, and no noise. The function should return the same parameters that you know to be present. Then escalate to two known sinusoidals with known parameters. Once you are happy with that case, start adding noise. And so on. Rune
Reply by ●January 7, 20082008-01-07
>On 7 Jan, 12:41, "jfaller30" <jfalle...@hotmail.com> wrote: >> Hi, >> >> I am trying to implement the 1D DMUSIC algorithm described in: >> >> [1] =A0 =A0 Y. Li, J. Razavilar, and K. Liu, "A Super-ResolutionParameter=> >> Estimation Algorithm for Multi- Dimensional NMR Spectroscopy,"University>> of Maryland, College Park, MD 1995. > >Not the easiest reference to get hold of... do you know if >a copy is available on-line? > >> I am not sure what is wrong with my MATLAB code. > >You need to test your code with ideal signals. Start >with one damped sinusoidal with known parameters, and >no noise. The function should return the same parameters >that you know to be present. Then escalate to two known >sinusoidals with known parameters. Once you are happy >with that case, start adding noise. > >And so on. > >Rune >The URL to get a copy of the reference is: https://drum.umd.edu/dspace/handle/1903/5631. The reference has some examples at the end. I implemented example 1: function dmusicExp() N = 24; J = 12; K = 2; SNR = 40; w = wgn(N+1,1,SNR); theta = 0.1; s1 = -0.2+j*2*pi*0.42; s2 = -0.2+j*2*pi*(0.42+theta); y = 0; for i=1:N+1 n=i-1; y(i) = exp(s1*n) + exp(s2*n) + w(i); % y(i) = exp(s1*n) + exp(s2*n); end [f s] = dmusic(y,K,J); John
Reply by ●January 7, 20082008-01-07
>>On 7 Jan, 12:41, "jfaller30" <jfalle...@hotmail.com> wrote: >>> Hi, >>> >>> I am trying to implement the 1D DMUSIC algorithm described in: >>> >>> [1] =A0 =A0 Y. Li, J. Razavilar, and K. Liu, "A Super-Resolution >Parameter= >> >>> Estimation Algorithm for Multi- Dimensional NMR Spectroscopy," >University >>> of Maryland, College Park, MD 1995. >> >>Not the easiest reference to get hold of... do you know if >>a copy is available on-line? >> >>> I am not sure what is wrong with my MATLAB code. >> >>You need to test your code with ideal signals. Start >>with one damped sinusoidal with known parameters, and >>no noise. The function should return the same parameters >>that you know to be present. Then escalate to two known >>sinusoidals with known parameters. Once you are happy >>with that case, start adding noise. >> >>And so on. >> >>Rune >> > >The URL to get a copy of the reference is: > >https://drum.umd.edu/dspace/handle/1903/5631. > >The reference has some examples at the end. I implemented example 1: > >function dmusicExp() > N = 24; > J = 12; > K = 2; > SNR = 40; > w = wgn(N+1,1,SNR); > theta = 0.1; > s1 = -0.2+j*2*pi*0.42; > s2 = -0.2+j*2*pi*(0.42+theta); > > y = 0; > for i=1:N+1 > n=i-1; > y(i) = exp(s1*n) + exp(s2*n) + w(i); >% y(i) = exp(s1*n) + exp(s2*n); > end > > [f s] = dmusic(y,K,J); > >John >BTW Rune, I got the noiseless case to work. I am having trouble with the spectrum aspect equation 19. Another reference that you can check out is: http://www.emo.org.tr/resimler/ekler/46f5c776b38ad4b_ek.pdf The following is my code for the noiseless case: function [f s] = dmusicExp() N = 24; J = 12; K = 2; SNR = 40; w = wgn(N+1,1,SNR); theta = 0.1; s1 = -0.2+j*2*pi*0.42; s2 = -0.2+j*2*pi*(0.42+theta); y = 0; for i=1:N+1 n=i-1; y(i) = exp(s1*n) + exp(s2*n); end [f s] = dmusicIdeal(y,K,J); function [finS s] = dmusicIdeal(y,K,J) N = length(y)-1; if((N-J) < K || J < K) error('J has to be between K and N-K'); end % Prediction Matrix for ind1=1:J i = ind1-1; endPoint = N-J-1+i; A(1:J,i+1) = y(i+1:endPoint+1); end % Singular Value Decomposition [U D V] = svd(A); Vn = V(:,K+1:J); % Ideal, if no noise indSi = 1; indR = (0:J-1)'; si = 0; for B=[0 0.2] for w=[0 2*pi*0.42 2*pi*0.52 2*pi] si(indSi) = -B + j*w; indSi = indSi + 1; end end s = nchoosek(si,2); indFinS = 1; finS = zeros(1,2); for indS=1:size(s,1) for indK=1:K r(:,indK) = exp(indR.*s(indS,indK)); end if(roundn(norm((Vn.')*r),-4) == 0) finS = s(indS,:); end end John
Reply by ●January 7, 20082008-01-07
Rune, I got it to work!!! The following code is example 1 and the DMUSIC algorithm from reference [1]. All you have to do is run dmusicEx1(). [1] Y. Li, J. Razavilar, and K. Liu, "A Super-Resolution Parameter Estimation Algorithm for Multi- Dimensional NMR Spectroscopy," University of Maryland, College Park, MD 1995 https://drum.umd.edu/dspace/handle/1903/5631 ---------------------------------------------------- MATLAB Code: function dmusicEx1() % function dmusicEx1() % % 1D Damped MUSIC (DMUSIC) Algorithm Example % % This is the implementation of example 1 from reference [1] % % Ex: % dmusicEx1(); % % References: % [1] Y. Li, J. Razavilar, and K. Liu, "A Super-Resolution Parameter % Estimation Algorithm for Multi- Dimensional NMR Spectroscopy," % University of Maryland, College Park, MD 1995. % [2] E. Yılmaz and E. Dilaveroğlu, "A Performance Analysis for DMUSIC, % "presented at 5TH International Conference on Electrical and % Electronics Engineering (ELECO 2007), Bursa, Turkey, 2007. % % Coded by: Kenneth John Faller II January 07, 2008 % N = 24; J = 12; K = 2; SNR = 40; theta = 0.1; s1 = -0.2+j*2*pi*0.42; s2 = -0.2+j*2*pi*(0.42+theta); y = 0; for i=1:N+1 n=i-1; y(i) = exp(s1*n) + exp(s2*n); end y = awgn(y,SNR); dmusic(y,K,J); function dmusic(y,K,J) % function dmusic(y,K,J) % % 1D Damped MUSIC (DMUSIC) Algorithm % % y = signal % K = # of damped sinusoids % J = Subsample of prediction matrix % % Ex: % dmusic(y,2,100); % % References: % [1] Y. Li, J. Razavilar, and K. Liu, "A Super-Resolution Parameter % Estimation Algorithm for Multi- Dimensional NMR Spectroscopy," % University of Maryland, College Park, MD 1995. % [2] E. Yılmaz and E. Dilaveroğlu, "A Performance Analysis for DMUSIC, % "presented at 5TH International Conference on Electrical and % Electronics Engineering (ELECO 2007), Bursa, Turkey, 2007. % % Coded by: Kenneth John Faller II January 07, 2008 % N = length(y)-1; if((N-J) < K || J < K) error('J has to be between K and N-K'); end % Prediction Matrix for ind1=1:J i = ind1-1; endPoint = N-J-1+i; A(1:J,i+1) = y(i+1:endPoint+1); end % Singular Value Decomposition [U D V] = svd(A); Vn = V(:,K+1:J); % If noise is present indB = 1; indW = 1; for B=0:0.05:0.5 for w=0:0.01:2*pi s(indB,indW) = -B + j*w; r = 0; for indR=0:J-1 n=indR+1; r(n,1) = exp(indR*s(indB,indW)); end rt = r/norm(r); f(indB,indW) = real(1/((rt')*(conj(Vn)*(Vn.'))*rt)); indW = indW + 1; end indW = 1; indB = indB + 1; end contour(imag(s),real(s),f,20); xlabel('Frequency'); ylabel('Damping Factor');
Reply by ●January 7, 20082008-01-07
Rune Allnor <allnor@tele.ntnu.no> writes:> On 7 Jan, 12:41, "jfaller30" <jfalle...@hotmail.com> wrote: > > Hi, > > > > I am trying to implement the 1D DMUSIC algorithm described in: > > > > [1] � � Y. Li, J. Razavilar, and K. Liu, "A Super-Resolution Parameter > > Estimation Algorithm for Multi- Dimensional NMR Spectroscopy," University > > of Maryland, College Park, MD 1995. > > Not the easiest reference to get hold of... do you know if > a copy is available on-line?Citeseer has it: http://citeseer.ist.psu.edu/221906.html Get it here: http://citeseer.ist.psu.edu/cache/papers/cs/6848/ftp:zSzzSzdspserv.eng.umd.eduzSzpaperzSzliyezSzJMR.pdf/a-super-resolution-parameter.pdf Ciao, Peter K. -- "And he sees the vision splendid of the sunlit plains extended And at night the wondrous glory of the everlasting stars."
