| Re: vector layer comparison in IDL [message #42751] |
Fri, 25 February 2005 15:16  |
Mark Hadfield
Messages: 783
Registered: May 1995
|
Senior Member |
|
|
yp wrote:
> Thanks for your suggestion. I did think of the same; to bring both
> layers to same grid. But I am stuck here... How should I go about
> interpolating the (x,y) pairs (No attributes are taken into account at
> the moment)? Here I have to interpolate the co-ordinate rather than
> data. If I go like interpolating all y's f(x) and all x's f(y) I'll end
> up with wrong interpolated data when there are more than one parallel
> lines.
So it's a problem of interpolating from one spatially irregular grid to
another, or perhaps of interpolating from both to a common grid. I have
posted on this in the past (to everybody's confusion) & will have a
think about it over the rest of the weekend. But first...
What is the shape of the grid elements in your vector layers? Triangles,
polygons? Or do your vector layers even have "grid elements". Are they
just lists of (x,y,data) with no implied spatial organization of x and
y? In this case it is normal to use Delaunay triangulation to create a
triagular grid linking the points.
Generally, the key to interpolating from one grid to another (say grid A
to grid B) is to locate the points in grid B relative to those in grid
A. In doing so you want to take advantage of regularity in the geometry
of the grids.
I *think* I understand your last point about parallel lines. Generally,
unless the grids are very simple, you'll have to consider x and y at the
same time.
IDL does have several functions, the most general one being GRIDDATA,
for interpolating from scattered data values. You may want to look at
the documentation.
--
Mark Hadfield "Ka puwaha te tai nei, Hoea tatou"
m.hadfield@niwa.co.nz
National Institute for Water and Atmospheric Research (NIWA)
|
|
|
|
| Re: vector layer comparison in IDL [message #42764 is a reply to message #42751] |
Fri, 25 February 2005 06:35  |
yp
Messages: 42
Registered: February 2005
|
Member |
|
|
Mark Hadfield wrote:
> yp wrote:
>> Mark Hadfield wrote:
>>
>>> yp wrote:
>>>
>>>> Is there a sensible way to compare two exported vector layers
>>>> and compute the residuals using IDL? I am sorry for asking this
>>>> simple question...
>>>
>>> Not at all. I, for one, have no idea what your simple question
>>> means. Is this some sort of GIS thing?
>>
>> Yes Mark, it is related to GIS. Then, you might advise to use some
>> GIS package. However, I have have 100s of *files* - these files are
>> basically ARC vector layers exported as ASCII files (lon, lat, ...,
>> attributes). For instance, I have 50 sites where the surveys have
>> been done in two different times (hence 100 ASCII tables); estimate
>> the absolute changes (residuals or differences). I believe that
there
>> should be some way of doing it in IDL...
>
> OK, so what do *you* think is the way to go about this and where do
you
> foresee difficulties? Can you get the data into IDL? Do the lon and
lat
> data differ between 2 surveys at the same site? If not, simple
> subtraction of the attribute data should be OK; if so, I guess you'll
> need to interpolate to a common grid first. Or is there more to it?
>
>
> --
> Mark Hadfield "Ka puwaha te tai nei, Hoea tatou"
> m.hadfield@niwa.co.nz
> National Institute for Water and Atmospheric Research (NIWA)
Hi Mark,
Thanks for your suggestion. I did think of the same; to bring both
layers to same grid. But I am stuck here... How should I go about
interpolating the (x,y) pairs (No attributes are taken into account at
the moment)? Here I have to interpolate the co-ordinate rather than
data. If I go like interpolating all y's f(x) and all x's f(y) I'll end
up with wrong interpolated data when there are more than one parallel
lines.
Any clue?
|
|
|
|
| Re: vector layer comparison in IDL [message #42809 is a reply to message #42764] |
Wed, 23 February 2005 12:06 |
Mark Hadfield
Messages: 783
Registered: May 1995
|
Senior Member |
|
|
yp wrote:
> Mark Hadfield wrote:
>
>> yp wrote:
>>
>>> Is there a sensible way to compare two exported vector layers
>>> and compute the residuals using IDL? I am sorry for asking this
>>> simple question...
>>
>> Not at all. I, for one, have no idea what your simple question
>>means. Is this some sort of GIS thing?
>
> Yes Mark, it is related to GIS. Then, you might advise to use some
> GIS package. However, I have have 100s of *files* - these files are
> basically ARC vector layers exported as ASCII files (lon, lat, ...,
> attributes). For instance, I have 50 sites where the surveys have
> been done in two different times (hence 100 ASCII tables); estimate
> the absolute changes (residuals or differences). I believe that there
> should be some way of doing it in IDL...
OK, so what do *you* think is the way to go about this and where do you
foresee difficulties? Can you get the data into IDL? Do the lon and lat
data differ between 2 surveys at the same site? If not, simple
subtraction of the attribute data should be OK; if so, I guess you'll
need to interpolate to a common grid first. Or is there more to it?
--
Mark Hadfield "Ka puwaha te tai nei, Hoea tatou"
m.hadfield@niwa.co.nz
National Institute for Water and Atmospheric Research (NIWA)
|
|
|
|
| Re: vector layer comparison in IDL [message #42824 is a reply to message #42809] |
Wed, 23 February 2005 02:23 |
yp
Messages: 42
Registered: February 2005
|
Member |
|
|
Mark Hadfield wrote:
> yp wrote:
>> Is there a sensible way to compare two exported vector layers and
>> compute the residuals using IDL? I am sorry for asking this simple
>> question...
>
> Not at all. I, for one, have no idea what your simple question means.
Is
> this some sort of GIS thing?
>
Yes Mark, it is related to GIS. Then, you might advise to use some GIS
package. However, I have have 100s of *files* - these files are
basically ARC vector layers exported as ASCII files (lon, lat, ...,
attributes). For instance, I have 50 sites where the surveys have been
done in two different times (hence 100 ASCII tables); estimate the
absolute changes (residuals or differences). I believe that there
should be some way of doing it in IDL...
|
|
|
|
| Re: vector layer comparison in IDL [message #42856 is a reply to message #42824] |
Thu, 03 March 2005 14:58 |
David Fanning
Messages: 11724
Registered: August 2001
|
Senior Member |
|
|
Mark Hadfield writes:
> Jim Pendleton of RSI has pointed out to me that there is an IDL function
> called PNT_LINE that does this.
Oh, no! *How* many hours did you say!? Some days
you just shouldn't get out of bed. :-)
Cheers,
David
--
David Fanning, Ph.D.
Fanning Software Consulting, Inc.
Coyote's Guide to IDL Programming: http://www.dfanning.com/
|
|
|
|
| Re: vector layer comparison in IDL [message #42857 is a reply to message #42824] |
Thu, 03 March 2005 14:42 |
Mark Hadfield
Messages: 783
Registered: May 1995
|
Senior Member |
|
|
Mark Hadfield wrote:
> function magic_distance_function, x, y, x0, y0, x1, y1
>
> compile_opt DEFINT32
> compile_opt STRICTARR
> compile_opt STRICTARRSUBS
> compile_opt LOGICAL_PREDICATE
>
> d0 = sqrt((x-x0)^2+(y-y0)^2)
> d1 = sqrt((x-x1)^2+(y-y1)^2)
> dp = 2.0*poly_area([x,x0,x1],[y,y0,y1]) / $
> sqrt((x1-x0)^2+(y1-y0)^2)
>
> return, dp > (d0 < d1)
>
> end
Jim Pendleton of RSI has pointed out to me that there is an IDL function
called PNT_LINE that does this.
--
Mark Hadfield "Ka puwaha te tai nei, Hoea tatou"
m.hadfield@niwa.co.nz
National Institute for Water and Atmospheric Research (NIWA)
|
|
|
|
| Re: vector layer comparison in IDL [message #42861 is a reply to message #42824] |
Thu, 03 March 2005 13:26 |
David Fanning
Messages: 11724
Registered: August 2001
|
Senior Member |
|
|
Mark Hadfield writes:
> Applying this to the surface of the earth is left as an exercise, as is
> testing.
Wow! I'm impressed. How much are you charging per hour
to do this. I have a couple of projects. :-)
Cheers,
David
--
David Fanning, Ph.D.
Fanning Software Consulting, Inc.
Coyote's Guide to IDL Programming: http://www.dfanning.com/
|
|
|
|
| Re: vector layer comparison in IDL [message #42862 is a reply to message #42824] |
Thu, 03 March 2005 13:08 |
Mark Hadfield
Messages: 783
Registered: May 1995
|
Senior Member |
|
|
yp wrote:
> Mark Hadfield wrote:
>> Are you interested, perhaps, in the maximum distance between
>> the two "coastlines"?
>
> You guessed it right. I got to calculate the absolute error
> (distance) assuming that one is *true* and the other is
> *estimated* from other sources.
Ok, so we have 2 polylines, the "true" one defined by vectors
xt and yt, both dimensioned [n], and the "estimated" one
defined by vectors xe and ye, both dimensioned [m]. We want to
calculate the maximum distance of the "estimated polyline" from
the "true" one. Here's a naive approach:
dmax = 0.
for i=0,m-1 do begin
for j=1,n-1 do begin
d = magic_distance_function(xe[i], ye[i], $
xt[j-1], yt[j-1], xt[j], yt[j])
if j eq 1 then dmin = d else dim = dmin < d
endfor
dmax = dmax > dmin
endfor
(I think I've got the logic right there. We're looping through the
"estimated" vertices, for each one calculating the distance to the
closest line segment on the "true" polyline.)
Now you might object that this solution is a little vague in the line
where magic_distance_function is invoked. That's a good point. So here's
an attempt at this function
function magic_distance_function, x, y, x0, y0, x1, y1
compile_opt DEFINT32
compile_opt STRICTARR
compile_opt STRICTARRSUBS
compile_opt LOGICAL_PREDICATE
d0 = sqrt((x-x0)^2+(y-y0)^2)
d1 = sqrt((x-x1)^2+(y-y1)^2)
dp = 2.0*poly_area([x,x0,x1],[y,y0,y1]) / $
sqrt((x1-x0)^2+(y1-y0)^2)
return, dp > (d0 < d1)
end
This uses three relevant distances between point [x,y] and line segment
[x0,y0] -> [x1,y1]. D0 and d1 are distances to the end points; dp is the
perpendicular distance between [x,y] and the line through [x0,y0] &
[x1,y1], extending the line as far as necessary. I'm too lazy to look up
the expression for perpendicular distance, so I've taken advantage of
the relationship between this distance and the area of the triangle
formed by the 3 points, calculated by IDL function POLY_AREA.
Applying this to the surface of the earth is left as an exercise, as is
testing.
You might want to think about the case where the "estimated" polyline
follows the "true" polyline closely, but has extra vertices at one end
or other. The above algorithm sees this as an error.
--
Mark Hadfield "Ka puwaha te tai nei, Hoea tatou"
m.hadfield@niwa.co.nz
National Institute for Water and Atmospheric Research (NIWA)
|
|
|
|
| Re: vector layer comparison in IDL [message #42881 is a reply to message #42824] |
Thu, 03 March 2005 03:00 |
yp
Messages: 42
Registered: February 2005
|
Member |
|
|
Mark Hadfield wrote:
> yp wrote:
>> Well, here I have two sets exported coordinated (X=lon, Y=Lat) of
the
>> same area at two different times. The task is to see the spatial
>> difference between the two layers (In this example I dont have any
>> closed polygon, but that should be kept in mind for future)
>
> Well, they do look very similar, but not identical...
>
> Have you given any thought to what the term "spatial difference"
could
> or should mean? What are you trying to achieve with your comparison?
> Once you work that out, there's no doubt you can use IDL to calculate
it
> for you. Are you interested, perhaps, in the maximum distance between
> the two "coastlines"?
>
> --
> Mark Hadfield "Ka puwaha te tai nei, Hoea tatou"
> m.hadfield@niwa.co.nz
> National Institute for Water and Atmospheric Research (NIWA)
Thanks Mark,
You guessed it right. I got to calculate the absolute error (distance)
assuming that one is *true* and the other is *estimated* from other
sources.
|
|
|
|
| Re: vector layer comparison in IDL [message #42923 is a reply to message #42824] |
Tue, 01 March 2005 12:55 |
Mark Hadfield
Messages: 783
Registered: May 1995
|
Senior Member |
|
|
yp wrote:
> Well, here I have two sets exported coordinated (X=lon, Y=Lat) of the
> same area at two different times. The task is to see the spatial
> difference between the two layers (In this example I dont have any
> closed polygon, but that should be kept in mind for future)
Well, they do look very similar, but not identical...
Have you given any thought to what the term "spatial difference" could
or should mean? What are you trying to achieve with your comparison?
Once you work that out, there's no doubt you can use IDL to calculate it
for you. Are you interested, perhaps, in the maximum distance between
the two "coastlines"?
--
Mark Hadfield "Ka puwaha te tai nei, Hoea tatou"
m.hadfield@niwa.co.nz
National Institute for Water and Atmospheric Research (NIWA)
|
|
|
|
| Re: vector layer comparison in IDL [message #42924 is a reply to message #42824] |
Tue, 01 March 2005 02:25 |
yp
Messages: 42
Registered: February 2005
|
Member |
|
|
Mark Hadfield wrote:
> yp wrote:
>> Hi Mark,
>> Thanks for all your suggestions. I guess, I cant use IDL GRIDDATA
>> function for this job since the 3rd dimension (z) is empty; i.e.,
here
>> we have only (x,y). So I am still scratching my forehead to about
>> interpolating the space. Let me read through the gridding section
(IDL
>> online refs are not excellent... I have few IDL handbooks written
by
>> David Fanning. Hope I'll find a solution..)
>
> What is it you actually want to do?
>
> --
> Mark Hadfield "Ka puwaha te tai nei, Hoea tatou"
> m.hadfield@niwa.co.nz
> National Institute for Water and Atmospheric Research (NIWA)
Well, here I have two sets exported coordinated (X=lon, Y=Lat) of the
same area at two different times. The task is to see the spatial
difference between the two layers (In this example I dont have any
closed polygon, but that should be kept in mind for future)
X1 Y1 X2 Y2
-4.226679575 50.31360946 -4.226832766 50.31355545
-4.226203672 50.31351466 -4.226254354 50.31269132
-4.226013311 50.31275628 -4.225242132 50.3125473
-4.225061507 50.31275628 -4.224374514 50.31168316
-4.224585605 50.31209269 -4.223073086 50.31096305
-4.224109702 50.3114291 -4.221627055 50.31168316
-4.222777176 50.31114471 -4.22090404 50.31240327
-4.22153983 50.3115239 -4.220614834 50.31312338
-4.220588025 50.31247188 -4.221048643 50.31398752
-4.220492845 50.31304067 -4.221627055 50.31427556
-4.221349469 50.31398865 -4.221627055 50.31427556
-4.22153983 50.31455744 -4.220614834 50.3145636
-4.22153983 50.31455744 -4.218734994 50.31470763
-4.220492845 50.31465224 -4.216565947 50.31557176
-4.219350679 50.31465224 -4.21483071 50.31643589
-4.218874777 50.31484183 -4.21295087 50.317156
-4.218494055 50.31493663 -4.210781824 50.31802013
-4.21744707 50.31531582 -4.209169115 50.31834138
-4.216114544 50.31597941 -4.209169115 50.31834138
-4.2153531 50.3160742 -4.208612778 50.3184522
-4.215638642 50.31673779 -4.206588335 50.3184522
-4.214496476 50.31654819 -4.203985479 50.31816415
-4.212973589 50.31730658 -4.200370402 50.31802013
-4.210594077 50.31806496 -4.197333737 50.31787611
-4.208880829 50.31844415 -4.194875485 50.31816415
-4.2070724 50.31844415 -4.192851042 50.31830818
-4.206120596 50.31853895 -4.191260408 50.31730002
-4.206120596 50.31806496 -4.189375556 50.31755032
-4.204788069 50.31806496 -4.189375556 50.31755032
-4.20278928 50.31844415 -4.189091361 50.31758807
-4.200314588 50.31797016 -4.188512949 50.31802013
-4.198315798 50.31797016 -4.189380568 50.31888426
-4.19660255 50.31815976 -4.18995898 50.32090057
-4.194603761 50.31834936 -4.189235964 50.32248481
-4.192985693 50.31825456 -4.190034148 50.32345645
-4.19241461 50.31787537 -4.190034148 50.32345645
-4.191653166 50.31778057 -4.190537392 50.32406905
-4.190225459 50.31778057 -4.192272629 50.32522123
-4.188702572 50.31787537 -4.195887706 50.32709352
-4.188702572 50.31787537 -4.198201356 50.32824569
-4.189273655 50.31825456 -4.200515005 50.32925385
-4.189178474 50.31863375 -4.20167183 50.33055004
-4.189273655 50.31929734
-4.189749557 50.31958173
-4.188988114 50.31977133
-4.189749557 50.32052971
-4.189559196 50.3209089
-4.189844738 50.32128809
-4.188892933 50.32176208
-4.189844738 50.32195168
-4.188892933 50.32252047
-4.189939918 50.32356324
-4.190199747 50.32378321
-4.190199747 50.32378321
-4.191843527 50.32517481
-4.193842317 50.32612279
-4.195460385 50.32621759
-4.197078452 50.32735516
-4.198315798 50.32830314
-4.199933866 50.32877713
-4.200885671 50.32981991
-4.201456753 50.33086269
|
|
|
|
| Re: vector layer comparison in IDL [message #42926 is a reply to message #42824] |
Mon, 28 February 2005 11:48 |
Mark Hadfield
Messages: 783
Registered: May 1995
|
Senior Member |
|
|
yp wrote:
> Hi Mark,
> Thanks for all your suggestions. I guess, I cant use IDL GRIDDATA
> function for this job since the 3rd dimension (z) is empty; i.e., here
> we have only (x,y). So I am still scratching my forehead to about
> interpolating the space. Let me read through the gridding section (IDL
> online refs are not excellent... I have few IDL handbooks written by
> David Fanning. Hope I'll find a solution..)
What is it you actually want to do?
--
Mark Hadfield "Ka puwaha te tai nei, Hoea tatou"
m.hadfield@niwa.co.nz
National Institute for Water and Atmospheric Research (NIWA)
|
|
|
|
| Re: vector layer comparison in IDL [message #42927 is a reply to message #42824] |
Mon, 28 February 2005 11:27 |
yp
Messages: 42
Registered: February 2005
|
Member |
|
|
Mark Hadfield wrote:
> My last post was rather hurried, which (oddly enough) caused me to
use a
> lot of words to say very little. What I should really have said
was...
>
> It sounds like you have a straightforward problem of 2D interpolation
> from data that are not on a rectangular grid. Use GRIDDATA. Read the
> section in the manual on gridding and check out the gridding demo.
>
>
> --
> Mark Hadfield "Ka puwaha te tai nei, Hoea tatou"
> m.hadfield@niwa.co.nz
> National Institute for Water and Atmospheric Research (NIWA)
Hi Mark,
Thanks for all your suggestions. I guess, I cant use IDL GRIDDATA
function for this job since the 3rd dimension (z) is empty; i.e., here
we have only (x,y). So I am still scratching my forehead to about
interpolating the space. Let me read through the gridding section (IDL
online refs are not excellent... I have few IDL handbooks written by
David Fanning. Hope I'll find a solution..)
|
|
|
|
| Re: vector layer comparison in IDL [message #42942 is a reply to message #42751] |
Sun, 27 February 2005 13:18 |
Mark Hadfield
Messages: 783
Registered: May 1995
|
Senior Member |
|
|
My last post was rather hurried, which (oddly enough) caused me to use a
lot of words to say very little. What I should really have said was...
It sounds like you have a straightforward problem of 2D interpolation
from data that are not on a rectangular grid. Use GRIDDATA. Read the
section in the manual on gridding and check out the gridding demo.
--
Mark Hadfield "Ka puwaha te tai nei, Hoea tatou"
m.hadfield@niwa.co.nz
National Institute for Water and Atmospheric Research (NIWA)
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
