Subscribe to this thread
Home - General / All posts - TileFilterVar
dyalsjas
157 post(s)
#17-May-20 01:16

Quick question to confirm my assumption

Does the TileFilterVar SQL function (and Variance template) provide the max minus min of values contained within the radius, filter shape, of the input image tile?

What is the difference between the TileFilterVar and TileFilterVarPop, and Variance / Variance Population templates?

tjhb
10,094 post(s)
#17-May-20 03:40

TileFilterVarPop should calculate the sum of squared differences between each value within the defined radius (weighted by the defined filter shape) and the average of all those values, the sum then being divided by the number of values (or sum of weights) including the centre value itself.

TileFilterVar (without Pop, meaning sample) should calculate almost exactly the same, except at the end, where it will divide by the number of values minus 1 (or if weights are used, then by the sum of weights divided by the number of values minus 1, I think).

The Pop version assumes we know all relevant neighbouring values; the sample version assumes that we don't (e.g. because NULLs are prevalent in the image, and the NULL values represent not true empty areas but unknowns).

That is simplistic. Anyone who knows statistics better (that is ~everyone), please correct it.

dyalsjas
157 post(s)
#17-May-20 14:25

Tim, thanks for confirming my assumption. It's been waaayyy too long since I thought hard about statistics.

I can probably still use variance instead of the cheetah flips needed to subtract the output of a TileFilterMin from the output of a TileFilterMax.

tjhb
10,094 post(s)
#18-May-20 03:38

I didn't think I had confirmed your assumption, but contradicted it. Perhaps we were not clear enough together?

There is no direct relation between TileFilterMin orTileFilterMax and TileFilterVar(Pop).

TileFilterVar(Pop) does not use Min or Max. It uses every value independently, and their combined average.

They can't be substituted or translated one for another.

tjhb
10,094 post(s)
#18-May-20 03:57

If you need the local range, then subtracting TileFilterMin from TileFilterMax is correct. Variance won't give a similar measure.

Take a tile like

1  1  1

1 19  1

1  1  1

The minimum is 1, the maximum 19, so the range is 18. The average is (8 + 19) / 9 = 3.

The population variance is ((1 - 3) ^ 2) * 8 + ((19 - 3) ^ 2) * 1 = 32 + 256 = 288.

The population standard deviation is the square root of 288, slightly under 17.

adamw


10,447 post(s)
#18-May-20 11:29

You forgot to divide the sum of squares by the number of measurements. So, the result is not 288, but 288 / 9 = 32.

--SQL9

VALUE @t TILE = StringJsonTile('[ 1, 1, 1, 1, 19, 1, 1, 1, 1 ]',

  3, 3, 1, false);

 

? TileJson(TileFilterVarPop(@t, 1, TileFilterDefSquare(1, 1)))

-- nvarchar: [

--  null, null, null,

--  null, 32, null,

--  null, null, null

-- ]

:-)

tjhb
10,094 post(s)
#18-May-20 21:03

Thanks! How silly of me.

Better working of the example too.

tjhb
10,094 post(s)
#17-May-20 10:49

So to answer

Does the TileFilterVar SQL function (and Variance template) provide the max minus min of values contained within the radius, filter shape, of the input image tile?

more directly, no.

adamw


10,447 post(s)
#18-May-20 11:13

As Tim says, TileFilterVar / VarPop compute sample variance and population variance.

You seem to need Max - Min, so you need to do two calls: TileFilterMax(...) - TileFilterMin(...). We have been thinking about adding functions to compute that directly - perhaps calling the measure Range.

dyalsjas
157 post(s)
#18-May-20 13:48

How can I tackle nesting these two functions (TileFilterMax and TileFilterMin) into a query?

I have a set of models that sum/subtract the outputs of two or more filters, either from the results of a filter or from the data the filter operates on.

For this one, I need to determine the range of values within a filter radius e.g. subtract TileFilterMin from TileFilterMax.

If nesting two filters into a query is easy, then perhaps a function/template for range isn't necessary, a suitable example for the help docs would suffice.

If there are performance benefits for having a specific function, then, yes, please.

I think "Range" is an easily understood description.

adamw


10,447 post(s)
#18-May-20 15:32

Well, wherever you are calling TileFilterMin or TileFilterMax, you can define a custom function called TileFilterRange and call that instead, with the same arguments:

--SQL9

FUNCTION TileFilterRange(@t TILE, @radius FLOAT64, @filter TILETILE AS

  TileFilterMax(@t, @radius, @filter) - TileFilterMin(@t, @radius, @filter)

END;

We will still likely add built-in variants of various XxxRange functions, but until we do, you can roll your own.

dyalsjas
157 post(s)
#21-May-20 16:41

So, I'm trying to roll a range function, and I'm starting with the Function description above, then trying to merge it into the generated SQL from the Maximum template. I also read the Process Images using Dual 3x3 Filters help topic.

I've attempted to update the SQL from the Maximum template as follows:

My source image is called [n30 UTM Source Image]

My output table is called [n30 Tiles Range]

My ouput image is called [n30 Tiles Range Image]

My desired radius is 1

My desired filter type is square

I suspect I've declared @t, @radius, and @filter incorrectly.

Here's my attempt at the query, it does not parse.

Thoughts? Ideas? Wry expressions at novice coding?

--SQL9

FUNCTION

TileFilterRange(@t TILE, @radius FLOAT64, @filter TILETILE AS TileFilterMax(@t, @radius, @filter) - TileFilterMin(@t, @radius, @filter)

END;

VALUE @t = [n30 UTM Source Image]

VALUE @radius = 1

Value @filter = 1

CREATE TABLE [n30 Tiles Range] (

  [mfd_id] INT64,

  [X] INT32,

  [Y] INT32,

  [Tile] TILE,

  INDEX [mfd_id_x] BTREE ([mfd_id]),

  INDEX [X_Y_Tile_x] RTREE ([X][Y][Tile] TILESIZE (128,128) TILETYPE INT16),

  PROPERTY 'FieldCoordSystem.Tile' '{ "Axes": "XYH", "Base": "World Geodetic 1984 (WGS84)", "CenterLat": 0, "CenterLon": -93, "Eccentricity": 0.08181919084262149, "FalseEasting": 500000, "LocalOffsetX": 500000, "LocalOffsetY": 3318785.352608442, "LocalScaleX": 27, "LocalScaleY": 27, "MajorAxis": 6378137, "Name": "Universal Transverse Mercator Zone 15 (N)", "ScaleX": 0.9996, "ScaleY": 0.9996, "System": "Transverse Mercator", "Unit": "Meter", "UnitScale": 1, "UnitShort": "m" }',

  PROPERTY 'FieldTileSize.Tile' '[ 128, 128 ]',

  PROPERTY 'FieldTileType.Tile' 'int16'

);

CREATE IMAGE [n30 Tiles Range Image] (

  PROPERTY 'Table' '[n30 Tiles Range]',

  PROPERTY 'FieldTile' 'Tile',

  PROPERTY 'FieldX' 'X',

  PROPERTY 'FieldY' 'Y',

  PROPERTY 'Rect' '[ 0, 0, 3574, 4122 ]'

);

PRAGMA ('progress.percentnext' = '100');

INSERT INTO [n30 Tiles Range] (

  [mfd_id][X][Y],

  [Tile]

SELECT

  [mfd_id][X][Y],

  CASTV ((TileRemoveBorder(TileFilterRange(TileCutBorder(@t, VectorMakeX2([X][Y]), 1), 1, TileFilterDefSquare(@radius, @filter)), 1)) AS INT16)

FROM [n30 UTM Source Image]

THREADS SystemCpuCount();

TABLE CALL TileUpdatePyramids([n30 Tiles Range Image]);

adamw


10,447 post(s)
#21-May-20 16:58

Well, the query does not parse because VALUE declarations miss type - should be VALUE @t TABLE = [n30 ...], VALUE @radius FLOAT64 = ..., VALUE @filter FLOAT64 = ... (the numeric types could also be INT32, the conversion is automatic). You should also end each VALUE statement with a semicolon - depending on further text the query engine might find that the statement ended without a semicolon, but it is better to use one, otherwise you might accidentally make an error and the query engine will execute something other from what you intended.

Frankly, since you are using each variable exactly once after it is declared, you could have just used the values you assigned to them instead.

dyalsjas
157 post(s)
#21-May-20 18:01

Wry facial expressions at novice coding.... yep, well deserved

Corrected the declared values and it works like a champ.

tjhb
10,094 post(s)
#21-May-20 19:05

Frankly, since you are using each variable exactly once after it is declared, you could have just used the values you assigned to them instead.

All the same, in my opinion two good reasons to declare them at the top are to make the given parameters more visible (they are not constants, but current assumptions, which can be nicely isolated), and to facilitate code re-use.

I suppose that it only one reason really, but I like doing it that way too.

A slightly different spin: I like using VALUE, and sometimes FUNCTION, even when they are not necessary, to leave breadcrumbs for a future self, to help understand what I did and to make arbitrary changes easier.

dyalsjas
157 post(s)
#22-May-20 04:09

Adam (and Tim),

Thanks for the guidance, your patience for noob errors is remarkable. I've included a working example of the query (with source data from the inport_dem_SDTS project in the Manifold.net examples.

I've added a number of questions about declaring VALUES for Image or Table names in the range query. I didn't see an easy answer in the user manual.

Finally, I asked earlier this week about subtracting the values in the Tiles of one band from the values in the Tiles of another band.

To tie my interest to the Range function, I would like to combine the output of Range (filter radius 1) and Range (filter radius 4) by dividing the Tile values in each output by 2 and summing the results into a single band. See the Example Combine Query for step 1.

Right now my only Manifold option is to export one of the bands as a set of points and performing a Drawing to Image join.

Attachments:
2020_0521_Range_Example.mxb

adamw


10,447 post(s)
#22-May-20 09:14

Regarding questions about values, can one declare the input / output image / table as a value so that whenever you want to switch to use a new image, you only have to change a single declaration - mostly yes:

If what we are talking about is having a query with a couple of VALUE declarations followed by statements using those values, and about manually editing these VALUE declarations before rerunning the query for each new image - then sure, this is fine. You cannot do CREATE TABLE @t ... / ALTER TABLE @t ... / DROP TABLE @t - these statements require the names of components be spelled out explicitly, they cannot take them as parameters, that's why I say 'mostly yes' instead of just 'yes' above. But you can do INSERT INTO @t ..., SELECT ... FROM @t ..., etc.

You can also create a query which will take components as parameters and use EXECUTE to invoke it. The attached MXB shows how - run the Range Filter query, it will create the output components and then run the Range Filter Parameterized query which will fill them.

Range Filter Parameterized:

--SQL9

 

-- accepts parameters:

--   @input TABLE - input image

--   @output TABLE - output image

--   @radius FLOAT64 - filter radius

--   @center FLOAT64 - center value in square filter

 

FUNCTION TileFilterRange(@IN TILE, @RADIUS FLOAT64, @FILTER TILETILE AS

  TileFilterMax(@IN, @RADIUS, @FILTER) - TileFilterMin(@IN, @RADIUS, @FILTER)

END;

 

INSERT INTO @output (

  [mfd_id][X][Y],

  [Tile]

SELECT

  [mfd_id][X][Y],

  CASTV ((TileRemoveBorder(TileFilterRange(

    TileCutBorder(@input, VectorMakeX2([X][Y]), @radius),

    @radius,

    TileFilterDefSquare(@radius, @center)

  ), @radius)) AS FLOAT32)

FROM @input

THREADS SystemCpuCount();

 

TABLE CALL TileUpdatePyramids(@output);

Range Filter:

--SQL9

 

-- create output image

CREATE TABLE [Range Tiles] (

  [mfd_id] INT64,

  [X] INT32,

  [Y] INT32,

  [Tile] TILE,

  INDEX [mfd_id_x] BTREE ([mfd_id]),

  INDEX [X_Y_Tile_x] RTREE ([X][Y][Tile] TILESIZE (128,128)

    TILETYPE FLOAT32),

  PROPERTY 'FieldCoordSystem.Tile'

    ComponentCoordSystem([1475 ELEVATION Raster]),

  PROPERTY 'FieldTileSize.Tile' '[ 128, 128 ]',

  PROPERTY 'FieldTileType.Tile' 'float32'

);

CREATE IMAGE [Range Image] (

  PROPERTY 'Table' '[Range Tiles]',

  PROPERTY 'FieldTile' 'Tile',

  PROPERTY 'FieldX' 'X',

  PROPERTY 'FieldY' 'Y',

  PROPERTY 'Rect' ComponentProperty([1475 ELEVATION Raster]'Rect')

);

 

-- fill output image

EXECUTE WITH (

  @input TABLE = [1475 ELEVATION Raster],

  @output TABLE = [Range Image],

  @radius FLOAT64 = 1,

  @center FLOAT64 = 1

[Range Filter Parameterized];

Is there a way to programmatically copy the image bounds from one image to the other - sure, see Range Filter above. It copies the coordinate system as well. You can copy any property.

Regarding the last question, if I am reading it correctly you want to combine the channels of two images derived from the same original image, correct? If so, you don't have to go through conversions of pixels to points and back - you can just join one of the images to the other, creating a new channel, and then going through the result subtracting one channel from the other. The join can be performed by the Join dialog. Alternatively, since both source images are derived from the same original image, you can just join the tables on X and Y (img1 INNER JOIN img2 ON img1.x=img2.x AND img1.y=img2.y) and subtract the tile for one image from the tile for another. If the images were in different coordinate systems you could have first reprojected one of the images to the other using one of the CoordConvertTileSet functions, and then proceeded with the same join on X and Y.

PS: Sorry if my post earlier felt wry, I didn't intend it that way. These queries require quite a bit of figuring out, no question. I'll help in any way I can.

Attachments:
2020_0521_Range_Example-v2.mxb

dyalsjas
157 post(s)
#22-May-20 19:35

Adam,

Thanks for highlighting some of the key things that cannot be declared in VALUE statements.

I was easily able to use join dialog to add the output channel of one image as a new channel for another image, resulting in a two band image.

I can't seem to find how to add the pixel values from a channel to the pixel values in another channel.

For this process I'm not contemplating using a filter, that has already been done.

Does the TileChannelsConcat function add the value of pixels/tiles/channel to the pixels/tiles/channel of another band?

Is there any plan to add a set of "TileChannel Math" functions?

There are an several geospatial analysis models I use that add or subtract the pixel values of one channel to or from another.

adamw


10,447 post(s)
#23-May-20 08:48

We do have functions to perform channel math.

You can extract a specific channel from a tile using TileChannel. That takes a tile and the number of the channel you want to extract, and returns a tile with that channel.

You can copy values from a specific channel in one tile to a specific channel in a different tile using TileChannelCopy. That takes two tiles, the number of the channel to read in one tile and the number of the channel to write in the other tile, and returns the modified tile.

You can rearrange or drop channels in a tile using TileChannels. That takes a tile and either a single numeric value or a vector of 2, 3 or 4 values, which specify the desired order for the channels. The output is a tile with as many channels as the number of values you passed. For example, this will recompose a BGR tile into RGB: TileChannels(tile, VectorMakeX3(2, 1, 0))

You can merge channels together using TileChannelsConcat. That takes two tiles and merges them into a new tile, first taking all channels from the first tile and then taking all channels from the second tile. The total number of channels in the result is limited by 4.

Single-channel tiles can be added / subtracted / and so on, using operators like + / - and functions like TileExp, TileFloor, etc. If you have a tile with two channels and want to create a sum of those channels, you can do it like this: TileChannel(tile, 0) + TileChannel(tile, 1). The result of this expression is a tile with a single channel, with each value being the sum of channel 0 and channel 1 in the respective pixel.

Hope this helps.

Dimitri


7,413 post(s)
#23-May-20 10:30

See examples in this topic.

dyalsjas
157 post(s)
#23-May-20 15:16

Thanks to you and Adam. I'll dig into this. I guess I've been overthinking the TileChannel functions.

On a side note, earlier (probably a year ago) I asked about raster morphology functions. Now I realize can probably do this with iterations of StringJsonTile functions.

Here's a website that describes how to compose the filters for morphology.

https://homepages.inf.ed.ac.uk/rbf/HIPR2/hipr_top.htm

adamw


10,447 post(s)
#23-May-20 16:41

Just in case, we are planning to provide built-in functions for morphology operations in the future - they aren't yet there, although yes, you can build custom ones (I would use a script: create a function that takes and returns a tile implementing the operation using a script, then call that function from SQL, there are examples of how to do this on the forum).

tjhb
10,094 post(s)
#18-May-20 21:21

If nesting two filters into a query is easy, ...

I think the hard thing here is that it is so easy.

It can take a little while to fully grasp that tile filter functions and tile arithmetic do much the same thing: they operate on whole tiles at once (at least in logic; often they will be broken up into smaller chunks inside the GPU or CPU), returning a tile as their result.

The significant difference is that tile filter functions take neighbours into account (and thus normally return blanks around the edges, unless the tile is expanded first), while tile arithmetic is just a straight-through "overlay" where cells in matching positions are pairs of operands (and nothing special happens at the edges).

The "stack" of tile operations, whether filter or arithmetic or a combination of both, can be as high as you like. And I'm pretty sure that when GPGPU is used, the whole stack will be processed as a single block of code, with only one transport operation to and from CPU. (Unless you introduce a conditional--often that can be replaced by more arithmetic.) So not just easy, but super fast.

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