[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: I NEED HELP: Differences between MATLAB and Octave
From: |
Paul Kienzle |
Subject: |
Re: I NEED HELP: Differences between MATLAB and Octave |
Date: |
Tue, 18 Jun 2002 02:24:59 -0400 |
On Fri, Jun 07, 2002 at 02:25:17PM -0500, Marco Boni wrote:
> The simple code I've attached to this mail runs very fast under Matlab,
> while it is executed very slow under Octave.
>
> Why? Where I've failed!
For-loops in octave are embarrassingly slow. You have nested
forloops (one in your code, and one in hist). Both are probably
removable with a bit of cleverness.
For example, in hist:
for i = 1:n-1
cutoff (i) = (x (i) + x (i + 1)) / 2;
endfor
is simply translated into
cutoff = ( x(1:n-1) + x(2:n) ) / 2;
The next loop in hist.m is more tricky:
freq(1) = sum(y < cutoff (1));
for i = 2:n-1
freq (i) = sum (y >= cutoff (i-1) & y < cutoff (i));
endfor
freq(n) = sum (y >= cutoff (n-1));
You can get part way there with sort:
[s, idx] = sort ( [cutoff(:); y(:)] );
pt = [0; idx > n; 0];
This gives you a sequence of 0's and 1's where the zeros represent
bin boundaries and the ones represent bin contents. Because sort
is 'stable' all y values which are identical to an x value will
appear after the x value, so the semi-open [x-1, x) logic will
be preserved. The leading and trailing zeros capture all the y's
which are not otherwise in a bin.
Using cumulative sum, you can turn the ones into frequency counts
chist = cumsum(pt);
But you only need the counts at the bin boundaries:
chist = s(find(diff(chist) == 0));
Now you have a 'cumulative histogram'. Differentiate it and you
should have the histogram:
freq = diff(chist);
I may messed up something along the way, but hopefully you get
the idea. Could you fix it, test it and submit it to bug-octave?
Since you are dealing with really large numbers of bins (up to 2^20 if
I'm reading your code correctly), fixing hist should address most of
the speed problems.
You can also remove a number of statements outside of your for loop,
but that probably won't gain you very much.
Paul Kienzle
address@hidden
>
>
>
> --
> Marco Boni
> Viale delle Magnolie, 4
> 50142 - Firenze
> Cell.: +39 335 6079353
> Home : +39 055 7398033
> Email: address@hidden
> -----------------------------------------------------
> My machine is an i686 and running Red Hat Linux 7.2
> Why don't you try it?
> -----------------------------------------------------
>
>
> %entropia di un segnale random
>
> closeplot;
> clear all;
> k=16;
> nt=2^k;
> a=floor(rand(1,nt)+.5);
> dt=k+4;
> for n=1:dt
> dimb=n;
> numb=floor(nt/dimb);
> dim=dimb*numb;
> aa=reshape(a(1:dim),dimb,numb);
> x=(1:1:n)-1;
> xx=2.^x;
> ris1=xx*aa+1;
> ris2=hist(ris1,2^n);
> ris3=ris2(find(ris2));
> %
> size(ris3,2)/size(ris2,2);
> n
> p=ris3/numb;
> ris(n)=-sum(p.*log(p));
> end;
> entr=diff(ris);
> entr(1)=log(2);
> plot(entr);hold;plot(entr,'*');
> plot(log(2)*ones(1,dt),'r');grid
> input ("Pick a number, any number! ")
-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.org
How to fund new projects: http://www.octave.org/funding.html
Subscription information: http://www.octave.org/archive.html
-------------------------------------------------------------