Subscribe to this thread
Home - General / All posts - M9 Modifying pixels in a raster using a vector features as analysis masks
danb


1,734 post(s)
#15-May-19 05:16

I can’t remember if I have asked this before but have wanted to be able to do it for a while in M9. I have a vector feature which I want to transpose onto an image to modify or extract the overlain pixels.

Thinking of ways to do this, I know there was forum mention of a method similar to M8’s 'no interpolation' method when pasting a vector feature to a grid, but I haven’t seen anything similar in the current range of interpolators. We also now have features for performing zonal type functions collecting data under a vector feature from an overlain raster, but as far as I am aware, not the other way round.

With this in mind, is it currently possible to extract the pixels under a vector feature or interact with them such as setting their value in SQL, or is this a feature request?

Making (irregular) selections in images and modifying the pixel values in some way is a very powerful feature of M8 which I would love to see in M9.


Landsystems Ltd ... Know your land | www.landsystems.co.nz

Dimitri


5,892 post(s)
#15-May-19 05:48

If you don't mean the raster to vector interaction as shown in Example: Transfer DEM Terrain Heights to Areas in a Drawing then that would be a feature request.

tjhb

9,159 post(s)
#15-May-19 08:09

Dan means the other way around: use geometry to look up underlying pixels, then do something with the pixels.

The "something" could be (a) do raster maths, (b) create a new raster from those pixels, (c) modify an existing raster (the source raster or another) only within the scope of those pixels. (Other things?)

I tend to think selecting pixels is a red herring. It is often the most convenient method of doing those things in 8, but is not necessary in 9.

Anyway, all of these things are possible already in 9, just not as easy as they could be.

An exact example might be good, e.g. a Manifold 8 workflow with simple data. Then we can agree on how it can currently be done in 9, and how it might be made easier/more approachable.

danb


1,734 post(s)
#15-May-19 21:16

Thanks Tim. I will put together some simple data for a use case.


Landsystems Ltd ... Know your land | www.landsystems.co.nz

danb


1,734 post(s)
#16-May-19 01:00

Find attached a M8 test project. It contains main three components [VECTOR], [RASTER] and

[1 UPDATE SELECTED PIXELS].

[VECTOR] and [RASTER] are subsets taken from larger extent components in a M9 project.

[1 UPDATE SELECTED PIXELS] is the way I would approach this in M8 i.e. select all objects in [VECTOR], transfer this selection to [RASTER] and then run this query.

The M9 project creates the surface as a subset of a lidar grid including all pixels with heights <=1.8m.and all other pixels set to invisible. The vector feature is the result of a ‘Trace Areas’ transform operating on the same grid with some polygons such as those along the river removed from the resulting drawing.

For this example, all I wish to do is use Manifold 9 to update the pixels under the vector polygons with a value of 99. As Tim suggests however, I am also interested in more general cases of interactions between geom and pixels as outlined.

Attachments:
M8 FORUM.map


Landsystems Ltd ... Know your land | www.landsystems.co.nz

tjhb

9,159 post(s)
#17-May-19 01:13

Cool data!

tjhb

9,159 post(s)
#17-May-19 02:20

Here are the steps I am thinking of, with some effort at optimization.

  • Reproject areas to CS of image
  • Subdivide areas into convex parts, with relevant attribute; result to temporary table T with RTREE on part
  • Make tile boxes for image; result to temporary table U with RTREE on box
  • Filter boxes in U overlapping some part in T
  • Split filtered tiles into pixels, drawing pixel boxes; result to temporary table V with RTREE on box
  • Filter pixel boxes in V overlapping some part in T; result to temporary table W with BTREEDUP on XY tile index
  • Adjust value of joined pixels in W, applying attribute from matched part
  • Left join tiles in U to adjusted pixels in W (join on tile index)
  • Recombine joined tiles with coalesce: adjusted pixel if any, else original
  • Substitute recombined tiles for original tiles as applicable
tgazzard
138 post(s)
#13-Jan-20 00:14

I would like to update the pixel values in a raster using the values in an area vector drawing. Is this now possible using a transform or sql approach in m9?

I have looked through the manual and looked at the latest video on closest rasters which does some interaction with vector points. However, it isn't obvious to me how this would work in m9. This is something that I regularly undertook in m8. Any help appreciated.

tjhb

9,159 post(s)
#13-Jan-20 00:22

It’s not built-in yet Tim.

I think this is Dan’s top outstanding feature request, mine too.

Besides that, it’s definitely time I had a proper go at writing as good a solution as I can.

Overdue. No excuses...

tgazzard
138 post(s)
#13-Jan-20 00:27

Thanks for confirmation that there isn't a feature (yet) that does this.

If you do write a solution, I would love to see this.

dyalsjas120 post(s)
#13-Jan-20 00:32

Concur, this is a needed capability when working with GIS and images.

Would also like to transform image pixels to points (and back) e.g. elevation image to x,y,z points and the reverse, x,y,z points to elevation image.

Kriging approaches the x,y,z points to elevation image, but direct translation (without interpolation) is sometimes necessary.

Dimitri


5,892 post(s)
#13-Jan-20 04:57

You can do the (and back) part now, transferring pixel values from an image to points in a drawing. See this example.

tjhb

9,159 post(s)
#13-Jan-20 05:30

Yes. That is the reverse, the "and back" part.

The missing part is to transfer data from geometry (of any kind) to coincident pixels.

dyalsjas120 post(s)
#13-Jan-20 11:14

I reviewed this example, but I'm unable to find an easy way to to create a drawing of points coincident with each pixel in an image. This effort is non-trivial with high resolution DEMs,

Dimitri


5,892 post(s)
#13-Jan-20 12:35

Ah, I see. I misunderstood your interest. I thought that given a raster image and a layer of points above that, you wanted to transfer the value of each pixel under a given point into a field in that point, for example, the way Height values from the raster in the Transform: Closest Raster topic were transferred into the eight points in the Points drawing.

Creating a dense grid of points to overlay a raster is not that difficult. It's tedious, but can be done remarkably quickly given a series of copy/paste/shift and repeat operations, doubling the size of the point swaths that are copied/pasted each time. There have been threads on that.

But, could it be there are easier ways to do what you want than to transform a raster into points?

adamw


8,971 post(s)
#14-Jan-20 16:18

Just wondering, do you need to do this for a couple of points or for massive amounts of points (are we talking about altering some points on a raster or basically re-creating it)?

Both things can be done using a script, but if it is the latter, then maybe there's a different approach you can take and it would be interesting to know why do you want to go from raster to point-for-each-pixel and then back to raster - the first guess is that you want to do something in the middle and cannot find a suitable combination of raster functions to do that in raster (and if so, maybe what we need is a couple more raster functions).

In general, we do have 90% of the machinery to go from vector to raster and will likely do the remaining 10% fairly soon. We have that on the radar for sure.

dyalsjas120 post(s)
#15-Jan-20 01:49

I am trying to recreate a model that compares surface variability factors across ranges of thirty and three hundred meter diameters. This model can be GPU accelerated when using 10m to 30m elevation models, but I haven’t been able to identify a software that will allow me to develop the same model to work on high resolution elevation models like LiDAR data sets.

CUDA raster acceleration seems to drop off after an analysis mask exceeds a modest radius. This model would have 15 and 150 pixel radius masks respectively when using 1m DEMs.

I thought I might try a brute force approach by converting an elevation model to x,y,x points and using a combination of buffers and queries. The underlying math capability already seems built into Manifold.

I also thought I might explore if a combination of file export and import methods might allow me to get to a suitable set of points.

The final result of the model is effectively a raster image of classes similar to a Land Use Land Cover image. I’m sure I could export a point file from Manifold into something I could kludge back into an image using GDAL or ESRI or something.

I'll once again ask for an option to select units (degrees or percentage) when calculating a slope image.

adamw


8,971 post(s)
#15-Jan-20 08:21

Thanks, that helps.

So, the back and forth between raster and vector is indeed fairly massive, but that's because traditional approaches in the raster seem too slow (150 pixel windows are quite slow, that's true), you want to experiment using geometry.

Let's see what you can do right now:

You can convert a raster into a drawing - that's a simple query. (See TileToValues functions, say if you want an example query.)

You can then analyze that drawing using transforms in 9. Or using a script.

You can also export the drawing to, say, SHP, or perhaps export the table to CSV with X-Y-Z values, and process that data in other software, then import the result back into 9.

After you have the result, if it is just pixels and the pixels are still in their integer positions with just the values changed, it is fairly simple to put them back into a raster using a script. Again, say if you want example code.

If the result is areas or lines, then perhaps wait until we do vector to raster as a built-in tool.

PS: Request for degrees or percentages in slopes noted - I'll bump the wishlist item.

dyalsjas120 post(s)
#16-Feb-20 00:43

Adam,

It's taken a while for me to get back to this posting.

I've reviewed the various topics here in the user forum that included examples of how to create a table from tile values using the various tiletovalue functions.

I haven't been able to get a working query.

Can you (or any of the other forum members) assist with an example query of how to create a table of x, y, and z values from a single band image. The image I'm currently working with is an image of slope values.

I'm in the middle a trying an x.y.x file export, but it's surprisingly slow.

v/r

jd

adamw


8,971 post(s)
#17-Feb-20 08:52

Sure.

Given an image named 'german_alps', with single-channel tiles of the standard 128x128 size, this will produce the table with x - y - value for each pixel:

--SQL9

VALUE @tilesizex INT32 = 128;

VALUE @tilesizey INT32 = 128;

SELECT tilex * @tilesizex + x AS x, tiley * @tilesizey + y AS y, value FROM (

  SELECT x AS tilex, y AS tiley, SPLIT CALL TileToValues(tile)

  FROM [german_alps]

);

The inner query splits all tiles into pixels and keeps tile X and tile Y for each pixel. The outer query combines tile X / Y and pixel X / Y within a tile into pixel X / Y within an image, and keeps pixel value.

You can convert the outer query into SELECT ... INTO to put this into a table. If you do that, it might be worth it to add MFD_ID / MFD_ID_X to that table to make it editable - or you can create a BTREE index on X / Y (the combinations of which are unique), that will make the table editable as well.

dyalsjas120 post(s)
#17-Feb-20 19:21

I've tried this query, but I'm not sure the image I'm working with has standard tile sizes.

The tile description shows <tile, 256x16, float32>

When running the query above, the value column in the query result shows NULL.

I've tried updating the tile sizes to 256, but still receive NULL within the value column of the query result.

tjhb

9,159 post(s)
#17-Feb-20 21:06

text

text

When running the query above, the value column in the query result shows NULL.

That might be normal...

(1) How far are you scrolling down? If you are scrolling all the way to the bottom, remember that that means you get the "first" 50,000 records, which here means the "first" 50,000 pixels. ("First" here is arbitrary in theory, but they will usually be drawn from the tiles with the lowest indexes, though this could be disrupted.)

(2) Could the displayed pixels all be <null>, and this still be normal? Yes, if you have a large number of empty tiles. That is, tiles that are present, but contain only invisible pixels. That can be the case with images imported/linked from other formats, including Manifold 8. Manifold 9 is normally more efficient, since it normally only writes tiles which have some visible pixel(s).

To test more definitely, just add a filter to Adam's outer query.

SELECT...

FROM

(

SELECT...

)

WHERE [Value] IS NOT NULL;

That will return up to 50,000 records showing visible pixels... or no records at all, if all pixels were invisible.

adamw


8,971 post(s)
#18-Feb-20 07:58

I've tried updating the tile sizes to 256, but still receive NULL within the value column of the query result.

Maybe that's what you did, but just in case - if your tiles are 256x16, set @tilesizex to 256 and @tilesizey to 16, do not set both to 256.

As Tim says, NULL values represent invisible pixels. If you don't want them, you can filter them out with WHERE.

dyalsjas120 post(s)
#22-Mar-20 17:55

I've just revisited this query, (including the filter to exclude NULL values), and the output is what I expected.

My next question is: How do I have the X and Y values reflect the projected coordinates of each pixel in the image instead of the pixel coordinates?

danb


1,734 post(s)
#15-Jan-20 04:36

I think this is Dan’s top outstanding feature request, mine too.

Right on the money Tim. Definitely my number one currently.

I have been away for a while but what a treat to come back to. Two great builds.

I am very much looking forward to really getting to grips with the distance tools and the las library looks like a game changer for me and the way we tend to store las data (primarily to cater to legacy GIS).

The improvements to day to day usability such as 'rename related' are just great.


Landsystems Ltd ... Know your land | www.landsystems.co.nz

Manifold User Community Use Agreement Copyright (C) 2007-2019 Manifold Software Limited. All rights reserved.