Error in storing RGB image pixel-values at point locations using rasterio and geopandas

Question

I have a shapefile showing the locations of some cities and I am trying to extract pixel values from a RGB image at the location of points and store it in new columns in the point shapefile. I can print the extracted pixels values, but once I want to store them in new columns in the shapefile, I receive an error. can anyone help me to solve this using rasterio and geopandas?

Here is the code:

import geopandas as gpd
import rasterio as rio
nbands=3
cities=gpd.read_file('path to the point shapefile/cities.shp')
fields=['band1','band2','band3']
cities[fields]=0

for index, row in cities.iterrows():
    cities = row['cities']
    x= row['geometry'].x
    y= row['geometry'].y
    
    raster = rio.open('path to RGB image/stacked_RGB.tif') 
    raster_data= raster.read()
    for i in range(nbands):
        band = raster_data[:,:,i]
        rowIndex,colIndex = raster.index(x,y)
        print(cities + ': '+ str(raster_data[i][rowIndex,colIndex])) #up to here the code works fine
        cities[['band1','band2','band3']].loc[index]=raster_data[i][rowIndex,colIndex]

error:

See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy
  field_samples[['band1','band2','band3']].iloc[index]=raster_data[i][rowIndex,colIndex]
/run/nvme/job_12710378/tmp/ipykernel_144365/871553810.py:14: SettingWithCopyWarning: 
A value is trying to be set on a copy of a slice from a DataFrame

Answer 1

You're overcomplicating things. Just use the rasterio dataset.sample method.

E.g

import geopandas as gpd
import rasterio as rio

cities = gpd.read_file('cities.shp')
stacked_RGB = 'stacked_RGB.tif'

with rio.open(stacked_RGB) as raster: # Only open raster once, not every loop iteration
                                      # and definitely don't read() the raster every iteration
    fields = [f'band{i+1}' for i in range(raster.count)]  # calculate output fields instead of hardcoding
    cities[fields] = list(raster.sample(zip(cities['geometry'].x, cities['geometry'].y)))
    print(cities.head())

Here's an alternate implementation from the geopandas docs.