Equalize image stack - CLIJ2_multiplyImageStackWithScalars

Hi @haesleinhuepf or one other who is familiar with CLIJ2.

I develop a IJm macro and want to equalize/normalize gray values in a stack. I want to have the average gray value in an ROI to be the same in every slice of the stack. I already have the average gray values in the ROI and now hold the derived correction factor in an array.

Next step is to use

Ext.CLIJ2_multiplyImageStackWithScalars(Image source, Image destination, Array scalars);

My question is: how to push a linear (1D) array to the GPU. In the CLIJ2 docs I see:

Ext.CLIJ2_pushArray(Image destination, Array input, Number width, Number height, Number depth);

Should I just use my 1D array name (correction_vector) for input and then set width, height, depth respectively to 1, 1, corection_vector.lenght ? (the last number is of course also the nr of slices in the stack).

Or is there another CLIJ routine to get an array into the GPU?

My confusion comes from the fact that multiplyImageStackWithScalars requires an Array, whereas pushArray delivers an image.

Thanks
Fred

`

1 Like

Hey @Fred_ckx,

Yes. Also in such cases it’s always worth a try to run a simple macro to see if the method does what one thinks it should do:

run("CLIJ2 Macro Extensions", "cl_device");

array = newArray(1, 2, 3, 4);

Ext.CLIJ2_pushArray(image, array, 4, 1, 1);

Ext.CLIJ2_print(image);

image

I also often use pushResultsTableColumn to get a column from the results table to GPU memory.

Very interesting, indeed! There are not many functions taking an array as input. If it hurts, I’m happy to make an additional function that takes an image instead of an array.

Last but not least, there is an alternative: I recently figured out that you can run mathematical operations on images with different sizes. CLIJ applies clamp-to-edges and thus, you can multiply an image stack with a vector in z - thats an image with width=1, height=1, depth>1.

If you try out this macro:


run("CLIJ2 Macro Extensions", "cl_device");

array = newArray(
	1, 1, 
	2, 2, 
	3, 3, 
	4, 4);

width = 2;
height = 1;
depth = 4;
Ext.CLIJ2_pushArray(image, array, width, height, depth);

print("Image before:");
Ext.CLIJ2_print(image);

factors = newArray(8,4,2,1);
Ext.CLIJ2_pushArray(factor_image, factors, 1, 1, 4);

// multiply images of different size
Ext.CLIJ2_multiplyImages(image, factor_image, result);

print("Image after:");
Ext.CLIJ2_print(result);

It prints that:

Let me know if that helps!

Cheers,
Robert

1 Like

Thanks @haesleinhuepf !

I will see what is in this for me :slight_smile:

Of course I tried to do it in the way that I suggested in the following way:

	//get measured mean gray values)
	ref_Value = Table.getColumn("Mean");
	target_Value = 85;   //must later be set in a suitable way

	pre_Factor = newArray(ref_Value.length);
	for (ii = 0; ii < ref_Value.length; ii++){
		pre_Factor[ii] = target_Value/ref_Value[ii];
		};
	//Now correct the stack 
    selectImage("frames");
    inputStack = getTitle();
    Ext.CLIJ2_push(inputStack);

    
	//  Ext.CLIJ2_pushArray(Image destination, Array input, Number width, Number height, Number depth);
	Ext.CLIJ2_pushArray(factorCLIJ2, pre_Factor, 1, 1, pre_Factor.length);

	//	Ext.CLIJ2_multiplyImageStackWithScalars(Image source, Image destination, Array scalars);
	Ext.CLIJ2_multiplyImageStackWithScalars(inputStack, equalizedStack, factorCLIJ2);
	Ext.CLIJ2_pull(equalizedStack);

This code throws an error “array expected” for factorCLIJ. Apparently it is not considered as an array.

If I don’t overlook other things… the main difference between my code and the first code in your reply is that you use multiplyImages and I use multiplyImageStackWithScalars.

Fred

1 Like

Yes: The parameters for pushArray are “Image destination, Array input”. Thus, pushArray converts an array into an image. In the next line you pass an image where you should pass an array.

Try this:

Ext.CLIJ2_multiplyImageStackWithScalars(inputStack, equalizedStack, pre_Factor);

Does that work?

1 Like

Maybe… :slight_smile: It probably does. The error is gone, but I get a different one in a different place which causes that I cant judge the result.

So, in you last suggestion pre_Factor is a variable in the CPU and not in the GPU… if that is correct it would probably help if you add such a remark to the docs :-).

Orrrrrr: [EDIT] maybe everything that isn’t an image, such as an array, is in the CPU???

I’ll notify you about the result when I can collect it.

thanks
Fred

1 Like

Great idea! I’ll add that. In the meantime: There are no arrays in clij. Arrays are a concept of ImageJ macro. When pushing them to the GPU, they become 1-D images…

1 Like

Haha that matches with my remark that I added in an EDIT of my previous message…

1 Like

Hmmm, I see strange behavior. Whatever I put in the pre_Factor array, even zero values in all elements, the equalizedStack is identical to inputStack…

EDIT: no actually your suggestion works:

I am probably hustling with variable names / references.

Fred

1 Like

@haesleinhuepf, sorry for the confusion:

It does not work. The following code should give an all black stack as I have filled the pre_Factor array with all zeroes:

    run("CLIJ2 Macro Extensions", "cl_device=");
	Ext.CLIJ2_clear();
    ref_Value = newArray(35);
	pre_Factor = newArray(ref_Value.length);
	
	for (ii = 0; ii < ref_Value.length; ii++){
		pre_Factor[ii] = 0; 					//set factor to zero for clearest result				//target_Value/ref_Value[ii];
		};
	
	//Now correct the stack 
    selectImage("frames");
    inputStack = getTitle();
    Ext.CLIJ2_push(inputStack);
	Ext.CLIJ2_multiplyImageStackWithScalars(inputStack, equalizedStack, pre_Factor);
	Ext.CLIJ2_pull(equalizedStack);
	rename("equalized frames");

The resulting output stack is identical to the input stack.

Very puzzling…

Fred

1 Like

Hey @Fred_ckx,

it’s concerning. But we’ll figure that out.

Could you please run this little macro:

run("Close All");
newImage("test", "8-bit ramp", 100, 100, 10);

run("CLIJ2 Macro Extensions", "cl_device=");
Ext.CLIJ2_clear();

ref_Value = newArray(35);
pre_Factor = newArray(ref_Value.length);

for (ii = 0; ii < ref_Value.length; ii++){
	pre_Factor[ii] = 0; 					//set factor to zero for clearest result				//target_Value/ref_Value[ii];
	};

//Now correct the stack 
inputStack = getTitle();
Ext.CLIJ2_push(inputStack);
Ext.CLIJ2_multiplyImageStackWithScalars(inputStack, equalizedStack, pre_Factor);
Ext.CLIJ2_pull(equalizedStack);
rename("equalized frames");

On my machine, two windows open:

image image

If I change the image to 32 bit and the factor to -0.5:

run("Close All");
newImage("test", "32-bit ramp", 100, 100, 10);

run("CLIJ2 Macro Extensions", "cl_device=");
Ext.CLIJ2_clear();

ref_Value = newArray(35);
pre_Factor = newArray(ref_Value.length);

for (ii = 0; ii < ref_Value.length; ii++){
	pre_Factor[ii] = -0.5; 					//set factor to zero for clearest result				//target_Value/ref_Value[ii];
	};

//Now correct the stack 
inputStack = getTitle();
Ext.CLIJ2_push(inputStack);
Ext.CLIJ2_multiplyImageStackWithScalars(inputStack, equalizedStack, pre_Factor);
Ext.CLIJ2_pull(equalizedStack);
rename("equalized frames");

image image

Do you see the same?

1 Like

Yes, in both cases. But that doesn’t make it less puzzling for me :slight_smile:

But wait!

What I want to do is to change the average gray scale in each slice of a stack. What you show is to change the intensity profile within an image.

I want to change each pixel in a slice by the same factor. The factor is different for each slice.

In the end the dynamic range within each slice is the same as before the treatment. Only the average gray value changes in each slice. The correction factors per slice have been measured before on a small Reference ROI .

Fred

@haesleinhuepf,

OK, now I get it…

The correction factors where all the same in you example. When I do this:

run("Close All");
newImage("test", "8-bit ramp", 100, 100, 35);

run("CLIJ2 Macro Extensions", "cl_device=");
Ext.CLIJ2_clear();

ref_Value = newArray(35);
pre_Factor = newArray(ref_Value.length);

for (ii = 0; ii < ref_Value.length; ii++){
	pre_Factor[ii] = ii/35; 	}				//set factor to zero for clearest result				//target_Value/ref_Value[ii];
	

//Now correct the stack 
inputStack = getTitle();
Ext.CLIJ2_push(inputStack);
Ext.CLIJ2_multiplyImageStackWithScalars(inputStack, equalizedStack, pre_Factor);
Ext.CLIJ2_pull(equalizedStack);
rename("equalized frames");

I can actually see that each factor now infact is the multplying factor over z !

This does not directly tell me why my original code is not working.

I will gave it some further inspection…

Maybe a clue is in the fact that my images are 32 bits (RGB). Your stack is 8 bit or 32 bit gray LUT. Here comes my lack of knowledge on how images are coded?

Fred

Can you make a minimal macro, that is executable and that reproduces the issue? Please share it together with the data set. It’s otherwise hard to guess what’s happening. I can’t see your screen :wink:

I’ll do that. That will be tomorrow.

I also have to prepare my dataset as I cannot share the original.

What is the best way to send you a (part of) the image stack?

1 Like

In order to reproduce this issue, a 10x10 pixel image with 10 frames should be enough, right? You can share that as tiff here in the forum :slightly_smiling_face:

OK. I have less time tomorrow. So it will be probably tomorrow night before I post something.

Thanks

Hi @haesleinhuepf,

I think I found a clue. When I run you small test macro on my stack, the console window shows the following:

My stack is 32 (RGB).

When I create several test ramps similar to yours, only 8,16,32 bit run well and RGB fails.

Fred

1 Like

Hey @Fred_ckx,

clij doesn’t support multi-channel data for technical reasons and also because most users don’t want to apply operations to channels. For example: Blurring between channels introduces cross-talk. Background subtraction between channels is hard to justify if you think about wavelengths and how images are acquired. Segmentation-based analysis typically means you segment one channel to analyze another. So typically channels are treated independently. There are not many scenarios, where you really want to process multi-channel data as such.

Can you process your data channel by channel?

Edit: Furthermore, do you need the information in the channels? Or would your intended analysis also work on a greyscale image?

Cheers,
Robert

Hi @haesleinhuepf,

I think that I need the color information for separation. So I need to think of a way to keep that.

Does this mean that I have to create a stack that contains the additional dimension c with 3 subpages R, G and B?

Fred

There are at least two ways for processing multi-channel data: As we discussed here you can push a channel stack. Or alternatively, as showin in this example by pushing a stack per timepoint. In the second case, channels are then represented as slices. You can then copy a slice to process one channel independent from another.

As pointed out above that you want to multiply slices (time points?) with values, I would suggest a strategy like this:

Stack.setChannel(1);
Ext.CLIJ2_pushCurrentZStack(image);
result1 = process(image);

Stack.setChannel(2);
Ext.CLIJ2_pushCurrentZStack(image);
result2 = process(image);

Stack.setChannel(3);
Ext.CLIJ2_pushCurrentZStack(image);
result3 = process(image);

Ext.CLIJ2_pull(result1);
Ext.CLIJ2_pull(result2);
Ext.CLIJ2_pull(result3);

run("Merge Channels...", "c1=" + result1 + " c2=" + result2 + " c3=" + result3 + " create");

I hope that helps

Cheers,
Robert