Re: [OctDev] imfilter works in the opposite way it's supposed to do

[Carnë Draug]

On 17 April 2012 07:56, Sunghyun Cho <address@hidden> wrote:
> 2012/4/16 Jordi Gutiérrez Hermoso <address@hidden>
>>
>> 2012/4/16 Carnë Draug <address@hidden>:
>> > On 16 April 2012 02:04, Sunghyun Cho <address@hidden> wrote:
>> >> Thanks for the reply, Carnë.
>> >> I tracked down the source of the bug, and found that conv2 works
>> >> incorrectly, when the shape parameter is 'valid'.
>> >> imfilter uses filter2, and filter2 calls conv2, which causes the
>> >> problem.
>> >> It looks like there is already a bug report about this.
>> >
>> > OK. Just for sake of completeness, I'll add that the bug report in
>> > question is #34893 and is marked to be fixed for the 3.8 releases
>> > https://savannah.gnu.org/bugs/?func=detailitem&item_id=34893
>>
>> There is still no fix for that bug. I tried to attack it, but it looks
>> like the bug is in the Fortran code, and I'm having trouble
>> understanding that code. If you care to help, it looks like the bug is
>> somewhere here in the *conv2.f functions:
>>
>>
>>  http://hg.savannah.gnu.org/hgweb/octave/file/a459d42bb0a6/libcruft/blas-xtra
>>
>> Thanks,
>> - Jordi G. H.
>
>
>
> I think I figured out how to fix it.
>
> The inner convolution routines have some bug.
> Following is a part of dconv2.f
>
> 65 integer ma,na,mb,nb
> 66 double precision a(ma,na),b(mb,nb)
> 67 double precision c(ma-mb+1,na-nb+1)
> 68 integer i,j,k
> 69 external daxpy
> 70 do k = 1,na-nb+1
> 71 do j = 1,nb
> 72 do i = 1,mb
> 73 call daxpy(ma-mb+1,b(i,j),a(mb+1-i,k+j-1),1,c(1,k),1)
> 74 end do
> 75 end do
> 76 end do
> 77 end subroutine
>
> Line 73 has a bug.
> The column indices of a and b moves in the same direction in the current
> implementation.
> It should be reversed. It might be fixed like below:
>
> call daxpy(ma-mb+1,b(i,j),a(mb+1-i,k+nb-j),1,c(1,k),1)
>
> I hope someone to be able to double-check this.

I have no knowledge of fortran but made this change for testing and
does not solve the problem. I still get the following

> conv2([1:5;1:5],[1:2],'valid')
ans =

    5    8   11   14
    5    8   11   14

Seems I should be getting [4 7 10 13; 4 7 10 13]

Also, this change adds 4 tests to make check (all in conv2, the
remaining 17 still succeed).

Carnë