> 
img1 := Read(cat(kernelopts(homedir),"/Documents/circshift/13ydvt0.jpg")): img2 := Read(cat(kernelopts(homedir),"/Documents/circshift/30rkxvt.jpg")):

> 
# Scale img1 to dimensions of img2 # (padding with zero seems not appropriate to this circular shift # and absolutedifference scheme. img1 := Scale(img1,1..184,1..184):

> 
A:=Reshape(img1,[1..184*184]):

> 
B:=Reshape(img2,[1..184*184]):

> 
p:=proc(N::integer, lo::integer, hi::integer, k::integer, a::Array(datatype=float[8],order=C_order), b::Array(datatype=float[8],order=C_order), c::Array(datatype=float[8],order=C_order)) local cbest::float,tot::float,i::integer,kps::integer, s::integer,sbest::integer,Nsq::integer; sbest:=0; cbest:=1.0e6; Nsq:=N*N; for s from lo to hi do tot:=0.0: for i from 1 to Nsq do tot := tot+abs(b[i]a[irem(Nsq1+is,Nsq)+1]); end do; c[k+s]:=tot; if c[k+s]<cbest then sbest:=s; cbest:=c[k+s]; end if; end do; return sbest; end proc: cp:=Compiler:Compile(p):

> 
C:=Array(1..401,datatype=float[8],order=C_order): SBest:=CodeTools:Usage( trunc(evalhf(p(184,200,200,201,A,B,C))) ); C[201+SBest];

memory used=0.92KiB, alloc change=0 bytes, cpu time=4.90s, real time=4.91s, gc time=0ns
> 
# The Compiled version is about 10 times faster than the evalhf'ed version C:=Array(1..401,datatype=float[8],order=C_order): SBest:=CodeTools:Usage( cp(184,200,200,201,A,B,C) ); C[201+SBest];

memory used=0.84KiB, alloc change=0 bytes, cpu time=468.00ms, real time=469.00ms, gc time=0ns
> 
M:=184^2/2; C:=Array(1..2*M+1,datatype=float[8],order=C_order): SBest:=CodeTools:Usage( cp(184,M,M,M+1,A,B,C) ); C[M+1+SBest];

memory used=0.81KiB, alloc change=0 bytes, cpu time=39.73s, real time=39.80s, gc time=0ns
> 
final:=CircularShift(A,SBest): imgfinal:=Reshape(final,[184,184]):

> 
View([img1,img2,imgfinal]); #Embed([img1,img2,imgfinal]);

> 
plots:pointplot(<<seq(i,i=M..M)>Vector[column](C)>, gridlines=false);

