"...how can I add the latitude and longitude of each station on my code"
Since you are using the shaperead argument pair {'UseGeoCoords',true},
data = shaperead('tl_2019_us_county.shp','UseGeoCoords',true, . . .);
data will be a structure array including fields lon and lat. To add additional latitude and longitude coordinates,
S(end+1).Lon = 32.52701484484009; % adds an additional structure to array
S(end).Lat = -111.07177734375; % adds the lat value to the new structure
% Add any other values to the new structure...
If the goal is to add those coordinates to the map, you can do that without altering the shaperead data.
geoshow(data);
hold on
plot(-111.07177734375, 32.52701484484009,'bo')
