Why does this program return complex numbers?

I am using this code

http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/doc/voicebox/distispf.html

to calculate the Itakura-Saito distance between two waveforms. When I used to it on several thousand different FFTs of the data it worked fine, however using it on the raw data produced results like "NaN -1.4918e+000i".

Looking at the code it seemed the NaNs were coming from any zeros in the input data as it divides by it. So I coded it to replace any zeros with 0.00001. This still leaves the complex numbers.

It's more down to application of theory then any coding errors but if someone could explain the reasoning behind it I'd be grateful.

Answer 1

[Babak]

I think the complex numbers come from the function log used in your code at line 0064.

Note that the input of log function should be positive real number or it will give out a complex output! Check to see if the argument of log positive or not. In other words, if

r=pf1(1:nx,:)./pf2(1:nx,:);

is always positive.