Release of G'MIC 3.0

@David_Tschumperle
It might be interesting to know that `ilaplacian` benefits from better accuracy provided by the math parser. Looking back at my old experiments, I thought I’d test a basic “seamless blend” with the divergence calculated as follows:

``````gcd_seamless_example :
sp sp
+r[0] 1,1,1,100%,2 avg={^} rm.
ri. .. to_rgba.
l. W={w} r2dx {round(0.25*w)},2 r2dx \$W,0,0,0.5,0.5 endl
split_opacity. n. 0,1 erode. 3
f[0] "begin(boundary=1;
fx0(x,y)=(j(x,y)-i);
fx1(x,y)=(j(#1,x,y)-i#1);
bx0(x,y)=(i-j(x,y));
bx1(x,y)=(i#1-j(#1,x,y));
);
lerp(fx0(1,0),fx1(1,0),i#2) -
lerp(bx0(-1,0),bx1(-1,0),j(#2,-1)) +
lerp(fx0(0,1),fx1(0,1),i#2) -
lerp(bx0(0,-1),bx1(0,-1),j(#2,0,-1))"
k[0] ilaplacian 0 +fc. \$avg +[-2,-1] c 0,255
``````

I don’t know how it affects speed, but the image quality does improve I think.

Interesting. I’ve always been a bit disappointed by `ilaplacian`, as it sometimes renders weird things. It’s actually a bit hard to manage the boundary conditions correctly, and this command probably requires some update.
Anyway, I’ve not tested (yet) your example, but I will !

Indeed, I’m beginning to wonder where the improvement actually lies… my assumption about precision could be nothing to do with it. Perhaps worth looking at regardless

Edit: yes, it looks like boundaries are the culprit. Using `expand_xy` with neumann boundary makes problems disappear as well. I see `blend_seamless` already deals with it. Sorry for the distraction!

Update:
I get it now. The reason is boundary conditions for the second gradient/difference. The first gradient should be neumann, but the second should be dirichlet. This illustrates the problem:
`gmic 5,1,1,1,x +g x,1 g. x,-1 f[0] "begin(boundary=1);j(1)-i" f[0] "begin(boundary=0);i-j(-1)"`

Doing it all in a single pass within the math parser avoids the issue:
`gmic 5,1,1,1,x f "begin(boundary=1);j(1)-i-(i-j(-1))"`

Conclusion: if `gradient` allowed boundary = dirichlet then the problem would be solved.
I imagine you already know this though

1 Like

Seems I was still wrong. The boundary needs to be periodic, for both first and second differences:
`gmic 5,1,1,1,x f "begin(boundary=2);j(1)-i" f "begin(boundary=2);i-j(-1)" ilaplacian 0`

That suggests `gradient` should be recoded to allow this. And probably as a custom command rather than a native one (I don’t see why a custom command would be less performant in that case, after all, this is just about applying some 3x3 kernels)

Yes, I think ideally it should allow all the usual boundaries. Periodic boundary laplacian actually solves the seamless problem completely. No glitches. That can be seen by taking `gcd_seamless_example` and set boundary=2.

Proof:
`gmic sp +r 1,1,1,100%,2 avg={^} rm. f "begin(boundary=2);j(1)+j(-1)+j(0,1)+j(0,-1)-4*i" ilaplacian 0 +fc \\$avg add`

That’s probably because `ilaplacian` uses a FFT-based reconstruction, which implicitely considers periodic boundary conditions.

I’ll try to rewrite `gradient` as a custom command, to allow any kind of boundary conditions.

1 Like

Thanks, both of you.

1 Like

I think I’m done with it.

It was not that easy, considering the various schemes `gradient` implement, as well as the fact that it can be invoked without arguments.
For those interested, that’s the current source code:

``````#@cli gradient : { x | y | z | c }...{ x | y | z | c },_scheme,_boundary_conditions : (no arg)
#@cli : Compute the gradient components (first derivatives) of selected images, along specified axes.
#@cli : (eq. to 'g').\n
#@cli : 'scheme' can be { -1=backward | 0=centered | 1=forward | 2=sobel | 3=rotation-invariant (default) | \
# 4=deriche | 5=vanvliet }.
#@cli : 'boundary_conditions' can be { 0=dirichlet | 1=neumann | 2=periodic | 3=mirror }.
#@cli : (no arg) compute all significant components.
#@cli : Default values: 'scheme=0' and 'boundary_conditions=1'.
v + _gradient \$axes,\$scheme,\$boundary v -
if \$is_noarg noarg fi

# Parse and check arguments.
is_noarg=0
l[]
if ['\$1']==0 axes,scheme,boundary=,0,1 is_noarg=1 # No 1st arg
else
is_axes={"s=['\$1'];fill(s,k,isin(s[k],'xyzc'));min(s)"}
if !\$is_axes axes,scheme,boundary=,0,1 is_noarg=1 # Invalid 1st arg
elif ['\$2']==0 axes,scheme,boundary=\$1,0,1 # No 2nd arg
else
is_scheme={"isint(\$2) && inrange(\$2,-1,5)"}
if !\$is_scheme axes,scheme,boundary=,0,1 is_noarg=1 # Invalid 2nd arg
elif ['\$3']==0 axes,scheme,boundary=\$1,\$2,1 # No 3rd arg
else
is_boundary={"isint(\$3) && inrange(\$3,0,3)"}
if !\$is_boundary axes,scheme,boundary=,0,1 is_noarg=1 # Invalid 3rd arg
else axes,scheme,boundary=\${1-3}
fi
fi
fi
fi
onfail axes,scheme,boundary=,0,1 is_noarg=1 # Invalid arguments
endl
u \$axes,\$scheme,\$boundary,\$is_noarg

# Main command.
axes,scheme,boundary=\${1=},\${2=},\${3=}

# Display log message.
if '\$axes'!=0 s_axes=" along axes '"\$axes"'" fi
s0,s1,s2,s3,s4,s5,s6=backward,centered,forward,sobel,rotation-invariant,deriche,vanvliet
s_scheme=\${s{\$scheme+1}}
s0,s1,s2,s3=dirichlet,neumann,periodic,mirror
s_boundary=\${s\$boundary}
e[0--3] "Compute gradient of image\$?"\$s_axes", with "\$s_scheme" scheme and "\$s_boundary" boundary conditions."
sx,sy,sz,sc=",",;,/,^

# Process images.
repeat \$! l[\$<] nm={n}
laxes=\$axes if '\$laxes'==0 if d>1 laxes=xyz else laxes=xy fi fi
('\$laxes')
repeat w
a={`i[#1,\$>]`} s=\${s\$a}
elif \$scheme==-1 # Backward
(-1\${s}1\${s}0) +correlate[0] .,\$boundary rm..
elif \$scheme==0 # Centered
(-0.5\${s}0\${s}0.5) +correlate[0] .,\$boundary rm..
elif \$scheme==1 # Forward
(0\${s}-1\${s}1) +correlate[0] .,\$boundary rm..
elif \$scheme==2 # Sobel
if _'\$a'==_'x' (-1,0,1;-2,0,2;-1,0,1)
elif _'\$a'==_'y' (-1,-2,-1;0,0,0;1,2,1)
else (-0.5\${s}0\${s}0.5)
fi
+correlate[0] .,\$boundary rm..
elif \$scheme==3 # Rotation-invariant
A,B={[0.25*(2-sqrt(2)),0.5*(sqrt(2)-1)]}
if _'\$a'==_'x' (-\$A,0,\$A;-\$B,0,\$B;-\$A,0,\$A)
elif _'\$a'==_'y' (-\$A,-\$B,-\$A;0,0,0;\$A,\$B,\$A)
else (-0.5\${s}0\${s}0.5)
fi
+correlate[0] .,\$boundary rm..
else # Deriche/Vanvliet
s4,s5=deriche,vanvliet com=\${s\$scheme}
if \$boundary<2 +\$com[0] 0,1,\$a,\$boundary
else
if _'\$a'==_'x' +r[0] {0,[w+2,h,d,s]},0,\$boundary,0.5 \$com. 0,1,\$a shrink_x. 1
elif _'\$a'==_'y' +r[0] {0,[w,h+2,d,s]},0,\$boundary,0,0.5 \$com. 0,1,\$a shrink_y. 1
elif _'\$a'==_'z' +r[0] {0,[w,h,d+2,s]},0,\$boundary,0,0,0.5 \$com. 0,1,\$a shrink_z. 1
fi
fi
fi
done
rm[0,1]
nm \$nm
endl done
u \$is_noaarg
``````

I’ve tried disabling the native command `gradient` so that the custom one can be used instead. So far, it seemed to work well.
I’ll build pre-release binaries with all that included.

1 Like

Also for the improvements :

• `gradient` can be called with `axis=c` now.
• And of course, all types of boundary conditions are managed now.

And cherry on the cake : this also made me discover a small bug in the `convolve`/`correlate` operators that led to a segfault when used kernel was mai,ly along Z-axis.

1 Like

That’s a good piece of work! I’m not surprised it took effort, having done my own basic equivalents using `convolve`. It’ll be tomorrow before I can test it myself, but I guess problems ought to show up quickly given how widely used `gradient` is.

I’m wondering now about the `laplacian` command possibly supporting periodic boundary… not important though, it should be simple enough to make with the new `gradient`!

Yes, and `hessian` too, that should be also converted as a custom command as well (should not be that hard, now that `gradient` has been done).

1 Like

Very good David, thank you.

Maybe you could allow deriche and vanvliet for “c” too with boundary >=2. Occasionally I use the spectral channel for time sequences.

That’s not possible no, because the algorithm does not allow it.
Deriche and Vanvliet filters use recursive filtering that limits the possible boundary conditions to neumann and dirichlet.
I could fake it by resizing the image first, then apply deriche or vanvliet with dirichlet boundaries, then crop the result. But as those are native commands, implementing an original algorithm that does not allow these bc, I’m not really sure this is a good idea
(you can easily write a custom command for that of course ).

Ok, I should only remind me, that it is not possible resp. I have to fake it by myself. I see in your implementation of _gradient, that you already do the faking by resize for “xyz” in the deriche and vanvliet case!

In fact it is a rare case having 3d temporal multispectral data.

Can we have a new command that is image, but with linear opacity? As in opacity adds rather than using alpha composition? It would solve my problem with Fragment Blur.

Choose `-value` as the opacity when drawing `J` on `I`, to do `I + J*value`

1 Like

Yes, but in that case, I know that the `sigma` is small enough so that I only need to add 1px-border to the images.
In the general case of `deriche` and `vanvliet` , the parameter `sigma` can be anything, and it’s kinda hard to determine the size of the needed borders. Actually, as the recursive filtering implies kernel with theoretically infinite support, the amount of borders to add could be really huge!

Ok, understood.
Thanks again, for all your work too.

• `hessian` has been converted to a custom command too, so with any kind of `boundary_conditions`.
• `sharpen` has been also converted (but without the choice for boundary conditions in this case, as it seems not really interesting).