1c Copyright (C) 2010-2021 The Octave Project Developers
2c
3c See the file COPYRIGHT.md in the top-level directory of this
4c distribution or <https://octave.org/copyright/>.
5c
6c This file is part of Octave.
7c
8c Octave is free software: you can redistribute it and/or modify it
9c under the terms of the GNU General Public License as published by
10c the Free Software Foundation, either version 3 of the License, or
11c (at your option) any later version.
12c
13c Octave is distributed in the hope that it will be useful, but
14c WITHOUT ANY WARRANTY; without even the implied warranty of
15c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
16c GNU General Public License for more details.
17c
18c You should have received a copy of the GNU General Public License
19c along with Octave; see the file COPYING.  If not, see
20c <https://www.gnu.org/licenses/>.
21c
22      subroutine zdconv2o(ma,na,a,mb,nb,b,c)
23c purpose:      a 2-dimensional outer additive convolution.
24c               equivalent to the following:
25c                 for i = 1:ma
26c                   for j = 1:na
27c                     c(i:i+mb-1,j:j+mb-1) += a(i,j)*b
28c                   endfor
29c                 endfor
30c arguments:
31c ma,na (in)    dimensions of a
32c a (in)        1st matrix
33c mb,nb (in)    dimensions of b
34c b (in)        2nd matrix
35c c (inout)     accumulator matrix, size (ma+mb-1, na+nb-1)
36c
37      integer ma,na,mb,nb
38      double complex a(ma,na)
39      double precision b(mb,nb)
40      double complex c(ma+mb-1,na+nb-1)
41      double complex btmp
42      integer i,j,k
43      external zaxpy
44      do k = 1,na
45        do j = 1,nb
46          do i = 1,mb
47            btmp = b(i,j)
48            call zaxpy(ma,btmp,a(1,k),1,c(i,j+k-1),1)
49          end do
50        end do
51      end do
52      end subroutine
53
54      subroutine zdconv2i(ma,na,a,mb,nb,b,c)
55c purpose:      a 2-dimensional inner additive convolution.
56c               equivalent to the following:
57c                 for i = 1:ma-mb+1
58c                   for j = 1:na-nb+1
59c                     c(i,j) = sum (sum (a(i:i+mb-1,j:j+nb-1) .* b))
60c                   endfor
61c                 endfor
62c arguments:
63c ma,na (in)    dimensions of a
64c a (in)        1st matrix
65c mb,nb (in)    dimensions of b
66c b (in)        2nd matrix
67c c (inout)     accumulator matrix, size (ma+mb-1, na+nb-1)
68c
69      integer ma,na,mb,nb
70      double complex a(ma,na)
71      double precision b(mb,nb)
72      double complex c(ma-mb+1,na-nb+1)
73      double complex btmp
74      integer i,j,k
75      external zaxpy
76      do k = 1,na-nb+1
77        do j = 1,nb
78          do i = 1,mb
79            btmp = b(i,j)
80            call zaxpy(ma-mb+1,btmp,a(mb+1-i,k+nb-j),1,c(1,k),1)
81          end do
82        end do
83      end do
84      end subroutine
85