<jupyter_start><jupyter_text>*Special library needed*- pip install git+https://github.com/wpgp/wpgpDownloadPy- pip install gdal - pip install pandasNote gdal may require install gdal OGR in your OS first (apt-get install) google it for your OSAlso included in the requirements.txt file<jupyter_code>from wpgpDownload.utils.wpcsv import ISO_LIST

for iso in ISO_LIST:
    print(iso)
iso_country = 'USA'
from wpgpDownload.utils.wpcsv import Product

products = Product(iso_country)  # Where instead of GRC it could be any valid ISO code.

#  to list all the products for GRC
for p in products:
    print('{}/{}\t{}\t{}'.format(p.idx, p.country_name,p.dataset_name,p.path))

from ftplib import FTP
import os

worldpopmanifest = "/assets/wpgpDatasets.csv"
worldpopftp = "ftp.worldpop.org.uk"
ftp = FTP(worldpopftp)
r = ftp.login()
print("login response : {}".format(r))
filename = '.' + worldpopmanifest
if not os.path.exists(os.path.dirname(filename)):
    try:
        os.makedirs(os.path.dirname(filename))
    except OSError as exc: # Guard against race condition
        if exc.errno != errno.EEXIST:
            raise

ftp.cwd(os.path.dirname(worldpopmanifest))
ftp.retrlines('LIST')
print(os.path.basename(worldpopmanifest))
            
with open(filename, 'wb') as f:
    ftp.retrbinary("RETR {}".format(os.path.basename(worldpopmanifest)), f.write)

ftp.quit()

    
import pandas as pd

df=pd.read_csv(filename)

df.head()
#mport geopandas as gpd
is_USA = df["ISO3"] == "USA"
USA_manifest = df[is_USA]
USA_manifest.head()<jupyter_output><empty_output><jupyter_text>*This step may take about 15 minutes to download the 15Gb tiff for whole US territory*<jupyter_code>list=['ppp_2010','ppp_2010']
#print(USA_manifest["Covariate"]=="ppp_2010")
print(USA_manifest[USA_manifest["Covariate"]=="ppp_2010"]["PathToRaster"])

rasterPath = USA_manifest[USA_manifest["Covariate"]=="ppp_2010"]["PathToRaster"].iloc[0]
print("rasterPath : {}".format(rasterPath))
print("rasterPath type : {}".format(type(rasterPath)))

worldpopftp = "ftp.worldpop.org.uk"
ftp = FTP(worldpopftp)
r = ftp.login()
print("login response : {}".format(r))
filename = './' + rasterPath
print("filename : {}".format(filename))
if not os.path.exists(os.path.dirname(filename)):
    try:
        os.makedirs(os.path.dirname(filename))
    except OSError as exc: # Guard against race condition
        if exc.errno != errno.EEXIST:
            raise
print("folder: {}".format(os.path.dirname(filename)))
ftp.cwd(os.path.dirname(filename))
ftp.retrlines('LIST')
print(os.path.basename(filename))

with open(filename, 'wb') as f:
    ftp.retrbinary("RETR {}".format(os.path.basename(filename)), f.write)

print("dl okay")
ftp.quit()


import gdal
import os
from gdalconst import GA_ReadOnly
geotiff="./GIS/Population/Global_2000_2020/2010/USA/usa_ppp_2010.tif"

buffer_deg=0.0027
xmax_MetroArea = -76.3881880661718 + buffer_deg
ymax_MetroArea = 39.5284754397557 + buffer_deg
xmin_MetroArea = -77.7356609923511 - buffer_deg
ymin_MetroArea = 38.3969644315365 - buffer_deg
#ulx> <uly> <lrx> <lry
extent = [xmin_MetroArea, ymax_MetroArea, xmax_MetroArea, ymin_MetroArea]
print(extent)
data = gdal.Open(geotiff, GA_ReadOnly)
newGeotiff = "{}/DC_{}".format(os.path.dirname(geotiff),os.path.basename(geotiff))
print("geotiff: {}".format(geotiff))
print("newGeotiff: {}".format(newGeotiff))
#geoTransform = data.GetGeoTransform()
#minx = geoTransform[0]
#maxy = geoTransform[3]
#maxx = minx + geoTransform[1] * data.RasterXSize
#miny = maxy + geoTransform[5] * data.RasterYSize
#https://gdal.org/programs/gdal_translate.html
data = None

from osgeo import gdal

data = gdal.Open(geotiff)
data = gdal.Translate(newGeotiff, data, projWin = extent)
data = None
print("done")
import pandas as pd
xls_template_file = "1km_t0.csv"

df_template=pd.read_csv(xls_template_file)

df_template.head()

#https://gis.stackexchange.com/questions/221292/retrieve-pixel-value-with-geographic-coordinate-as-input-with-gdal/221430

from osgeo import gdal

driver = gdal.GetDriverByName('GTiff')
filename = "./GIS/Population/Global_2000_2020/2010/USA/DC_usa_ppp_2010.tif" #newGeotiff #path to raster
dataset = gdal.Open(filename)
band = dataset.GetRasterBand(1)

cols = dataset.RasterXSize
rows = dataset.RasterYSize

transform = dataset.GetGeoTransform()

xOrigin = transform[0]
yOrigin = transform[3]
pixelWidth = transform[1]
pixelHeight = -transform[5]

data = band.ReadAsArray(0, 0, cols, rows)
print("...")
coordinates = df_template[['x', 'y']].apply(tuple, axis=1).tolist()
#get list of 25k coordinates tuples as [(-77.7356609923511, 39.5284754397557), (-77.7266778395099, 39.5284754397557), ...

#import pdb; pdb. set_trace()

#points_list = [ (355278.165927, 4473095.13829), (355978.319525, 4472871.11636) ] #list of X,Y coordinates
pop = []
for point in coordinates:
    col = int((point[0] - xOrigin) / pixelWidth)
    row = int((yOrigin - point[1] ) / pixelHeight)
    #print ("{},{}={}".format(row, col, data[row][col]))
    pop += [data[row][col]]

#import pdb; pdb. set_trace()

df_template["Population"]=pop

df_template.head()
#for idx, line in df_template.iterrows():
#    print(line)
#    print("{}) {}/{} coordinates".format(idx, line["x"], line["y"]))
#    die()<jupyter_output>...
{'fraction_non_alphanumeric': 0.15180446820698873, 'fraction_numerical': 0.04983769333587932, 'mean_word_length': 5.2431466030989275, 'pattern_counts': {'":': 0, '<': 11, '<?xml version=': 0, '>': 11, 'https://': 3, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '2', 'dir': 'jupyter-structured-clean-dedup', 'id': '295118', 'license': 'MIT', 'n_tokens_mistral': 2143, 'n_tokens_neox': 1941, 'n_words': 389, 'path': 'data/worldpop/exploring_population2010.ipynb', 'repo_name': 'loudog40oz/Air-Pollution', 'size': '5237'}