matrix problems

[jonnysukes]

I'm trying to run a matlab script for a class I'm taking now, but I can't get
it to work. I already went through the script and changed all the syntax to
follow octave notations. The problem I'm having is trying to run the
operation:

C =D\f

I've gone through all the checks in the code and I found that the matrix D
has entries of +-Inf and so it can't run the operation. I'm almost positive
that's the problem because thats where the code stops running. But the
entries in D aren't actually Inf, they're just bigger than double precision
supports. Is there anyway around that problem, or an operation similar to \
that works on infinity? Here's the code as of now, I tried fixing the
problem by just replacing Inf with 2*2**1022, but it doesnt seem to be
working well. I also tried using the symbolic package, but it doesn't seem
to let you index matrices. Thanks for the help!

function nmse_error = refine_prony_hildebrand_error(ri, z); %(ri, t, f)

t = z{1}; 
                                
f = z{2};

N = length(t);
n = length(ri)/2;


delta = t(2) - t(1);
x = t/delta;

a = ri(1:n) + i*ri((n+1):end); 

% compute mu

mu = exp(a);

% construct matrix equation for series coefficients

D = ones(N,n);
bigconstant = 2*2**1022;
for r = 1:N
    for s = 1:n
       b = mu(s)^(r-1);
        if (abs(real(b)) == Inf && abs(imag(b)) == Inf)
           D(r,s) = sign(real(b))*bigconstant + sign(imag(b))*i*bigconstant;
        elseif(abs(real(b)) == Inf)
           D(r,s) = sign(real(b))*bigconstant + i*imag(b);
        elseif(abs(imag(b)) == Inf)
           D(r,s) = real(b) + sign(imag(b))*i*bigconstant;
        else
            D(r,s) = mu(s)^(r-1);
        end
    end
end

% solve system

C = D\f; 

% compute Prony fit

f_prony = zeros(size(t));
for k = 1:length(a)
    f_prony  = f_prony + exp(a(k)*x)*C(k);
end


% compute error

nmse_error = this_nmse(f, f_prony)

plot(t, real(f), t, real(f_prony))
drawnow

%
% computes the Normalized Mean Square Error between two vectors
%
% e = nmse(f,g)
%

function e = this_nmse(f,g)
e = sum(sum((abs(f-g)).^2)/sum((abs(f)).^2)); 
-- 
View this message in context: 
http://www.nabble.com/matrix-problems-tp15292917p15292917.html
Sent from the Octave - General mailing list archive at Nabble.com.