# import raster image
data(orforest)
orforest <- terra::unwrap(orforest)
# find the 2nd order least squares polynomial surface
polyfit <- fitplane(orforest, order = 2)
# create raster of polyfit
x <- terra::setValues(orforest, polyfit)
# plot the fit
terra::plot(x)
Run the code above in your browser using DataLab