zigzag scan: Image processing function

[Muthiah Annamalai-3]

Hello there,
I have a function zigzag() which walks-off a (square)matrix and reads its elements in a zigzag fashion. It is useful for image processing.

Example:
octave  --eval "x=[1:4]'*[1:4];eval('x');tic;printf('%d\n',zigzagscan(x));toc()" -q
x =

    1    2    3    4
    2    4    6    8
    3    6    9   12
    4    8   12   16

1
2
2
3
4
3
4
6
6
4
8
9
8
12
12
16
Elapsed time is 0.177223 seconds.


However I dont know where I should commit it to in Octave-forge.
Someone please suggest the right place.

Code is attached.

Cheers
Muthu


[zigzag.m]

---

## Copyright (C) 2006, June, Muthiah Annamalai, <muthiah.annamalai@...>
##
## This program is free software; you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation; either version 2 of the License, or
## (at your option) any later version.
##
## This program is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
## GNU General Public License for more details.
##
## You should have received a copy of the GNU General Public License
## along with this program; if not, write to the Free Software
## Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
##

## -*- texinfo -*-
## @deftypefn {Function File} {} zigzag (@var{M})
##
## Returns the zig-zag scanned coefficients of given matrix @var{M}.
## Function of ZigZag scan is to accumulate the transform (DCT like)
## coefficients in the LF->HF in that order, from a square matrix.
##
## An exmaple of the use of @code{zigzag} is
## @example
## @group
##   x=[1:4]'*[1:4];
##   zigzag(x)
## @end group
## @end example
## @end deftypefn

% (C) 2006 June, Muthiah Annamalai, <muthiah.annamlai@...>
%
% For a simple case_ we can use a ZZ scan using many loops, to
% generate the indices.
% eg: case_ 4x4 we have ZZ scan indices being (in that order)
% (each of the rows have sums of 2:8) i.e 2:2S)
% (each row of t=2:2S has sum= (t-1) elements if_ t<S+2, and  sum=2S+1-t, elsewhere)
% (each row contains all the permutations summing up to t, in that order)
% (permutations are equidistant)
% (permutations are ordered as oddrows have descending in  x)
% (even rows have ascending in x)
% (if_ even row & non-multiple of 3, increase y index, else_ x index \
%   for_ first element of row.)
% (all odd rows have to increase x component first)
% (if_ only with first element, we can get the rest of the elements)
% (odd rows only have odd no of elements, even rows have only even no
%   of elements)
%
%1 (1,1),2
%2 (1,2) (2,1),3,x:y,y:x
%3 (3,1),(2,2),(1,3),4
%4 (1,4),(2,3),(3,2),(4,1),5
%5 (4,2),(3,3),(2,4),6
%6 (3,4),(4,3),7
%7 (4,4),8
%
%
%1 (1,1),2
%2 (1,2) (2,1),3,x:y,y:x
%3 (3,1),(2,2),(1,3),4
%4 (1,4),(2,3),(3,2),(4,1),5
%5 (5,1),(4,2),(3,3),(2,4),(1,5),6
%
%6 (3,4),(4,3),7
%7 (4,4),8
%
function rval=zigzag(mtrx)

  if (nargin < 1)
    error('usage: zigzag_opt(square-matrix)')
  else
    S=issquare(mtrx);
    if(S <= 0)
      error('usage: zigzag_opt(square-matrix)')
    end
  end
  x=1;
  y=1;

  itrx=1;
  rval=zeros(1,S*S);
  rval(itrx)=[mtrx(x,y)];
  iseven=true;
 
  for itr=2:2*S-1
    %[x,y]
    %itr
    iseven=~iseven;

    %predicting the first element
    %is the trick!
    if(iseven)
      if(itr > S)
        x=x+1;
        y=y;
      else
        x=x;
        y=y+1;
      end
    else
      %odd row
      if(itr <= S)
        y=y;
        x=x+1;
      else
        y=y+1;
        x=x;
      end
    end
    %[x,y]

    %compute nelem for_ this row.
    if(itr <S+1)
      nelem=itr;
    else
      nelem=2*S-itr;
    end
    rval(itrx+1)=mtrx(x,y);

    %compute target sum
    tgt_sum=itr;    
   
    %generate the elements
    if(iseven)
      %even rows ascend in X,
      %and have only even no of elements
      for itr2=2:nelem
        x=x+1;
        y=y-1;
        %[x,y]
        rval(itrx+itr2)=mtrx(x,y);
      end
    else
      %odd rows ascend in Y,
      for itr2=2:nelem
        y=y+1;
        x=x-1;
        %[x,y]
        rval(itrx+itr2)=mtrx(x,y);
      end
    end
    itrx=itrx+itr2;
  end
end

---

_______________________________________________
Octave-sources mailing list
Octave-sources@...
https://www.cae.wisc.edu/mailman/listinfo/octave-sources

[Fredrik Bülow-2]

Hi Muthiah!

I took the liberty of writing my own quicker version of zigzag.  I've
also added support for mxn matrices (where m!=n).

Code is attached.

Cheers
Fredrik Bülow

On 10/5/06, Muthiah Annamalai <muthuspost@...> wrote:

>
>  Hello there,
>  I have a function zigzag() which walks-off a (square)matrix and reads its
> elements in a zigzag fashion. It is useful for image processing.
>
>  Example:
>  octave  --eval
> "x=[1:4]'*[1:4];eval('x');tic;printf('%d\n',zigzagscan(x));toc()"
> -q
>  x =
>
>      1    2    3    4
>      2    4    6    8
>      3    6    9   12
>      4    8   12   16
>
>  1
>  2
>  2
>  3
>  4
>  3
>  4
>  6
>  6
>  4
>  8
>  9
>  8
>  12
>  12
>  16
>  Elapsed time is 0.177223 seconds.
>
>
>  However I dont know where I should commit it to in Octave-forge.
>  Someone please suggest the right place.
>
>  Code is attached.
>
>  Cheers
>  Muthu
>
>
> _______________________________________________
> Octave-sources mailing list
> Octave-sources@...
> https://www.cae.wisc.edu/mailman/listinfo/octave-sources
>
>
>
>

[my_zigzag.m]

---

## Copyright (C) 2006, October, Fredrik Bulow, <fredrik.bulow@...>
##
## This program is free software; you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation; either version 2 of the License, or
## (at your option) any later version.
##
## This program is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
## GNU General Public License for more details.
##
## You should have received a copy of the GNU General Public License
## along with this program; if not, write to the Free Software
## Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
##



function rval = my_zigzag(mtrx)
  n=size(mtrx);

  if(issquare(mtrx)) #Square matrix (quick case)

    ##We create a matrix of the same size as mtrx where odd elements are
    ##1, others 0.
    odd=kron(ones(ceil(n/2)),eye(2))((1:n(1)),(1:n(2)));

    ##We transpose even elements only.
    mtrx = mtrx.*odd + (mtrx.*(1-odd))';

    ##Now we mirror the matrix. The desired vector is now the
    ##concatenation of the diagonals.
    mtrx=mtrx(:,1+size(mtrx,2)-(1:size(mtrx,2)));

    ##Picking out the diagonals.
    rval  = [];
    for i = n(2)-1:-1:1-n(1)
      rval=[rval diag(mtrx,i)'];
    endfor

  else #Not square (Slow cases)
    mtrx=mtrx(:,1+size(mtrx,2)-(1:size(mtrx,2)));

    ##Picking out the diagonals and reversing odd ones manually.
    rval  = [];
    for i = n(2)-1:-1:1-n(1)
      new = diag(mtrx,i);
      if(floor(i/2)==i/2) ##Even?
        rval=[rval new'];
      else                ##Odd!
        rval=[rval new((1+length(new))-(1:length(new)))'];
      endif
    endfor
  endif
endfunction

---

_______________________________________________
Octave-sources mailing list
Octave-sources@...
https://www.cae.wisc.edu/mailman/listinfo/octave-sources

[Muthiah Annamalai-3]

Hello Fredrik,
I think your code does zigzag on a transposed matrix than the one I am
speaking about (I have in my zigzag function).

Your code does a zigzag of matrix = reshape(1:9,3,3) as follows,
##   1   4- 7
##   | /  /
##   2   5  8
##     /   /|
##   3 - 6  9
but I would expect the output to be as [1 4 2 3 5 7 8 6 9]; instead.
This is the expected behaviour of a DCT frequency coeffients walk-order.

I have also added some documentation to your code. Please see that,
and also adjust the code.

One way  I did adjust your code was to add one extra line

 mtrx=mtrx'

before the if statement. Your code however is not modified, as I figured
you may have another way to do the same without a transpose.

-Muthu


On Fri, 2006-10-06 at 16:40 +0200, Fredrik Bülow wrote:

> Hi Muthiah!
>
> I took the liberty of writing my own quicker version of zigzag.  I've
> also added support for mxn matrices (where m!=n).
>
> Code is attached.
>
> Cheers
> Fredrik Bülow
>
> On 10/5/06, Muthiah Annamalai <muthuspost@...> wrote:
> >
> >  Hello there,
> >  I have a function zigzag() which walks-off a (square)matrix and reads its
> > elements in a zigzag fashion. It is useful for image processing.
> >
> >  Example:
> >  octave  --eval
> > "x=[1:4]'*[1:4];eval('x');tic;printf('%d\n',zigzagscan(x));toc()"
> > -q
> >  x =
> >
> >      1    2    3    4
> >      2    4    6    8
> >      3    6    9   12
> >      4    8   12   16
> >
> >  1
> >  2
> >  2
> >  3
> >  4
> >  3
> >  4
> >  6
> >  6
> >  4
> >  8
> >  9
> >  8
> >  12
> >  12
> >  16
> >  Elapsed time is 0.177223 seconds.
> >
> >
> >  However I dont know where I should commit it to in Octave-forge.
> >  Someone please suggest the right place.
> >
> >  Code is attached.
> >
> >  Cheers
> >  Muthu
> >
> >
> > _______________________________________________
> > Octave-sources mailing list
> > Octave-sources@...
> > https://www.cae.wisc.edu/mailman/listinfo/octave-sources
> >
> >
> >
> >

[zigzag.m]

---

## Copyright (C) 2006, October, Fredrik Bulow, <fredrik.bulow@...>
##
## This program is free software; you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation; either version 2 of the License, or
## (at your option) any later version.
##
## This program is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
## GNU General Public License for more details.
##
## You should have received a copy of the GNU General Public License
## along with this program; if not, write to the Free Software
## Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
##

## -*- texinfo -*-
## @deftypefn {Function File} {} zigzag (@var{mtrx})
## Returns zigzag walk-off of the elements of @var{mtrx}.
## Essentially it walks the matrix in a Z-fashion.
##  
## mat =
##   1   4   7
##   2   5   8
##   3   6   9
## then zigzag(mat) gives the output,
## [1   2   4   7   5   3   6   8   9], by walking as
## shown in the figure from pt 1 in that order of output.
## The argument @var{mtrx} should be a matrix
##
## @example
## @group
## mat = reshape(1:9,3,3);
## zigzag(mat)
## ans =[1   2   4   7   5   3   6   8   9]
##
## @end group
## @end example
##
## @end deftypefn
##

## Author:   Fredrik Bulow, <fredrik.bulow@...>

function rval = zigzag(mtrx)
  if nargin < 1
    error('usage: zigzag(matrix); see help zigzag');
  end
  n=size(mtrx);
 
  if(issquare(mtrx)) #Square matrix (quick case)

    ##We create a matrix of the same size as mtrx where odd elements are
    ##1, others 0.
    odd=kron(ones(ceil(n/2)),eye(2))((1:n(1)),(1:n(2)));

    ##We transpose even elements only.
    mtrx = mtrx.*odd + (mtrx.*(1-odd))';

    ##Now we mirror the matrix. The desired vector is now the
    ##concatenation of the diagonals.
    mtrx=mtrx(:,1+size(mtrx,2)-(1:size(mtrx,2)));

    ##Picking out the diagonals.
    rval  = [];
    for i = n(2)-1:-1:1-n(1)
      rval=[rval diag(mtrx,i)'];
    endfor

  else #Not square (Slow cases)
    mtrx=mtrx(:,1+size(mtrx,2)-(1:size(mtrx,2)));

    ##Picking out the diagonals and reversing odd ones manually.
    rval  = [];
    for i = n(2)-1:-1:1-n(1)
      new = diag(mtrx,i);
      if(floor(i/2)==i/2) ##Even?
        rval=[rval new'];
      else                ##Odd!
        rval=[rval new((1+length(new))-(1:length(new)))'];
      endif
    endfor
  endif
endfunction

---

_______________________________________________
Octave-sources mailing list
Octave-sources@...
https://www.cae.wisc.edu/mailman/listinfo/octave-sources

[Fredrik Bülow-2]

Hi again!

I've done the following.


1. Changed the "argument check" to nargin != 1. I was thinking about
   putting the whole function in the else part of this if clause. I
   know that the error statements exits the function so I guess one
   might as well keep it the way it is. I don't know what the standard
   way of doing these things within the octave community is.

2. It's perfectly fine to implement the new ordering with mtrx=mtrx'
   but in order to show off I did it by changing odd and even operations
   around instead. Since both your and my zigzag behave the same way
   I also renamed the new function to zagzig.

3. Cleaned the square case a bit by replacing n=size(mtrx) with
   n=length(mtrx).

One thing I thought of... If you have lots of matrices A1, A2, A3 and
they all are of the same dimensions, then zagzig can be used "once and
for all" with massive CPU savings. First make a "index mapping vector"
like this:

n=size(A1);
i=zagzig(reshape(1:prod(n),n(1),n(2)));

Thereafter the zagzigs are given by

A1(i), A2(i) and A3(i).

I suspect that this is about an order of magnitude faster than doing
zagzig(A) for each matrix individually.

The perfect way of doing zagzig would be finding out a quick way of
generating the i vector without for-loops and if-statements all
together. I suspect that this can be done easily but I can't figure
out how.

Code attached.

Cheers
Fredrik

On 10/7/06, Muthiah Annamalai <muthuspost@...> wrote:

> Hello Fredrik,
> I think your code does zigzag on a transposed matrix than the one I am
> speaking about (I have in my zigzag function).
>
> Your code does a zigzag of matrix = reshape(1:9,3,3) as follows,
> ##   1   4- 7
> ##   | /  /
> ##   2   5  8
> ##     /   /|
> ##   3 - 6  9
> but I would expect the output to be as [1 4 2 3 5 7 8 6 9]; instead.
> This is the expected behaviour of a DCT frequency coeffients walk-order.
>
> I have also added some documentation to your code. Please see that,
> and also adjust the code.
>
> One way  I did adjust your code was to add one extra line
>
>  mtrx=mtrx'
>
> before the if statement. Your code however is not modified, as I figured
> you may have another way to do the same without a transpose.
>
> -Muthu
>
>
> On Fri, 2006-10-06 at 16:40 +0200, Fredrik Bülow wrote:
> > Hi Muthiah!
> >
> > I took the liberty of writing my own quicker version of zigzag.  I've
> > also added support for mxn matrices (where m!=n).
> >
> > Code is attached.
> >
> > Cheers
> > Fredrik Bülow
> >
> > On 10/5/06, Muthiah Annamalai <muthuspost@...> wrote:
> > >
> > >  Hello there,
> > >  I have a function zigzag() which walks-off a (square)matrix and reads its
> > > elements in a zigzag fashion. It is useful for image processing.
> > >
> > >  Example:
> > >  octave  --eval
> > > "x=[1:4]'*[1:4];eval('x');tic;printf('%d\n',zigzagscan(x));toc()"
> > > -q
> > >  x =
> > >
> > >      1    2    3    4
> > >      2    4    6    8
> > >      3    6    9   12
> > >      4    8   12   16
> > >
> > >  1
> > >  2
> > >  2
> > >  3
> > >  4
> > >  3
> > >  4
> > >  6
> > >  6
> > >  4
> > >  8
> > >  9
> > >  8
> > >  12
> > >  12
> > >  16
> > >  Elapsed time is 0.177223 seconds.
> > >
> > >
> > >  However I dont know where I should commit it to in Octave-forge.
> > >  Someone please suggest the right place.
> > >
> > >  Code is attached.
> > >
> > >  Cheers
> > >  Muthu
> > >
> > >
> > > _______________________________________________
> > > Octave-sources mailing list
> > > Octave-sources@...
> > > https://www.cae.wisc.edu/mailman/listinfo/octave-sources
> > >
> > >
> > >
> > >
>
>
>

[zagzig.m]

---

## Copyright (C) 2006, October, Fredrik Bulow, <fredrik.bulow@...>
##
## This program is free software; you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation; either version 2 of the License, or
## (at your option) any later version.
##
## This program is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
## GNU General Public License for more details.
##
## You should have received a copy of the GNU General Public License
## along with this program; if not, write to the Free Software
## Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
##

## -*- texinfo -*-
## @deftypefn {Function File} {} zigzag (@var{mtrx})
## Returns zigzag walk-off of the elements of @var{mtrx}.
## Essentially it walks the matrix in a Z-fashion.
##  
## mat =
##   1   4   7
##   2   5   8
##   3   6   9
## then zigzag(mat) gives the output,
## [1   2   4   7   5   3   6   8   9], by walking as
## shown in the figure from pt 1 in that order of output.
## The argument @var{mtrx} should be a matrix
##
## @example
## @group
## mat = reshape(1:9,3,3);
## zigzag(mat)
## ans =[1   2   4   7   5   3   6   8   9]
##
## @end group
## @end example
##
## @end deftypefn
##

## Author:   Fredrik Bulow, <fredrik.bulow@...>

function rval = zagzig(mtrx)

  if nargin != 1 #Checking arguments.
    error('usage: zigzag(matrix); see help zigzag');
  endif

  if issquare(mtrx) #Square matrix (quick case)
    n=length(mtrx);
    ##We create a matrix of the same size as mtrx where odd elements are
    ##1, others 0.
    odd=kron(ones(n,n),eye(2))((1:n),(1:n));

    ##We transpose even elements only.
    mtrx = (mtrx.*odd)' + (mtrx.*(1-odd));

    ##Now we mirror the matrix. The desired vector is now the
    ##concatenation of the diagonals.
    mtrx=mtrx(:,1+size(mtrx,2)-(1:size(mtrx,2)));

    ##Picking out the diagonals.
    rval  = [];
    for i = n-1:-1:1-n
      rval=[rval diag(mtrx,i)'];
    endfor

  else #Not square (Slow cases)
    n=size(mtrx);
    mtrx=mtrx(:,1+size(mtrx,2)-(1:size(mtrx,2)));

    ##Picking out the diagonals and reversing odd ones manually.
    rval  = [];
    for i = n(2)-1:-1:1-n(1)
      new = diag(mtrx,i);
      if floor(i/2)==i/2 ##Even?
        rval=[rval new((1+length(new))-(1:length(new)))'];
      else                ##Odd!
        rval=[rval new'];
      endif
    endfor
  endif
endfunction

---

_______________________________________________
Octave-sources mailing list
Octave-sources@...
https://www.cae.wisc.edu/mailman/listinfo/octave-sources