[GIS] How to estimate distance from points to line feature with v.distance

grass

I'm trying to use the v.distance module in GRASS 6.4.3 to estimate the distance from a series of point features (animal tracking locations) to the nearest point along a line feature (a coastline). I load in the two shapefiles (both UTM/WGS84),create a new column in the attribute table of the point features to take the estimated distances and start v.distance. I set the 'from' vector map to be the point layer and the 'to' vector map to be the coastline, tick "distance to the nearest point on the 'to' feature along that linear feature", and set the column to upload to be my new distance column.

The script runs ok but the estimated values in the attribute table are complete nonsense. At least half the values are identical, even though the animal was progressively moving away from the coastline and therefore distance should have been increasing. When the animal was 2000km from the coast the distance estimate is 166m, but right nect to the coast it is often >3000m! Funnily though, if I select the option to create an output vector map containing lines connecting the nearest elements, the lines perfectly connect each point to the nearest part of the coastline and are all different lengths.

What is going wrong here?! The frustrating thing is that I ran it once when testing and it worked perfectly, but all subsequent attempts have resulted in the same problems. I can't remember if I used any different settings on the first go. I'm sure it's something simple but it is driving me nuts trying to work out what it is! I know it may be hard to comment without seeing the data, but i'm really hoping someone has run into this before or has some idea where my mistake might be.

All help greatly appreciated! Many thanks,

Best Answer

If you're selecting "distance to the nearest point on the 'to' feature along that feature" then you're asking for the distance along that line feature from it's start. If you want the distance from the points then just choose "minimum distance to nearest feature". The command at the bottom of the v.distance window should include: "upload=dist" and NOT "upload=to_along".

Related Solutions

[GIS] How to calculate distances in a point sequence

I return to this issue because it is very similar to How do I find vector line bearing in QGIS or GRASS? and it can be solved with Python in the same way:

1) Haversine distance

One can find lots of scripts by searching Haversine distance with Python on the Internet and I choose one of them in Haversine Formula in Python (Bearing and Distance between two GPS points)

def haversine(lon1, lat1, lon2, lat2):
    """
    Calculate the great circle distance between two points 
    on the earth (specified in decimal degrees)
    """
    # convert decimal degrees to radians 
    lon1, lat1, lon2, lat2 = map(math.radians, [lon1, lat1, lon2, lat2])
    # haversine formula 
    dlon = lon2 - lon1 
    dlat = lat2 - lat1 
    a = math.sin(dlat/2)**2 + math.cos(lat1) * math.cos(lat2) * math.sin(dlon/2)**2
    c = 2 * math.asin(math.sqrt(a)) 
    km = 6367 * c
    return km

We have a series of lines (points) in the file that must be treated in pairs (point1 - point2) to calculate the distance. For this we will use a simple iterator from Most pythonic way to get the previous element

def offset(iterable):
    prev = None
    for elem in iterable:
        yield prev, elem
        prev = elem

Now it is possible to read the file (example of Kerrie) in pairs of lines/points

import csv
with open('testhavers.csv', 'rb') as f:
   reader = csv.DictReader(f)
   for  pair in offset(reader):
       print pair

 (None, {'LAT': '10.08527', 'LON': '124.59833', 'ID': '1', 'TIME': '21:24:37'})
 ({'LAT': '10.08527', 'LON': '124.59833', 'ID': '1', 'TIME': '21:24:37'},
 {'LAT':    '10.08523', 'LON': '124.59830', 'ID': '2', 'TIME': '21:25:07'})
 ({'LAT': '10.08523', 'LON': '124.59830', 'ID': '2', 'TIME': '21:25:07'}, 
 {'LAT': '10.08526', 'LON': '124.59832', 'ID': '3', 'TIME': '21:25:37'})
 ({'LAT': '10.08526', 'LON': '124.59832', 'ID': '3', 'TIME': '21:25:37'}, 
 {'LAT':    '10.08526', 'LON': '124.59831', 'ID': '4', 'TIME': '21:26:07'})

Then create a shapefile containing the original fields of the csv file and a new field for the distance with the Python modules Shapely and Fiona of Sean Gillies:

import fiona
from shapely.geometry import Point, mapping
# creation of the schema of the shapefile (geometry and fields)
schema = { 'geometry': 'Point', 'properties':{'ID': 'int', 'LAT':'float', 'LON':'float', 'TIME':'str','distance' : 'float'}}
# creation of the shapefile:
with fiona.collection("result.shp", "w", "ESRI Shapefile", schema) as output:
    # reading the csv file
    with open('testhavers.csv', 'rb') as f:
       reader = csv.DictReader(f)
       # we need here to eliminate the first pair of point with None
       for i, pair in enumerate(offset(reader)):
            if i == 0: (pair with None)
                # writing of the point geometry and the attributes
                point = Point(float(pair[1]['LON']), float(pair[1]['LAT']))
                dist = 0 # None
                output.write({'properties': {'ID':int(pair[1]['ID']),'LAT':float(pair[1]['LAT']),'LON':float(pair[1]['LON']), 'TIME':pair[1]['TIME'],'distance': dist},'geometry': mapping(point)})
             else:
                # writing of the point geometry and the attributes
                point = Point(float(pair[1]['LON']), float(pair[1]['LAT']))
                # Haversine distance between pairs of points
                dist = haversine(float(pair[0]['LON']), float(pair[0]['LAT']), float(pair[1]['LON']),float(pair[1]['LAT']))
                output.write({'properties': {'ID':int(pair[1]['ID']),'LAT':float(pair[1]['LAT']),'LON':float(pair[1]['LON']), 'TIME':pair[1]['TIME'],'distance': dist},'geometry': mapping(point)})

and the result: enter image description here

It is also possible to do it with PyQGIS but it is more complex than Fiona which uses simple dictionaries to create shapefiles.

You can use another function to calculate the Haversine distance (Why is law of cosines more preferable than haversine when calculating distance between two latitude-longitude points?) without any problem, only the distance calculation changes, not the process of creating the shapefile.

[GIS] Finding line with minimum distance to point and update attribute of points

A simple solution for your question is using QGIS tool Distance to Nearest Hub. However, this tool does not bring into account your particular way of calculating distance using a 90 degree angle, as you specified. It simply measures the distance from each input point to the middle of each line feature.

Here is an example of how the tool can be used in a similar context:

Input:

Parameters (Hub Layer name attribute is the attribute from the line layer you wish to keep, in your case 'id'):

Output (same points as input, with added attribute):

As with any other processing algorithms, make sure both your layers are using the same coordinate systems (in your case you had a different CRS for each input layer). If the distance field (Hub Dist) is of importance to you too (not just which feature is closest), then make sure you are also using projected coordinate systems (not 4326).
Hope this helps!

Best Answer

Related Solutions

[GIS] How to calculate distances in a point sequence

[GIS] Finding line with minimum distance to point and update attribute of points

Related Question