topic Re: Find 3d distance between 2 sets of points in 3D Mapping Questions

Find 3d distance between 2 sets of points

SB5 — Mon, 05 Dec 2016 10:21:16 GMT

Hi,

I'd like to find the 3d distance between all points in two point feature classes. One is 12 000 building points, the other 22 000 point locations of a construction machine (along a line), constructing many km of tunnel. My aim is to be able to select the building points within x meters when the construction machine is at a specific location. How can I do this without running Near 3d several thousand times?

I am using arcgis 10.3 with basic license. I am not very familiar with scripting.

Re: Find 3d distance between 2 sets of points

XanderBakker — Mon, 05 Dec 2016 12:37:50 GMT

The Near 3D—Help | ArcGIS for Desktop tool will get you the distance to the nearest object in the other layer, not all objects. The Generate Near Table—Help | ArcGIS for Desktop will calculate the distance between all objects in two featureclasses, but not in 3D.

I suppose you would have to script this, but this means evaluating 12000 x 22000 = 264M combinations. That will be a pretty long process. Also to work with the results will be slow, since your table will have 264M records.

To calculate the 3D distance between is not very complex. Just take the square root from (dx² + dy² + dz²).

Re: Find 3d distance between 2 sets of points

SB5 — Mon, 05 Dec 2016 13:17:00 GMT

Thanks Xander. I guess the easiest solution is to run Near 3D after specifying one location of the construction machine in the first feature class. My aim was to have this available for my colleagues to query from a table, but there would be too much data.

Re: Find 3d distance between 2 sets of points

XanderBakker — Mon, 05 Dec 2016 18:10:20 GMT

Just wondering... is there a lot of difference in Z between the building points and the construction machine points?

The fact that you will limit the search distance can limit a lot the number of results. Also if you apply a "Relate" and define your attribute indexes correctly, this might yield be a workable solution.

Is it possible to share your data, to see what's can be achieved?

Kind regards, Xander

Re: Find 3d distance between 2 sets of points

DarrenWiens2 — Mon, 05 Dec 2016 19:11:00 GMT

You could use Generate Near Table to calculate the 2d distances ('c' below), then Spatial Join to transfer start/end height attribute (difference is 'd' below), then calculate the 3d distance to a new field ('e' below).

Re: Find 3d distance between 2 sets of points

SB5 — Wed, 07 Dec 2016 09:15:59 GMT

The Z difference is between 5 and 100 meters.

Limiting the search radius would at the most cut the number of building points in two, leaving me with 6000 points.The building points are all within 400 meters of the line that the construction machine points run along, and I would need to select building points within 150-250 meters from each construction machine point.

I unfortunately cannot share the data.

Thanks again Xander.

Re: Find 3d distance between 2 sets of points

SB5 — Wed, 07 Dec 2016 09:16:50 GMT

Thanks, this seems interesting. I unfortunately only have the basic license and need advanced to use Generate Near Table.

Re: Find 3d distance between 2 sets of points

XanderBakker — Wed, 07 Dec 2016 12:20:08 GMT

The search radius that I was referring to is for searching the buildings for each point on the tunnel line. As you mention that you only want to consider buildings up to 250 m from a tunnel point, this will mean that you will not have to generate the 264M combinations.

I will come back to this and see if I can think of an easy way to generate a result using a Basic license.

Re: Find 3d distance between 2 sets of points

SB5 — Wed, 07 Dec 2016 13:17:06 GMT

Very much appreciated.

I am in the process of trying out the following solution, but not quite sure about it:

I ran 3D buffer tool for each construction machine point and got 22 000 multipatch 200 meter buffer spheres. I then ran Inside 3D to get a table with all the building points within each buffer sphere/multipatch (which had to be joined to the two point tables to preserve all my information in one table). In the resulting table (only 1 million rows) I have xyz for each construction machine point (=the center of the 3D buffer multipatch) and xyz for each building point, can I calculate 3D distance using the coordinates?

Re: Find 3d distance between 2 sets of points

SB5 — Wed, 07 Dec 2016 15:15:57 GMT

Hi again,

As mentioned in a comment to Xander above, I used 3D Buffer and Inside 3D plus some table joining to end up with a table with xyz coordinates for both building points and construction machine points. Can you tell me how to calculate 3D distance in a new field, either using excel or in the attribute table?

This is an example (with very many zeros for some reason)

POINT_X_building	POINT_Y_building	POINT_Z_building	POINT_X_constructionmachineno2850	POINT_Y_constructionmachineno2850	POINT_Z_constructionmachineno2850
599017.980859000000000	6640584.964650000400000	123.000000000000000	599075.213999999690000	6640639.921700000800000	-6.040000000000000
598973.976436999980000	6640633.025369999900000	121.950400000000000	599075.213999999690000	6640639.921700000800000	-6.040000000000000
598973.976436999980000	6640633.025369999900000	121.950400000000000	599075.213999999690000	6640639.921700000800000	-6.040000000000000
599111.963903000000000	6640589.985279999700000	126.000000000000000	599075.213999999690000	6640639.921700000800000	-6.040000000000000
599090.217955000000000	6640602.671959999900000	126.000000000000000	599075.213999999690000	6640639.921700000800000	-6.040000000000000

Thanks so much in advance.

Re: Find 3d distance between 2 sets of points

XanderBakker — Wed, 07 Dec 2016 18:38:07 GMT

Add a new field (double) and calculate the field using this formula:

Sqr (([POINT_X_building]-[POINT_X_constructionmachineno2850]^2 +  [POINT_Y_building]-[POINT_Y_constructionmachineno2850]^2 + [POINT_Z_building]-[POINT_Z_constructionmachineno2850]^2 ))‍

Re: Find 3d distance between 2 sets of points

XanderBakker — Sun, 12 Dec 2021 10:49:53 GMT

As an alternative I generated a script that will calculate the 3D distance between two sets of 3D point featureclasses.

The result is a table that contains the OID of the buildings, the OID of the tunnel points and the 3D distance:

By defining two relates in ArcMap you will be able to do what you are looking for.

Relate 1: based on the tunnel points, relating to the new table.

Relate 2: based on the table relating to the building points:

Make sure to define some attribute indexes in the table to increase the performance:

How does it work. In ArcMap you select a tunnel point and open the attribute table of the tunnel point featureclass:

Click on the relate button and select the relate:

This will open the table with the related combinations selected:

Then select the relate to propagate the selection to the building points:

... which will open the selected building points:

The map will show the selection:

The script I used for this is listed below. In my case I only had 6.5M combinations and it took 2 minutes to process the combinations on an average computer. This will probably take all night processing the 264M combinations.

Script:

import math
import arcpy
from time import localtime, strftime

def main():
    import os
    ShowProgress("Start process...")

    fc_buildings = r'D:\Xander\GeoNet\Generate3DNearTable\gdb\data.gdb\building_points_3D'
    fc_tunnel_points = r'D:\Xander\GeoNet\Generate3DNearTable\gdb\data.gdb\tunnel_points_3D'

    tbl_out = r'D:\Xander\GeoNet\Generate3DNearTable\gdb\data.gdb\Near3D_table_v01'
    fld_building_id = 'OID_Building'
    fld_tunnel_point_id = 'OID_TunnelPoint'
    fld_dist3D = 'Distance3D'

    # max distance
    max_dist = 250

    # sr = arcpy.Describe(fc_buildings).spatialReference

    # create lists of buildings and tunnel points
    ShowProgress("Create lists of buildings and tunnel points...")
    flds = ('OID@', 'SHAPE@')
    lst_building = [[r[0], r[1]] for r in arcpy.da.SearchCursor(fc_buildings, flds)]
    lst_tunnelpoints = [[r[0], r[1]] for r in arcpy.da.SearchCursor(fc_tunnel_points, flds)]

    # create output table for relate
    ShowProgress("Create output table...")
    tbl_ws, tbl_name = os.path.split(tbl_out)
    arcpy.CreateTable_management(tbl_ws, tbl_name)
    ShowProgress("Add fields to output table...")
    AddField(tbl_out, fld_building_id, "LONG", 8)
    AddField(tbl_out, fld_tunnel_point_id, "LONG", 8)
    AddField(tbl_out, fld_dist3D, "DOUBLE", 8)


    cnt_total = len(lst_building) * len(lst_tunnelpoints)
    flds = (fld_building_id, fld_tunnel_point_id, fld_dist3D)

    ShowProgress("Start creating combinations...")
    with arcpy.da.InsertCursor(tbl_out, flds) as curs:

        cnt = 0
        cnt_found = 0
        for building_data in lst_building:
            oid_building = building_data[0]
            pntg_building = building_data[1]
            for tunnel_data in lst_tunnelpoints:
                oid_tunnel_point = tunnel_data[0]
                pntg_tunnel_point = tunnel_data[1]
                cnt += 1
                if cnt % 10000 == 0:
                    ShowProgressPercentage("Processing combination: {0},  found: {1}".format(cnt, cnt_found), cnt, cnt_total)

                dist3D = get3Ddistance(pntg_building, pntg_tunnel_point)
                if dist3D <= max_dist:
                    cnt_found += 1
                    row = (oid_building, oid_tunnel_point, dist3D, )
                    curs.insertRow(row)

    ShowProgressPercentage("Ready processing {0} combinations,  found: {1}".format(cnt, cnt_found), cnt, cnt_total)


def AddField(fc, fld_name, fld_type, fld_length):
    if len(arcpy.ListFields(fc, fld_name)) == 0:
        arcpy.AddField_management(fc, fld_name, fld_type, None, None, fld_length)


def get3Ddistance(pntg1, pntg2):
    pnt1 = pntg1.firstPoint
    pnt2 = pntg2.firstPoint
    distance = math.sqrt((pnt1.X-pnt2.X)**2 + (pnt1.Y-pnt2.Y)**2 + (pnt1.Z-pnt2.Z)**2)
    return distance


def ShowProgressPercentage(txt, cnt, cnt_tot):
    current_time = strftime("%H:%M:%S", time.localtime())
    if cnt_tot is None:
        print "{0}   - {1}".format(current_time, txt)
    else:
        percentage = float(cnt) * 100.0 / float(cnt_tot)
        print "{0}   - {1} ({2}%)".format(current_time, txt, round(percentage, 2))

def ShowProgress(txt):
    # current_time = strftime("%Y-%m-%d %H:%M:%S", gmtime())
    current_time = strftime("%H:%M:%S", time.localtime())
    print "{0} - {1}".format(current_time, txt)

if __name__ == '__main__':
    main()
‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

Re: Find 3d distance between 2 sets of points

SB5 — Thu, 08 Dec 2016 11:58:18 GMT

Thank you so much for spending time on the script, Xander. I ran/started the script with the two point sets and it got to 20% in two hours before I realised that I wanted to split up the point sets, so I just recently started over. I don't think the process will finish before my long weekend starts, but it looks like it's going very well! I appreciate this a lot.

Re: Find 3d distance between 2 sets of points

SB5 — Thu, 08 Dec 2016 13:38:06 GMT

I have a small problem: Since the objectid in my tunnel point fc is called "objectid" and the new table's objectid is also called "objectid", the tunnel point objectid is not included and cannot be related to the original tunnel point file afterwards. Is there a way of changing the objectid, either choosing another field or renaming it?

Re: Find 3d distance between 2 sets of points

SB5 — Thu, 08 Dec 2016 13:43:25 GMT

...or can another identifying field be added from the tunnel point attribute table? Each tunnel point has a unique number that is used by the contractors to track the machines, so it needs to be there (either joined or added now)

Re: Find 3d distance between 2 sets of points

SB5 — Thu, 08 Dec 2016 15:49:23 GMT

(Sorry for spamming)

Nevermind, I just changed what field to include from the tunnel point attribute table, thought it had to be the objectid field. The script works great, I'm very grateful!

Re: Find 3d distance between 2 sets of points

SB5 — Thu, 08 Dec 2016 15:58:51 GMT

I'm back. Do the fields have to be indexing fields / objectid fields? I used the unique tunnel identifying field now, as mentioned above, and the numbers were changed during the process. So the solution was not as simple as I had hoped. I'm back to wondering if it is possible to change the objectid field name or being able to add a few more fields from the tunnel point attribute table.

Thanks in advance.

Re: Find 3d distance between 2 sets of points

XanderBakker — Thu, 08 Dec 2016 22:00:37 GMT

The idea is was to create a table that contains both ObjectID values. Those are read from the input featureclasses. The input ObjectID values won't change and allow you to make the connection between the two featureclasses. The output table contains this relation as would a relationshipclass (which would be available if you would have access to a Standard or Advanced license). It os possible to use a different field as identifier. I used ObjectID since this is standard available in a featureclass. However, this could require to adapt the type of field that is created in the output table.

The alternative could be to join the identifier you want to the output table, using the ObjectID field of the featureclass and the OID_* field in the output table as join fields.

I would suggest to run the process as is and follow the instructions described earlier to see if it works. The advise to create attribute indexes on the fields you use in the relates is to speed up the visualization in ArcMap.

Re: Find 3d distance between 2 sets of points

DanPatterson_Retired — Thu, 08 Dec 2016 23:15:15 GMT

Xander... do you still have that point file? Not the ridiculously large one, but the option within a selection constraint. I want to do some tests and I don't have a real 3D point file that large with real world coordinates to compare to. I am testing e_dist (from previous blogs) in 3D. A post here or email to me would be great

Re: Find 3d distance between 2 sets of points

XanderBakker — Fri, 09 Dec 2016 03:44:53 GMT

Hi Dan Patterson‌, I have attached a ZIPped GDB that contains the two 3D points files that I used and some intermediate data I used to create the data. The two point featureclasses I used are marked in red:

I simply drew a line, extracted a point every 10 meters, buffered the line, and generated random points in the buffer polygon, created Z fields with random values and converted the 2D points files to 3D using the attribute.

Very interested to see your results.

Kind regards, Xander