Re: Oct-file version of x>0

[Keith Goodman]

On 5/23/06, Dmitri A. Sergatskov <address@hidden> wrote:
> On 5/23/06, Keith Goodman <address@hidden> wrote:
> > While porting a m-file to an oct-file I got stuck on the line
> >
> > y = x'*(x > 0);
> >
> > where x is Nx1.
> >
> > How do I port x>0?
> >
>
> I cannot help you with your oct file, but the expression gets
> parsed correctly on octave 2.9.5:
>
> octave:1> x=rand(10,1)
> x =
>
>   0.323223
>   0.372934
>   0.021102
>   0.876365
>   0.498951
>   0.183033
>   0.863289
>   0.186314
>   0.469338
>   0.622294
>
> octave:2> y = x'*(x > 0)
> y = 4.4168

Plus it is faster than sum(x.*(x > 0)) or sum(x(x > 0)).

Since I don't know how to do x>0 in an oct-file, I resorted to a trick:

x'*(x > 0) = sum(x + abs(x))/2

#include <oct.h>
DEFUN_DLD(abc, args, nargout,"y = x'*(x > 0) where x is Nx1")
{
 octave_value retval;
 const Matrix x = args(0).matrix_value();
 return retval = (x + x.abs()).column_sum()/2.0;
}

I also tried brute force:

#include <oct.h>
DEFUN_DLD(abc, args, nargout,"y = x'*(x > 0) where x is Nx1")
{
 octave_value retval;
 const Matrix x = args(0).matrix_value();
 const int n = x.rows();
 double s = 0.0;
 for (int i=0; i < n; i++) {
   s = s + (x(i) > 0)*x(i);
 }
 return retval = s;
}