arcpy - Extracting value of raster for specific points (with X and Y coord) and write to column in Attribute Table?

I have specific points in shp file and list of rasters in the same coord system and want to extract values from rasters for each point and place these values in column in Attribute Table of shp file of points.

Rasters are separate files.

Is it possible to go from vector to raster?

Shp file pf points:

infc = 'samplepoints1.shp'

code just a start:

import arcpy

infile = r"C:\Users\anna_pavlenko\Desktop\SILVER\biscuit_dnbr.tif"
#create raster object

ftc = arcpy.sa.Raster(infile)
#rows and columns
ftc.height 
ftc.width 
#cell size
cellSize = ftc.meanCellHeight
print 'The cell  size is ', cellSize, 'm'

rasterList = r"C:\Users\biscuit_dnbr.tif"
rstArray = arcpy.RasterToNumPyArray(rasterList)


rows, cols = rstArray.shape
print "number of rows: " ,  rows
print "number of columns: ", cols

for rowNum in xrange(rows):
    for colNum in xrange(cols):
        value = rstArray.item(rowNum, colNum)

If I use

inRasterList = [["biscuit_dnbr.tif", "foothills_dnbr.tif", "dNBR1"] ,     ["biscuit_rdnbr.tif", "foothills_rdnbr.tif", "RdNBR1"],["silver_dnbr.tif", "pony_complex_dnbr.tif", "dNBR2"], ["silver_rdnbr.tif", "pony_complex_rdnbr.tif", "RdNBR2"]] 
#  ["foothills_dnbr.tif","dNBR1"], "foothills_rdnbr.tif", "RdNBR1"],
#  ["pony_complex_dnbr.tif", "dNBR2"], "pony_complex_rdnbr.tif", "RdNBR2"]]

ExtractMultiValuesToPoints(inPointFeatures, inRasterList, "NONE")

This function create values for 1 place (polygon) and extract value for point inside this polygon only and place them in columns: "dNBR1", "dNBR2", "RdNBR1", "RdNBR2".
For another set of rasters it create 0 -zeros for first points and add new columns : "dNBR1_1", "dNBR2_1", "RdNBR1_1", "RdNBR2_" each time.

Answer

If you don't have the spatial analyst extension or want to use Numpy, you can convert your x and y point coordinates to the row and column that is close to that point. The Extract Values tools have many more options, and produce slightly different results than those in the code below. I'm not advocating one method over another. I also don't know if this would be any faster.

import numpy as np

Get information about your raster file:

rast = arcpy.Raster(r"somePath")
desc = arcpy.Describe(rast)

Get the upper left corner of the raster and spatial reference system (you can ignore the spatial reference if all your points and rasters are the same).

ulx = desc.Extent.XMin
uly = desc.Extent.YMax
#Spatialrefernce

sr = desc.spatialReference

Convert the raster to a numpy array

rstArray  = arcpy.RasterToNumPyArray(rast)

Access your Feature Class using an update cursor to add the value. I'm using the projectAs tool to make sure the point is in the same coordinate system as the raster. Calculate the row and col based on how far the point is from the upper left corner. This can very depending on rounding. Then pass this into your numpy array.

with arcpy.da.UpdateCursor("NHDPoint",["SHAPE@","value"]) as uc:
    for row in uc:
        pnt = row[0].projectAs(sr)
        #assuming every point falls to the left and below uperleft corner

        deltaX = pnt.centroid.X - ulx
        deltaY = lly- pnt.centroid.Y
        arow = int(deltaY/rast.meanCellHeight)
        acol = int(deltaX/rast.meanCellWidth)
        row[1] = rstArray[arow,acol]
        uc.updateRow(row)

This code assumes that the points will fall to the right of upper left corner, and below it. It also just used the default rounding when converting the float to an integer.

Another option, which might be closer to what the Extract tools is to use a distance to the center of the cell from your points to determine which value to use.