Best methods for extracting raster cell values to points from Mosaic Dataset

[adamwade1]

I have roughly 50 mosaic datasets with 50 corresponding point feature classes. My goal is to append the raster cell value (elevation) to the overlying point. The mosaic datasets are all from a single band rasters (elevation).  I've been trying to iterate the process in arcpy, but I've been having mixed results. I'm running ArcGIS Desktop 10.6.1, and all the feature classes and mosaic datasets are in the same projection. When I run my script I'll get the results I want for the first few iterations, but then the remaining point feature classes contain extracted values that don't exist within any of the rasters. 

I'm not sure on the best method to extract values from a mosaic dataset. The one way that seems to work is to make a a mosaic layer from the mosaic dataset. I've been creating the mosaic layer in the same geodatabase where the point feature classes are stored and then running extraction. 

import arcpy
from arcpy import env
import os

#create a list of all folders in the directory
directory = "path to folders with rasters and mosaic datasets"
folders = os.listdir(directory)

#create list of the point feature classes
env.workspace = r'path to point feature classes'
files = arcpy.ListFeatureClasses()


point_path = "path to point feature classes" #for use in creating the mosaic layer
arcpy.CheckOutExtension("Spatial")
for i, j in zip(folders, files):
    in_workspace = i+"path to each geodatabase"
    in_mosaic_dataset_name = "mosaic dataset name"
    path_to_raster = in_workspace+in_mosaic_dataset_name
    outlayername = point_path+in_mosaic_dataset_name
    arcpy.MakeMosaicLayer_management(path_to_raster,outlayername)
    arcpy.sa.ExtractMultiValuesToPoints(j, outlayername)‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

I've run this both in the python window of ArcMap and in spyder, and in both instances I run into problems. The first few iterations work fine and then the raster cell values that are extracted do not exist in any of the rasters. Is there a better way to extract raster cell values from a mosaic dataset?

Cheers!

[DavidPike]

Maybe values from secondary rasters are being applied? When you create the mosaic layer you could use a query. There also might be overlap operations not being controlled?

[adamwade1]

It seems to happen starting with the third iteration. The first two iterations give me expected values, and then the problem starts. If I rerun the code starting at the previous third iteration I'll again get the expected values for two more iterations and then it goes bad again. So it's not a problem with rasters. How would I check for overlap operations?

[DavidPike]

Looking at the help it seems like the extract multi value to point append the raster values to the input points? Negating any iteration. Unless I'm looking at old syntax.

The script seems to be trying to create point outputs for each raster if I'm reading it right? 

Apologies if I'm reading it wrong.

[adamwade1]

The script does an extraction for about 50 different point featureclasses over 50 different mosaic datasets. It does exactly what I need for two iterations and then doesn't. I've had similar problems doing the same operation using the toolbox. I'll run the tool using a point feature class and mosaic dataset and get completely unexpected values. I'll run it again with the same data and it will work fine. It must be something with reading the mosaic datasets. It's driving me nuts. 

[DavidPike]

Haha! What about running it through the source rasters of your mosaic datasets? 

Also the OID field of your points could be getting changed if that's giving misleading info.

[adamwade1]

There are hundreds of rasters in each mosaic, so I would end up with hundreds of new new columns in each point feature class. I only want one column with all the values. I could concatenate all the columns into one, and afterward delete the previous columns but that's a serious amount of extra time, and defeats the purpose of using a mosaic dataset. The processes obviously works, but runs into a problem during the iteration. Just not sure what's causing the problem. 

[adamwade1 · Accepted solution]

I figured out a way to get this to work. For some reason, and I have no idea why, if you run the ExtractMultiValuesToPoints from a mosaic dataset you end up with values not found in the raster cells. I'm still not sure the source of the values. However, if you run it again with the exact same data you will get the desired results. I modified my above script to so that the ExtractMultiValuesToPoints uses the mosaic dataset twice. Also, I determined there is no need to make the mosaic layer, the tool works from the mosaic dataset alone. This solution will give you two new columns in the attribute table; one with the unknown values and the second with the desired values. Depending on how many files you are running you can either manually delete the extra column from each table or just modify the script to delete the first column of unknown values.  

import arcpy
from arcpy import env
import os

#create a list of all folders in the directory
directory = "path to folders with rasters and mosaic datasets"
folders = os.listdir(directory)

#create list of the point feature classes
env.workspace = r'path to point feature classes'
files = arcpy.ListFeatureClasses()

#iterate the point files and mosaic datasets
arcpy.CheckOutExtension("Spatial")
for i, j in zip(folders, files):
    in_workspace = i+"path to each geodatabase"
    in_mosaic_dataset_name = "mosaic dataset name"
    path_to_raster = in_workspace+in_mosaic_dataset_name
    rasterlist = [path_to_raster, path_to_raster] 
    arcpy.sa.ExtractMultiValuesToPoints(j, outlayername)‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

If anyone knows why the initial run gives you the unknown values, or where those values come from, please let me know. 

Cheers!