<jupyter_start><jupyter_text>Demo of CW-Tiler - Atlanta - V3This is a walk through for creation of the Atlantas Dataset.  This demo will walk through the use of the tiler to create utm tiles with the SpaceNet Data Repository.  We will be taking advantage of cloud optimized geotiffs.  For more information about SpaceNet visit https://spacenetchallenge.github.io/<jupyter_code># Import base tools

## Note, for mac osx compatability import something from shapely.geometry before importing fiona or geopandas
## https://github.com/Toblerity/Shapely/issues/553  * Import shapely before rasterio or fioana
from shapely import geometry
import rasterio
import random
from cw_tiler import main
from cw_tiler import utils
from cw_tiler import vector_utils
import numpy as np
import geopandas as gpd
import os
from tqdm import tqdm
# Setting Certificate Location for Ubuntu/Mac OS locations (Rasterio looks for certs in centos locations)  This is important for s3 access  (s3 must be configured using aws cli / boto3)
os.environ['CURL_CA_BUNDLE']='/etc/ssl/certs/ca-certificates.crt'
from stac_tools.stac_item import spacenetStacItem
import glob
import os<jupyter_output><empty_output><jupyter_text>Create List of Images for Processing<jupyter_code>from stac_tools.stac_item import spacenetStacItem
import glob
import os
baseLocation = "/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/srcData/rasterData/"
workorder_list = glob.glob(os.path.join(baseLocation, '**', "*.tif"), recursive=True)

spacenetImageList = []
imageDictList = []
stac_list = []
from shapely.geometry import shape 

for workorder_item in workorder_list:
    item_basename = os.path.splitext(os.path.basename(workorder_item))[0]
    item_folder = os.path.dirname(workorder_item)
    raster_path = workorder_item
    image_id = os.path.basename(os.path.dirname(item_folder))
    ## Note adding extra insertion
    image_id = image_id.replace('nadir', 'nadir_')
    xml_path = os.path.join(item_folder, item_basename+".XML")
    imd_path = os.path.join(item_folder, item_basename+".IMD")
    

    
    item = spacenetStacItem(raster_path, "DigitalGlobe", "CC 4.0 BY-SA", 'ID1', assetDict=[], imdPath=xml_path)
    #print(shape(item.geometry).area/(1000*1000))
    imageDict = {'image_id': image_id, 
                 'image_path': workorder_item,
                 'xml_path': xml_path, 'imd_path': imd_path,
                'image_type': item_basename.rsplit('_',1)[1],
                'image_area':shape(item.geometry).area/(1000*1000),
                'geometry': shape(item.geometry)}
    imageDictList.append(imageDict)
    
    
    
    stac_list.append(item)

print(len(imageDictList))
import geopandas as gpd
imageDF = gpd.GeoDataFrame(imageDictList)


imageDF.head()
imageDF.loc[imageDF['image_type']=='assembly','image_type']='Pan-Sharpen'
imageDF[imageDF['image_id']=='rasterData']['image_path'].values

imageDF.to_file("/raid/nfs/workingDir/cw-tiler/Atlanta_data.shp")<jupyter_output><empty_output><jupyter_text>Create Vector File Atlanta Dataset for tiling extent<jupyter_code># open s3 Location
### This Shapefile will be fed to cw-tiler to create tiles.  
tile_index_shp = '/raid/nfs/workingDir/cw-tiler/AOI_6_Atlanta_Cells_v1.shp'

# Each grid starting point will be spaced 400m apart
stride_size_meters = 450

# Each grid cell will be 400m on a side
cell_size_meters   = 450

# Specify the number of pixels in a tile cell_size_meters/tile_size_pixels == Pixel_Size_Meters
tile_size_pixels   = 900

# Generate list of cells to read from utm_bounds 
cells_list = main.calculate_analysis_grid(aoi_geom.bounds, stride_size_meters=stride_size_meters, cell_size_meters=cell_size_meters)

## 
cells_list = cells_list[0]
    
# Convert to Shapefile
#shapely.geometry.box(minx, miny, maxx, maxy, ccw=True)
from shapely import geometry
geom_list = [geometry.box(*cell_selection) for cell_selection in cells_list]
cell_gdf = gpd.GeoDataFrame(geometry=geom_list)
cell_gdf.crs = utm_crs
cell_gdf.to_file(tile_index_shp)<jupyter_output><empty_output><jupyter_text>Read List of Images into GeoDataframe<jupyter_code>import geopandas as gpd

imageDF = gpd.read_file("/raid/nfs/workingDir/cw-tiler/tmpGeoFiles/Atlanta_data.shp")
imageDF.columns
print("Total Images = {}".format(imageDF.shape))
imageDF.head()<jupyter_output>Total Images = (81, 7)<jupyter_text>Create JSON Work order Files for Processing Raw Imagery into Tiles.<jupyter_code>### 
# Feed index file to chip all other images into.  The Shape File has a set of cells.  CW-Tiler will clip imagery to each cell extent.  
tile_index_shp = '/raid/nfs/workingDir/cw-tiler/tmpGeoFiles/AOI_6_Atlanta_Cells_v1.shp'

# Pick location for Cells
results_data_base_location = '/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav4/'

## Location to write Command Script to:
command_script_location = '/raid/nfs/workingDir/cw-tiler/Atlanta_AOI_parallelScript_v1.sh'

import json
commandScript = []
## Python promt calling cw_tiler_AOI_6_Atlanta.py to process
pythonCommand = 'python cw_tiler_AOI_6_Atlanta.py {}'
## Iterate through each image to create the work order command.  This includes the creation of a json file with the details of the work order.  
for idx, image_row in imageDF.iterrows():
        workorder = image_row.to_dict()
        workorder['geometry']=workorder['geometry'].to_wkt()
        if workorder['image_type'] == 'Pan-Sharpen':
            ## Pan Sharpen tile should bee 900 x 900px
            workorder['tile_pixel_size'] = 900
        elif workorder['image_type'] == 'PAN':
            ## Pan  tile should bee 900 x 900px
            workorder['tile_pixel_size'] = 900
        elif workorder['image_type'] == 'MS':
            ## Pan  tile should bee 225 x 225px
            workorder['tile_pixel_size'] = 225
            
        workorder['tile_index_shp'] = tile_index_shp
        workorder['results_data_base_location'] = results_data_base_location
        
        ## Worker Orders should work on a _cog for quicker tiling.  This is not required
        workorder['image_path'] = workorder['image_path'].replace('.tif', '_cog.tif')
        
        ## Save Work order
        workorder_location = os.path.join(results_data_base_location, "{}_{}_v1.json".format(workorder['image_id'], workorder['image_type']))
        
        #print(workorder_location)
        with open(workorder_location, 'w') as fp:
            json.dump(workorder, fp, indent=4)
            
        commandScript.append(pythonCommand.format(workorder_location))

with open(command_script_location, 'w') as fp:
    
    fp.writelines(["{}\n".format(item)  for item in commandScript])
        
 
## Create Train /Test / Split
import os
import json
import glob
tile_index_shp = '/raid/nfs/workingDir/cw-tiler/AOI_6_Atlanta_Cells_v1.shp'
results_data_base_location = '/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav4/'

results = glob.glob(os.path.join(results_data_base_location, '**', '*.shp'), recursive=True)
## Each Imagery Workfile generates a ShapeFile with results.  
## Process Vector File
import geopandas as gpd
import pandas as pd
import numpy as np
baseIndexFile = gpd.read_file(tile_index_shp)
baseIndexFile1 = baseIndexFile.copy()

## Analyze Results to check all tiles are made.  
for result in results:
    if os.path.basename(result) == 'Pan-Sharpenworkorder.shp':
        tmpFile = gpd.read_file(result)
        idName = os.path.basename(os.path.dirname(result))
        tmpFile[idName] = tmpFile['result']
        tmpFile = tmpFile.drop(columns=['result', 'geometry'])
        baseIndexFile = pd.merge(baseIndexFile, tmpFile, how='outer', on='FID')
    if os.path.basename(result) == 'spacenet-buildingsworkorder.shp':
        tmpFile = gpd.read_file(result)
        #Index(['FID', 'object_cou', 'object_are', 'geometry'], dtype='object')
        idName = 'spacenet-buildings'
        tmpFile['spacenet-buildings_count'] = tmpFile['object_cou']
        tmpFile['spacenet-buildings_area'] = tmpFile['object_are']
        tmpFile = tmpFile.drop(columns=['object_cou', 'object_are','geometry'])
        baseIndexFile = pd.merge(baseIndexFile, tmpFile, how='outer', on='FID')
    
    if os.path.basename(result) == 'spacenet-tile-idworkorder.shp':
        tmpFile = gpd.read_file(result)
        #Index(['FID', 'object_cou', 'object_are', 'geometry'], dtype='object')
        tmpFile['label_suffix'] = tmpFile['tile_index']
        tmpFile['image_suffix'] = tmpFile['tile_ind_1']
        tmpFile = tmpFile.drop(columns=['tile_index', 'tile_ind_1','geometry'])
        baseIndexFile = pd.merge(baseIndexFile, tmpFile, how='outer', on='FID')
        

baseIndexFile.crs = baseIndexFile1.crs
#baseIndexFile.to_file('/raid/nfs/workingDir/cw-tiler/AOI_6_Atlanta_PAN_details.geojson', driver='GeoJSON')
baseIndexFile.head()
baseIndexFile['_TotalCount'] = baseIndexFile.drop(columns=['FID', 'geometry', 'spacenet-buildings_count', 'spacenet-buildings_area', 'label_suffix', 'image_suffix']).sum(numeric_only=True, axis=1)
#Split Train/Test/Validate
TrainSplit = .50
TestSplit  = .40
ValidateSplit = .10
totalTiles = np.sum(np.sum(baseIndexFile.drop(columns=['FID', 'geometry', "_TotalCount"])==0))
totalSpatialTiles = baseIndexFile.shape[0]
print("Total Tiles Spatially = {}".format(totalSpatialTiles))
print("Total Tiles Total = {}".format(totalTiles))
print("Train Tiles = {}, {}% of Tiles".format(totalSpatialTiles*TrainSplit, TrainSplit*100))
print("Test Tiles = {}, {}% of Tiles".format(totalSpatialTiles*TestSplit*1, TestSplit*100))
print("Validate Tiles = {}, {}% of Tiles".format(totalSpatialTiles*ValidateSplit, ValidateSplit*100))
## Assign all tiles to Train by default
baseIndexFile['Category'] = 'Train'
## Assign All Tiles that have a Total Count of missing Tiles >0 To Test  This is only neccessary because all imagery does not have the same extent.  \
## Therefore along the edges there are some nadir angles that don't have complete tiles.  These are passed to test because we only need one tile from each
baseIndexFile.loc[baseIndexFile['_TotalCount']>0, 'Category']='Test'
baseIndexFile['Category'].value_counts()
## Assign more tiles to Test based on the requires split
TestIndex = np.random.choice(baseIndexFile[baseIndexFile['Category']=='Train']['FID'].values, np.int(np.round(totalSpatialTiles*TestSplit - float(baseIndexFile['Category'].value_counts()['Test']))), replace=False)
print(len(TestIndex))    
baseIndexFile.iloc[TestIndex, baseIndexFile.columns.get_loc("Category")]='Test'
baseIndexFile['Category'].value_counts()
## Assign Tiles to Validate
TestIndex = np.random.choice(baseIndexFile[baseIndexFile['Category']=='Train']['FID'].values, int(totalSpatialTiles*ValidateSplit), replace=False)
print(len(TestIndex))    
baseIndexFile.iloc[TestIndex, baseIndexFile.columns.get_loc("Category")]='Validate'
baseIndexFile['Category'].value_counts()
#baseIndexFile.to_file("/raid/nfs/workingDir/cw-tiler/AOI_6_Atlanta_Train_Test_Validate_Split.geojson", driver='GeoJSON')
## Split Test into 3 categories based on Nadir Angle
TestIndex = np.random.choice(baseIndexFile[baseIndexFile['Category']=='Test']['FID'].values, int(totalSpatialTiles*TestSplit/3), replace=False)
print(len(TestIndex))    
baseIndexFile.iloc[TestIndex, baseIndexFile.columns.get_loc("Category")]='Test_Very-Off-Nadir'
TestIndex = np.random.choice(baseIndexFile[baseIndexFile['Category']=='Test']['FID'].values, int(totalSpatialTiles*TestSplit/3), replace=False)
print(len(TestIndex))    
baseIndexFile.iloc[TestIndex, baseIndexFile.columns.get_loc("Category")]='Test_Off-Nadir'
print(len(TestIndex))    
baseIndexFile.loc[baseIndexFile['Category']=='Test', "Category"]='Test_Nadir'
baseIndexFile['Category'].value_counts()




## Write Resulting geojson File
baseIndexFile.to_file("/raid/nfs/workingDir/cw-tiler/AOI_6_Atlanta_Train_Test_Validate_Split_V3.geojson", driver='GeoJSON')
baseIndexFile.groupby("Category")['spacenet-buildings_count', 'spacenet-buildings_area'].describe().reset_index()
### Read Final Tile Split details

final_Tile_Split = gpd.read_file("/raid/nfs/workingDir/cw-tiler/AOI_6_Atlanta_Train_Test_Validate_Split_V4.geojson", driver='GeoJSON')

print(final_Tile_Split['Category'].value_counts())
final_Tile_Split[final_Tile_Split['Category']=='Train'][['label_suffix', 'image_suffix']].head()
## Move files from initial processing into seperate Train/Test/Validate Split using mv command
import shutil
from tqdm import tqdm
results_data_base_location = '/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav2/'
processedData_path = '/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav3/'
commandPathList = []
for idx, row in tqdm(final_Tile_Split[final_Tile_Split['Category']=='Train'].iterrows()):
    results = glob.glob(os.path.join(results_data_base_location, '**', '*{}'.format(row['image_suffix'])), recursive=True)
    for result in results:
        newPath = os.path.join(processedData_path, 'spacenet_4_Train', result.replace('/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav2/', ''))
        if not os.path.exists(os.path.dirname(newPath)):
            os.makedirs(os.path.dirname(newPath))
        commandPath = "mv {} {}".format(result, newPath)
        commandPathList.append(commandPath)
        


## Makes sure folder exists
for pathName in tqdm(list(set([os.path.dirname(item.split(' ')[2]) for item in commandPathList]))):
    if not os.path.exists(pathName):
            os.makedirs(pathName)
        
        
with open('/raid/nfs/workingDir/cw-tiler/Atlanta_AOI_mv_Train_Data.sh', 'w') as fp:
    fp.writelines(["{}\n".format(item)  for item in commandPathList])
    
# To Execute command use sh ./Atlanta_AOI_mv_Train_Data.sh
  
## Work with clipping Labels.  
import shutil
from tqdm import tqdm
results_data_base_location = '/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav2/'
processedData_path = '/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav3/'
commandPathList = []
for idx, row in tqdm(final_Tile_Split[final_Tile_Split['Category']=='Train'].iterrows()):
    results = glob.glob(os.path.join(results_data_base_location, '**', '*{}'.format(row['label_suffix'])), recursive=True)
    for result in results:
        newPath = os.path.join(processedData_path, 'spacenet_4_Train', result.replace('/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav2/', ''))
        if not os.path.exists(os.path.dirname(newPath)):
            os.makedirs(os.path.dirname(newPath))
        commandPath = "mv {} {}".format(result, newPath)
        commandPathList.append(commandPath)
with open('/raid/nfs/workingDir/cw-tiler/Atlanta_AOI_mv_Train_Data_label.sh', 'w') as fp:
         fp.writelines(["{}\n".format(item)  for item in commandPathList]) 
for pathName in tqdm(list(set([os.path.dirname(item.split(' ')[2]) for item in commandPathList]))):
    if not os.path.exists(pathName):
            os.makedirs(pathName)
## Creation of Test Data Set

import pandas as pd
Off_nadir=25
Very_Off_nadir=40

columnList = [column for column in  final_Tile_Split.columns if column[0:3]=='Atl']
AOI_List_Dict = []
for column in columnList:
    nameParts = column.split('_')
    nadirAngle = int(nameParts[1].replace('nadir', ""))
    if nadirAngle <= 25:
        nadirClass = 'Nadir'
    elif nadirAngle < 40:
        nadirClass = 'Off-Nadir'
    else:
        nadirClass = 'Very-off-Nadir'
    
    AOI_List_Dict.append({'FolderName': column,
     'nadir_angle': nadirAngle,
    'nadir_class': nadirClass})
    
        
final_Tile_Split.columns[3][0:3]
aoi_Dict = pd.DataFrame(AOI_List_Dict)
aoi_Dict.head()
aoi_Dict['nadir_class'].value_counts()
NadirList = aoi_Dict[aoi_Dict['nadir_class']=='Nadir']['FolderName'].values
OffNadirList = aoi_Dict[aoi_Dict['nadir_class']=='Off-Nadir']['FolderName'].values
VeryOffNadirList = aoi_Dict[aoi_Dict['nadir_class']=='Very-off-Nadir']['FolderName'].values
final_Tile_Split = gpd.read_file("/raid/nfs/workingDir/cw-tiler/AOI_6_Atlanta_Train_Test_Validate_Split_V4.geojson", driver='GeoJSON')
final_Tile_Split['Category'].value_counts()
import shutil
from tqdm import tqdm
import glob
results_data_base_location = '/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav2/'
processedData_path = '/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav3/'
commandPathList = []

NadirList = aoi_Dict[aoi_Dict['nadir_class']=='Nadir']['FolderName'].values
OffNadirList = aoi_Dict[aoi_Dict['nadir_class']=='Off-Nadir']['FolderName'].values
VeryOffNadirList = aoi_Dict[aoi_Dict['nadir_class']=='Very-off-Nadir']['FolderName'].values

testCategoryList = ["Test_Nadir","Test_Off-Nadir", "Test_Very-Off-Nadir"]
nadirCategoryListList = [NadirList, OffNadirList, VeryOffNadirList]
rowList = []

for testCategory, nadirCategoryList in zip(testCategoryList, nadirCategoryListList):
    for idx, row in tqdm(final_Tile_Split[final_Tile_Split['Category']==testCategory].iterrows()):

        
        rowTmp = row[nadirCategoryList]
        baseID = np.random.choice(rowTmp[rowTmp>=450].index,1)[0]
        print(baseID)
        row['TestimageID'] = baseID


        results = glob.glob(os.path.join(results_data_base_location, baseID, '**', '*{}'.format(row['image_suffix'])), recursive=True)
        for result in results:
            newPath = os.path.join(processedData_path, 'spacenet_4_Test', result.replace('/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav2/', ''))
            if not os.path.exists(os.path.dirname(newPath)):
                os.makedirs(os.path.dirname(newPath))
            commandPath = "cp {} {}".format(result, newPath)
            commandPathList.append(commandPath)
            

        results = glob.glob(os.path.join(results_data_base_location, '**', '*{}'.format(row['label_suffix'])), recursive=True)
        for result in results:
            newPath = os.path.join(processedData_path, 'spacenet_4_Test', result.replace('/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav2/', ''))
            if not os.path.exists(os.path.dirname(newPath)):
                os.makedirs(os.path.dirname(newPath))
            commandPath = "cp {} {}".format(result, newPath)
            commandPathList.append(commandPath)

        rowList.append(row)

    
final_test_train_split_update = gpd.GeoDataFrame(rowList)
final_test_train_split_update.crs = final_Tile_Split.crs
#final_test_train_split_update.to_file("/raid/nfs/workingDir/cw-tiler/AOI_6_Atlanta_Train_Test_Validate_Split_V3p3.geojson", driver='GeoJSON')

        
with open('/raid/nfs/workingDir/cw-tiler/Atlanta_AOI_mv_Test_Data_v4.sh', 'w') as fp:
         fp.writelines(["{}\n".format(item)  for item in commandPathList]) 
        
for pathName in tqdm(list(set([os.path.dirname(item.split(' ')[2]) for item in commandPathList]))):
    if not os.path.exists(pathName):
            os.makedirs(pathName)
import shutil
from tqdm import tqdm
results_data_base_location = '/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav2/'
processedData_path = '/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav3/'
commandPathList = []


rowList = []


for idx, row in tqdm(final_Tile_Split[final_Tile_Split['Category']=='Validate'].iterrows()):
    results = glob.glob(os.path.join(results_data_base_location, '**', '*{}'.format(row['image_suffix'])), recursive=True)
    for result in results:
        newPath = os.path.join(processedData_path, 'spacenet_4_Validate', result.replace('/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav2/', ''))
        if not os.path.exists(os.path.dirname(newPath)):
            os.makedirs(os.path.dirname(newPath))
        commandPath = "mv {} {}".format(result, newPath)
        commandPathList.append(commandPath)

    results = glob.glob(os.path.join(results_data_base_location, '**', '*{}'.format(row['label_suffix'])), recursive=True)
    for result in results:
        newPath = os.path.join(processedData_path, 'spacenet_4_Validate', result.replace('/raid/nfs/data/Datasets/CosmiQ_General_Study/AOI_6_Atlanta/processedDatav2/', ''))
        if not os.path.exists(os.path.dirname(newPath)):
            os.makedirs(os.path.dirname(newPath))
        commandPath = "mv {} {}".format(result, newPath)
        commandPathList.append(commandPath)
         

        
with open('/raid/nfs/workingDir/cw-tiler/Atlanta_AOI_mv_Validate_Data_v2.sh', 'w') as fp:
         fp.writelines(["{}\n".format(item)  for item in commandPathList]) 
        
for pathName in tqdm(list(set([os.path.dirname(item.split(' ')[2]) for item in commandPathList]))):
    if not os.path.exists(pathName):
            os.makedirs(pathName)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13055408462908583, 'fraction_numerical': 0.00914648822991383, 'mean_word_length': 4.876661951909477, 'pattern_counts': {'":': 0, '<': 22, '<?xml version=': 0, '>': 23, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 4}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '159324', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 7432, 'n_tokens_neox': 7090, 'n_words': 1247, 'path': 'cw_tiler_Atlanta-V3.ipynb', 'repo_name': 'CosmiQ/cw-tiler', 'size': '20773'}