Borrar filtros
Borrar filtros

Regionprops handling of big images

17 visualizaciones (últimos 30 días)
Marcio Teixeira
Marcio Teixeira el 12 de Oct. de 2023
Comentada: Ashish Uthama el 16 de Oct. de 2023
Hello. I am processing a quite big image (120 GB) and need to find connected regions and calcularte some properties.
The "preprocessing" of the image works pretty fine. Here is the code:
bim = blockedImage("my_big_image.tif");
tmp = apply(bim,@(bs)processImgBlock(bs,filter_value));
And the function is like this
function bw = processImgBlock(bs,m)
bw = bs.Data == m; % masked image
bw = bwareaopen(bw,250);
bw = bwconncomp(bw);
bw = labelmatrix(bw);
end
So far so good. I get my label matrix ready to be processed. It is a 80k x 80k lines, single band image. According to documentation I will get as return a blocked image that I can retrieve using:
img = gather(tmp);
The issue is that if I run regionprops on img, Matlab does not process it due to memory. Note that I have 64GB in my computer.
So my next try that (quite obviously) failed was:
my_props = apply(tmp,@(bs)get_stats(bs,land_code));
% % % getstats is like follows
function s = get_stats(bs,land_code)
s = regionprops('table',bs.Data,'Area','Perimeter','Eccentricity','Orientation','Centroid');
end
because the function should return a blocked image and not a table.
Here comes the question(s):
  1. Is there any function that works similarly to apply that can perform the regionprops operations?
  2. If not, any suggestion on how to handle it? Splitting image in small parts would be the obvious choice but it could split a connected region in 2 or more subimages and I'd like to avoid it (if possible).
Thank you in advance and please reach me if you need more information or if you find my explanation unclear.

Respuesta aceptada

Walter Roberson
Walter Roberson el 12 de Oct. de 2023
Splitting image in small parts would be the obvious choice but it could split a connected region in 2 or more subimages
Any method that involves splitting the image into sub-images and processing the sub-images, runs the risk that you might miss connectivity that extends over sub-image boundaries.
blockedImage and apply and the older blockproc pass in partial images and assume the processing on the partial images can be handled independently. They do not transform the called function to be able to handle large arrays without running out of memory, and they have no idea how to synthesize the results of processing the sub-images together to form a complete result.
The sort of thing you would need would be to use regionprops with a tall array. Unfortunately, regionprops does not support tall arrays.
You are going to have to split the image up into sub-images, process that, get the results out, and then patch up the boundaries -- such as checking to see whether any connected components extend to the sub-image boundaries in adjacent sub-images, and if so then adjusting the information, such as re-labeling.
bw = bwconncomp(bw);
I am morally certain that is not going to return "correct" information when you use blocked images.
because the function should return a blocked image and not a table.
If you use the older blockproc() then it is not required that the output be the same size as the input block -- but if you use overlapping blocks, then make sure you use options to not trim the output. You have to return something numeric, but as long as you can make the output a consistent size then you could potentially pack data into a numeric array. This would not work very well for returning per-region information (unless you set a fixed maximum number to return), but it would work, for example, to return information about the largest region.

Más respuestas (3)

Image Analyst
Image Analyst el 12 de Oct. de 2023
How big is your bs.data binary image? I would start by trying to delete any variables the in workspace that you don't need anymore. Use the clear function. If you still don't have enough memory, try writing out the binary image to a disk file and then, in a separate program read in the binary file and get your properties. That should work since you won't have any other large variables hanging around in memory anymore.
  5 comentarios
Marcio Teixeira
Marcio Teixeira el 16 de Oct. de 2023
Of course! I have images of land use and land cover. The map consists of land uses and each one has a code e.g. forest = 10, pasture = 20 etc. When I apply the maks in my code I am measuring the desired region. I am interested in perimeter, area, the centroids and orientation of these areas.
The solution from @Walter Roberson was the one I was trying to avoid: split the image. If it is the only possible path, I will go for it. Now there is a suggestion from by @Ashish Uthama. I will try it.
About downsampling I am not so sure. My idea is to collect statistics on thse areas, so it would be nice to keep them in the original 30 x 30 m per pixel.
cheers!
Márcio
Ashish Uthama
Ashish Uthama el 16 de Oct. de 2023
Depending on how large these regions are, downsampling might be worth trying. i.e if they occupy say ~>5-10% each of the full image, the resolution lost in area/perimeter might be small when compared to the full image. (its likely also a function of the details in the mask edges, if its mostly uniform/straight then the loss might be even less).
If you only have a few well defined regions, you could consider using blockedImage's zero-copy crop like this example shows: crop and mask geoTiff Files. Adapt the idea here like so: 1)create a low res image of the mask (retain the WorldStart/WorldEnd of the orginal!) Use regionprops on this to get bounding boxes for each object, then iterate on that to crop the finest level one by one to recompute stats on the finest level individually. This will work as long as your objects are well localized, so a crop is actually much smaller than the full image.
(For questions like these, it really helps to share sample images so some 'non standard' ideas can be formed :D)

Iniciar sesión para comentar.


Ashish Uthama
Ashish Uthama el 13 de Oct. de 2023
Editada: Ashish Uthama el 13 de Oct. de 2023
@Marcio Teixeira - what is the maximum size of an individual 'region' that you are trying to analyze? If its small enough w.r.t to the image the approach outlined in this example: Detect Nuclei in Large Whole Slide Images Using Cellpose ought to work. It shows how to use regionprops on a blockedImage using apply.
It shows how to use blockedImage to process 100kx100k+ images while ensuring that objects on the border are correctly captured only once.
In short:
  1. Determine the maxium size of your individual object. This is the crucial part.
  2. Process large image in blocks with additional border size of max size/2 at least (use apply())
  3. In the processing function, reject objects whose centroids fall on the border region and only retain stats for objects whose centroids fall within the core block region
You ought to be able to modify the helpers in the linked example for your use case.
  2 comentarios
Marcio Teixeira
Marcio Teixeira el 16 de Oct. de 2023
Definetely I can determine the maximum size of my object. This step is OK. Unfortunatelly I don't have the medical image toolbox in my Matlab distribution (I use one of my University).
Nevertheless I found a workaround that I detailed in my answer below.
Thank you for the suggestion and the very interesting example you provided.
Cheers,
Márcio
Ashish Uthama
Ashish Uthama el 16 de Oct. de 2023
Glad you have something that works for you!
Point of clarification: You only need Medical if you want to use cellpose for segmentation. For just blockedImage+regionprops, you only need Image processing Toolbox. So you could just replace the cellpose call in the above example with whatever you currently use for segmentation to adapt it. Or you could just study the code and modify it for your use.

Iniciar sesión para comentar.


Marcio Teixeira
Marcio Teixeira el 16 de Oct. de 2023
Guys!
First of all thanks for your messages. I tried all your suggestions and I found one that is working nicely for me. After loading my image, I filtered it for my pixel of interest (let's say 20). With this binary image (bw) I do the following:
  1. Instead of removing the areas with less than 250 pixels I selected the top 100 using a1 = bwareafilt(bw,100). I came to 100 images by trial and error. This was the limit that I could use regionprops without crashing.
  2. I extracted the properties if interest and kept them in a table
  3. I can use the a1 as a filter to clear the top 100 areas from my original image bw. Then I apply step one again and get the next patches (101 to 200).
Of course I need to clear memory every step and read again the new filtered image. It is not the most elegant solution, but it is working for me.
Again thank you all for valuable inputs.

Productos


Versión

R2022b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by