Distance to value in raster instead of to shapefile

[ViriatoII]

Hey,

This package looks awesome! Is it possible to calculate distance to a pixel value in other raster?

For example, I have this binary raster of values 0/255. I'd like to calculate distance to 255 values:

It looks like it should be easy, since you transform any shapefile into a binary raster already, I just don't know what "rv_array and affine" should be then.

Cheers, Ricardo

[sgoodm]

Hey @ViriatoII,

Yes, you absolutely can use an existing raster. The rv_array and affine would be defined along the lines of:

import rasterio
with rasterio.open("your_raster.tif") as src:
    affine = src.transform()
    rv_array = src.read(1)

I'll add an example to the repo based on this approach for other folks to reference in the future. Thanks!

[ViriatoII]

That's great, thank you! 2 things:

• src.transform() gives a "not callable" error with those parenthesis, it works without them.

• Further down the line I get a math error:

  
  ~/.local/lib/python3.8/site-packages/distancerasters/main.py in _calculate_distance(self)
   97 
   98                     else:
  ---> 99                         km_min_dist = calc_haversine_distance(
  100                             convert_index_to_coords(cur_index, self.affine),
  101                             convert_index_to_coords(min_index, self.affine),

~/.local/lib/python3.8/site-packages/distancerasters/utils.py in calc_haversine_distance(p1, p2) 243 * math.sin(delta_lon / 2) * 2 244 ) --> 245 c = 2 math.atan2(math.sqrt(a), math.sqrt(1 - a)) 246 # km 247 d = radius * c

ValueError: math domain error

Just to make sure, does this affine and rv_array look good?

Affine: | 1.00, 0.00, 0.00| | 0.00, 1.00, 0.00| | 0.00, 0.00, 1.00| rv_array: [[0 0 0 ... 0 0 0] [0 0 0 ... 0 0 0] [0 0 0 ... 0 0 0] # There are values of 255 here but they are a minority ... [0 0 0 ... 0 0 0] [0 0 0 ... 0 0 0] [0 0 0 ... 0 0 0]]

The problem must stem from my data. It's not geographical but cellular data, the coordinates vary between 0 and 15,000. I guess I should replace **calc_haversine_distance** with my own function.

Cheers.
Ricardo

[sgoodm]

Ah, yes the current implementation will always apply haversine distance calculations if an affine is provided and assumes it is just geographic, unprojected data. I would definitely like to improve this approach to work with other types of data and be a bit more intuitive.

A quick fix that may work for your case is to not provide an affine or output_path (e.g., just dr.DistanceRaster(rv_array)) in order to avoid the block of code using the haversine calc. You would just need to save the resulting distance array yourself using the export_raster function since the DistanceRaster class requires an affine in order to export as a raster.

[ViriatoII]

Thank you for the help. The function runs now but unfortunately after 2 hours it still hasn't finished. Perhaps my raster is too big (195 MB)? Using a similar function in imageJ (based in Java, I think) takes less than a minute to run.

[sgoodm]

What are your raster dimensions? Finer-resolution data over large areas are likely to take some time, and unfortunately will definitely be slower based on Python than Java.

I would like to produce some detailed performance/timing metrics in the future, but for now the best I can say is that I have produced a global ~1km resolution (36000x18000) distance raster which also took around a few hours.

[ViriatoII]

Hey @sgoodm

I was trying with a file that had 13080x15696 pixels (my analysis is at pixel resolution). Now I'm using the function on subsets of the data and it works beautifully.

I work in spatial transcriptomics in intra- and inter-cellular space. So here you have RNA and their distance to the nucleus of the cell.

For anyone searching for a simple solution like me without Geographical projections, I found 3 other functions that work with rasters as numpy arrays, explained in detail here:

• cv2.distanceTransform()
• PyDIP.EuclideanDistanceTransform()
• scipy.ndimage.distance_transform_edt()

With that I close the issue, thank you!

Ricardo