[GIS] Comparing two TINs created using ArcGIS for Desktop

arcgis-desktopdemlidarspatial-analysttin

I have the following two files:

LiDAR las file1
Let's give it a name LiDAR las file2, contains the same file1 data, where number of points randomly removed.

So if file1 contains X number of points, file2 would contain Y number of points, where Y is less than X.

Now, I am using ArcGIS "LAS to Multipoint" to create from both files above their shape counterparts. Next I am creating two empty TIN(triangulated irregular network) files using ArcGIS. Then I edit both TINs, let's say tin1 and tin2, with the appropriate input feature classes, i.e. file1 to tin1 and file2 to tin2.

Is there any tool/package/library or any software that would let me calculate the RMSE(root mean square error) or any difference information that would give me information on the accuracy of the second file, i.e. file2, with respect to file1? I tried to use the TIN difference but it gives system error and cannot obtain the result. I would really appreciate if anyone could show how to obtain accuracy information and interpretation of the results obtained.

Best Answer

Approach

You either need to solve the system problem or perform a substitute or work-around comparison. Because you don't provide details about the error, I conclude you're looking for a work-around.

Solution

Consider converting both TINs to raster format, where a little map algebra does the rest. For example, the RMS measure of difference between the TINs is obtained by subtracting one grid from the other, squaring the result, and taking the square root of its average value. The answer depends on the cellsize (and grid origin) you use, but as the cellsize grows small, the answer will converge to the one you would obtain using the TINs directly.

Tools

There are many tools and libraries to perform such comparisons. They are the bread-and-butter operations of raster GIS packages like Spatial Analyst, IDRISI, and GRASS.

Accuracy

Note that no such procedure, whether applied to the TINs or some other representation, can estimate accuracy (which concerns agreement between the TIN and reality): it can only compare one data set to the other. Both could be hugely inaccurate, yet in close agreement.

Interpretation

There is nothing unusual about interpreting these results: one inspects both maps and the difference map, looks for outlying values in their differences, characterizes the statistical distribution of the differences, assesses the spatial correlation of the differences, identifies locations where clusters of high or low differences occur, etc., just as always.

Related Solutions

[GIS] Comparing two Digital Elevation Models (DEMs) from LAS files

How to compare two Digital Elevation Models (DEMs) in R.

#-------------------------------------------------------------------------
#Creating a reproducible example

library(raster)

  #simulating raster_1

  f = system.file("external/test.grd", package="raster")
  DEM_1 = raster(f)

  #simulating raster_2

  DEM_2 = DEM_1
  # replacing values from raster_1 to create a new raster sample (raster_2)
    DEM_2[(DEM_2>500 & DEM_2<900)] = 550
    DEM_2[(DEM_2>200 & DEM_2<300)] = 500

#-------------------------------------------------------------------------
# Comparison 1 (DEM_3 resulted from subtracting DEM_2 from DEM_1)

  DEM_3 = DEM_1 - DEM_2

    par(mfrow=c(1,3))

    plot(DEM_1, main = "DEM_1")
    plot(DEM_2, main = "DEM_2")
    plot(DEM_3, main = "DEM_3 = DEM_1 - DEM_2")

      dev.off()

enter image description here

#-------------------------------------------------------------------------
#Comparison 2 (histogram)

  hist(DEM_1, prob=T, main="DEM_1", xlab="")
  hist(DEM_2, prob=T, main="DEM_2", xlab="")
  hist(DEM_3, prob=T, main="DEM_3 = DEM_1 - DEM_2", xlab="")

    par(mfrow=c(1,1))

  standard_deviation = sd(c(as.matrix(DEM_3)),na.rm=T)

    dev.off()

enter image description here

#-------------------------------------------------------------------------
#comparison 3 (RMSE)

  library(hydroGOF)

  DEM_1_matrix = c(as.matrix(DEM_1))
  DEM_2_matrix = c(as.matrix(DEM_2))

  rmse = rmse(DEM_1_matrix,DEM_2_matrix)
  rmse
  [1] 135.3675 # this is the root mean squared error (RMSE) result.

See @whuber's answer on Comparing two TINs created using ArcGIS for Desktop? for a theoretical insight about this issue.

[GIS] ArcGIS Desktop:Transformations: Geographic vs Projected

The use of geographic vs projected is more up to you and your analysis needs. Based on what I've read so far in this question, it is not necessary for you to project into a UTM zone. However if you aren't consistent in picking the correct type of CRS, you will get alignment issues. Generally speaking (and there are exceptions) if your coordinate units are in degrees, you need a geographic coordinate system. If they are in feet, meters, or anything other than degrees, you need a projected coordinate system.

Part of what may confuse you is that NAD83 is actually a datum, and you can see it as part of the names for both geographic and projected coordinate systems. In ArcGIS's dialog they are grouped into folders. When I search for "North American 83" I get no results, so I'm not sure what options you would have been presented with to choose from. When you first see the dialog there should be two folders - one for Geographic Coordinate Systems and one for Projected Coordinate Systems. If you expand the GCS folder, and then the North America folder under that, you'll see the NAD83 that you want. Note there are several variants as well, such as CORS, which are updates/revisions to the original NAD83 definition. But you can ignore these unless you know you need to use one.

I would go back to your original tabular data and redo the conversion from CSV to points with the Add XY Event tool. In that tool, make sure you select the correct GCS system for the coordinates. Then export those layers to feature classes/shapefiles. At that point each of your four files (one and the three parts of two) should have the correctly defined coordinate system.

If you then add them one at a time to ArcGIS, whichever you add first will become the dataframe projection. All the ones after that will automatically be reprojected on-the-fly to that first CRS. Note that because you're changing datums you'll need to select a transformation as mkennedy mentions. This is separate from picking a CRS, though done via the same initial dialog box (look for the transformations button), and if you don't select a transformation at all or choose the wrong one, that can also result in a misalignment. If things don't line up at this stage, there is an error in which CRS you selected - either because something you selected was wrong or because what they're labeled as in the CSV is wrong.

If everything looks like it lines up correctly, you can then use the Project tool on the non-NAD83 layers to convert them to that CRS. Alternatively, if your dataframe is in NAD83 you can right-click on them and choose export, and in the dialog switch from the source CRS to the current dataframe CRS option. When that's done you should have all four feature classes/shapefiles in NAD 83 and you could merge them or add your xy coordinates for whatever.