Best practices for writing efficient user commands

Morning @afre!

Once upon a time, this was THE book to read (and I believe that it still is very valuable):

https://en.wikipedia.org/wiki/The_Unix_Programming_Environment

Have fun!
Claes in Lund, Sweden

Here, smoothness is optional because of the underscore. However, this presents a small inconsistency where

gmic sp lena my_command

works but

gmic sp lena my_command b 3

doesn’t; i.e., my_command must be followed by a comma.



How about best practices for efficient math? But I might be getting ahead of myself since I am still trying to grasp the content in [Mathematical expressions](http://gmic.eu/reference.shtml#section9).

I don’t see inconsistency here. After all, b could be a valid argument for my_command. Arguments of commands are not typed, so the G’MIC interpreter knows only there are arguments or not.
If arguments are expected, then forcing the user to put at least a comma seems to be logical.

Anyway, in my_command, you could additionally check if the provided argument seems to be correct, and if not, use the noarg command to tell no arguments have been provided.
Like this:

my_command : skip ${1=10}
   smoothness=$1
   if {!isval($1)} noarg smoothness=10 fi
   blur $smoothness

and in that case, a call like $ gmic sample lena my_command mirror x will apply my_command with default smoothness value set to 10.
This is not the ultimate solution to manage default parameters though, as the interpreter will try to subtitute a given argument of my_command anyway, with possible side effects giving unexpected results, as in:

foo : 
  sp lena
  my_command x=${"-echo foo u booh"}

that will display foo twice, as the item following my_command will be evaluated twice.

I read this more in the context of ‘how to do it quickly’. There was discussion a long time ago about various tricks from contributors, some of it made it to the wiki pages. Some of it likely doesn’t apply by now anyway so I wouldn’t heed it much.

The math parser (e.g. in the fill command) used to be avoided but it’s quite good now, which also removes two other problems; using loops over an image (like repeat… done) isn’t required and resorting to ‘abusing’ the 3D commands is no longer needed. If you’re writing filters from scratch and you really need speed I find sticking to ‘native’ commands (i.e. only running C++ code not G’MIC script) is best. Even then you can get a lot done by mixing that with the math parser. With larger files using shared buffers can make quite a large difference due to not having to copy data around.

This as well, partly motivated by my old system, which chokes even at native commands.


I am slowly depending more on native commands and the math parser, and aware of shared buffers. However, this depends on my knowledge and experience in G’MIC, math and programming, all of which is lacking :slight_smile:.

It might be best if I gave you an example. I have been working on a command called deets. What would you do to improve it? (If anyone decides to use deets, consider it to be under the license that G’MIC user commands are typically with credit to @afre.)

deets:
  r={[im,iM]} n 0,1 split_details 5
  +l[^0] n 0,1
    repeat $! l[$>]
      s$>={log10(ia#$>/sqrt(iv#$>))}
    endl done
    rm
  endl
  m={min($s0,$s1,$s2,$s3)}
  *[1] {$s0/$m} *[2] {$s1/$m} *[3] {$s2/$m} *[4] {$s3/$m}
  + n $r
1 Like

Without knowing much about what this is doing (something to do with SNR I guess given average/stddev?), the main thing I’d try to do is remove the ‘+l[^0]’ part. It’s the most expensive part other than split_details because it copies all the split images. Here’s one way without that (I added a parameter for number of levels):

deets : skip ${1=5}
  r={[im,iM]} split_details $1
  l[^0]
    $!,1,1,1,"(ia#x-im#x)/sqrt(iv#x)" log. /. {im}
    repeat {$!-1} -*[$>] {@$>} done rm.
  endl
  + n $r

The first step was to rearrange the maths to avoid having to normalize to (0,1). It turns out that’s just a case of subtracting image minimum from the average, because this is a ratio and variance is not affected by the minimum.

Then instead of using separate variables to store the values you can use a 5x1 ‘image’ which is very like a 1 dimension array. That way when fill loops over each point you can reference the relevant image using the x coord.

Finally now that we have an array the log part, minimum and division can all be done to the array as a whole. Then it’s a case of looping over the images to multiply by the relevant array value: {@$>}.

I’m really not sure what you’re doing maths wise here but ratios of logs amount to taking a log of the numerator with a base of the denominator. What I mean is you have log10 (a) / log10 (b) and that’s the same as logb (a) which also makes the base irrelevant in the fraction.

Incidentally I really don’t like usual log() syntax as if it were a function, I tend to treat it as an operator in my scribbles. I find it easier to write x|10 instead of log10 (x). After all it’s only repeated division in a way!

1 Like

@garagecoder Thanks, I will see what I can glean from your explanation. The command itself isn’t important. I just find myself looping and storing info all of the time, and wanted to see what you thought about how I structured it.

By the way I think you’re doing very well indeed, it took me an incredibly long time to get anywhere with G’MIC and I still feel like I know 10%. The main thing is having something that does what you want and it’s great for making something quickly. Don’t be discouraged! :slight_smile:

You’re doing fine @afre, indeed. First step is always to make something work as expected :slight_smile:
Then, it is just cosmetic.
Here is my proposal, from the version of @garagecoder, except I’m not using an image to store the log numbers, but a single variable s. There is almost no differences, I suspect @garagecoder’s version is even slightly faster by the way :slight_smile:

deets2 :
  r={[im,iM]} split_details 5
  l[^0]
    s={"s = vectorl(); for (k = 0, k<l, ++k, s[k] = log((ia#k - im#k)/sqrt(iv#k))); s/min(s)"}
    repeat $! *[$>] {arg(1+$>,$s)} done
  endl
  + n $r
1 Like

Also just for fun, this is how the command could be done using “only” code written in the math parser :

deets :
  eval "
    range = [im,iM];
    ext('split_details 5');
    s = vector4();
    for (k = 1, k<5, ++k, s[k-1] = log((ia#k - im#k)/sqrt(iv#k)));
    s/=min(s);
    for (k = 1, k<5, ++k, ext('*[',vtos(k),'] ',vtos(s[k-1])));
    ext('+ n ',vtos(range));
  "

Thanks for the feedback. I enjoy experimenting with decomposition :clown_face:. I have more questions but I might not understand your answers :stuck_out_tongue:.

About split_details

  • What do base_scale and detail_scale do? Should I just stick with a trous wavelets?
  • Are there quicker ways to do (almost) the same thing? I might do crazy things like nesting.

Bonus questions

  • How do I select a region of an image on which to target a series of commands? E.g., haar generates an image with sub-images but I would like to modify only one or two of them.
  • Some of the filtered images exhibit boundary issues. Certain native and stdlib commands have boundary_conditions built in. How would I go about giving my filters these capabilities?

The a trous wavelet basically uses dyadic decomposition scales, which is fine for most images, as long as you ask for a sufficient number of scales. The interest of the gaussian decomposition is you can more finely tune the beginning and ending scales. For instance, with only 3 scales, you can get a base scale that has very low frequencies and a higher scale that has very high frequencies (fine details) :

$ gmic sp parrots split_details 3,20,0.2

This is not possible with the a trous wavelets implementation (try split_details 3 to compare).
So if you are interested in manipulating only a single scale, the gaussian decomposition may be more adapted.

I don’t think so, as generating a pyramid of image scales requires only basic operations (blurring and subtraction, basically). Reducing the number of scales is probably the easiest way to speed up the operations.

I’ll suggest using the command crop (shortcut z) to retrieve the region to work on as a single image, then the command image (shortcut j) to put the result back in the haar decomposition.

Not sure I understand this question :slight_smile:
This would mean probably adding a boundary_conditions parameter to your command, and use it when invoking the builtin command used in your custom command ?

Certain commands (I don’t know how many) yield (slightly) different results each run. This might be related to the discrepancy in processing time. I normally ignore this but, now that I am trying to write my own commands and testing their effectiveness against existing ones, I would like to find out why and explore more consistent testing methods. E.g.,

gmic sp tiger +denoise_haar[0] , psnr psnr={i(0,1)} q
# yields 3 different results
# psnr=32.87042236328125
# psnr=32.921798706054688
# psnr=32.954963684082031

I am also curious about SNR, PSNR and SSIM. So far I have been calculating and comparing images using SNR and PSNR. Anyway, I just need more info on these metrics and how they ought to be applied in interpreting the effectiveness of a smoothing, sharpening, etc, filter.

Commands that make use of a random number generator have this behavior. This is the case for denoise_haar that shifts the image with random offsets before applying the haar transform for each iteration of the loop.

Definitely not an expert in understanding those metrics, so I’m skipping my turn here :slight_smile:

It would be great if G’MIC had a command for SSIM. I need to do more reading :slight_smile:. See:


Looks like ImageMagick recently got SSIM added to its -metric command. Unsure how close it is to Fred’s implementation or the original uwaterloo one. Maybe G’MIC could do better :slight_smile:: more robust implementations are available now. E.g., this tool contains a bunch of them.

@David_Tschumperle Any thoughts on my previous post? In the meantime, I thought of 2 more questions:

  1. What are the eqns that you use for MSE and PSNR? I can’t seem to calculate them on my own. Sorry, I am still not used to reading the code and doing the math. I am also figuring out a way to calculate a comparative SNR. Basically, I want to push out some interesting stats. (Entropy seems like another interesting one but I will stick with requesting SSIM {and maybe others listed on my links that you find interesting} for now :sunny:).

  2. I don’t know how to select certain pixels in an image. Taking an image with a Bayer pattern and separating it into its components would be a nice exercise. The components would either have holes or be compacted to half the dimensions of the original image. That knowledge would translate to fun ways to play with images.

Please forgive me if I spit out what you already know:

  1. The G’MIC versions match the usual definitions; MSE takes the difference between each matching pixel from two images and squares it, adds all those up and divides by number of pixels per image. Very like variance but comparing per pixel instead of to an average. You can see the code in the G’MIC std lib if you search for ‘psnr’. Normally the snr is used to see which of some images is closest to the original, so perhaps not very useful for comparing the numbers - just that one is higher/closer than the other.
    Edit: example of MSE could be gmic image.png ++ 2 - sqr e {is/whds}

  2. It really depends what you intend to do with them afterwards - there’s no real concept in G’MIC of ‘selection’ such as you might see in gimp, so you’re left with copying values around buffers. Spreading pixels in a pattern can be done lots of different ways in G’MIC (as usual!). I usually use warp, fill or permute. It’s even possible to use a 3D point cloud. The math parser is probably the most intuitive of those, but not necessarily fastest.

1 Like

Thanks for the insight.

That is what I have been doing. I have read somewhere that there could be a comparative SNR and the main reason that I ask is that PSNR as its name indicates is interested in an image’s peak rather than the average, so it got me thinking.

I would want to learn how to do it via fill, since I have yet to learn about warp, permute or pointcloud :blush:. This is what I would like to do, which of course could be expanded to other pixel herding possibilities:

image
Creative Commons Licence by @afre, using Excel :shushing_face:


New question: So far I have been using i in my fill calculations but how about I? I must admit that I struggle with anything that has to do with vectors at this point. However, in terms of thresholds, it might be useful. E.g., when I do a threshold on a color image using i, it does it for each channel, yielding a colored result. But say I wanted to replace all colors with negative values with nans or some other value, how would I do that?

Fill iterates over each pixel of the referenced image(s), so you can imagine the commands of your fill string running each time it moves to a new pixel. A simple example of removing pixels with fill:

gmic sp dog f "x%2?i:0"

which is saying “if the current x value modulo 2 is 1 then keep the current pixel, else replace with 0”.
More complex things can be done by referencing pixels of another image, for example by using i(#index,…).

Vectors at most basic can be good for operating on several values at once, e.g. if you multiply I you’re multiplying R, G, B each by the same. The usual vector functions are available too. I wouldn’t suggest replacing values with nan because it’s usually more of a problem to handle than just setting some extreme value (very large/very negative/very small), or keeping a mask in another image. Comparison operators on the image as a whole might help (eq/gt/lt…) as well as min/max. Maybe not the best way to learn, but most of what I use was gleaned by reading other filters which did a similar thing.

It’s quite ‘low level’ for this sort of thing and not a lot different from working with buffers in something like C. Loops, logic and maths!

I remember reading about modulo… I figured it out. E.g., to get R, I do

gmic 2,2,1,1 f 1 f (x+1)%2?i:0 f (y+1)%2?i:0

I couldn’t combine it into one fill command and had to handle x and y separately. Is there a way to do that? I suppose the alternative way to tackle this problem is to unfold the grid into a 1d vector, black out every other or few pixels, and then fold it back.


I am confused about ||.

image

E.g., I assumed the following would black out the Gs but it doesn’t:

gmic 2,2,1,1 f 1 f "(x+1)%2?i:0||(y+1)%2?i:0"

image

which is

image


As for I, multiplying I is the same as multiplying i. What I would like to know is how to replace RGB(A) pixels that have negative values. If I used i, it would only address the offending channels but not the pixel as a whole.