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 :slight_smile:

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 :slight_smile:

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:

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'.
#@cli : $ image.jpg gradient
#@cli : $$ https://gmic.eu/oldtutorial/_gradient
gradient :
  _gmic_s="$?" axes,scheme,boundary,is_noarg=${"_gradient_get_args $*"}
  v + _gradient $axes,$scheme,$boundary v -
  if $is_noarg noarg fi

# Parse and check arguments.
_gradient_get_args : skip "${1=},${2=},${3=}"
  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.
_gradient :
  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}
      if ${gradient_$a} [gradient_$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
      nm. gradient_$a
    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 :slight_smile: ).

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).