<jupyter_start><jupyter_text>In this notebook, we reproject all point clouds that have been clipped to lidar plots from the Washington State Plane coordinate reference system over to UTM. The WA DNR Lidar Portal stores all its lidar point clouds in this CRS, but we will be moving everything over to UTM.<jupyter_code>import geopandas as gpd
import pandas as pd
import os
import glob
import subprocess
import dask
from dask.distributed import Client, LocalCluster, progress
PLOTS = '../data/processed/blm_usfs_wadnr_plot_footprints.shp'
plots = gpd.read_file(PLOTS)
plots.info()
plots.crs
UTM_10 = '../data/external/utm_zone10_epsg4326.shp'
UTM_11 = '../data/external/utm_zone11_epsg4326.shp'
utm_10 = gpd.read_file(UTM_10)
utm_11 = gpd.read_file(UTM_11)
utm_10.info()
utm_10.crs
utm10_plots = gpd.sjoin(plots, utm_10)
utm11_plots = gpd.sjoin(plots, utm_11)

utm_lookup = pd.concat((utm10_plots, utm11_plots), 
                       axis=0, 
                       ignore_index=True)[['uuid', 'ZONE']].set_index('uuid')
utm_lookup['epsg'] = utm_lookup.ZONE.apply(lambda x: 6339 if x == '10' else 6340)
utm_lookup.head()
len(utm_lookup)
pd.unique(utm_lookup.epsg)
@dask.delayed
def reproject(infile, odir):
    basename = os.path.basename(infile)
    outfile = os.path.join(odir, basename)
    
    uuid = basename.split('_')[0]
    target_epsg = utm_lookup.loc[uuid]['epsg']
    
    proc = subprocess.run(['wine', '/storage/lidar/LAStools/bin/las2las.exe',
                           '-i', infile,
                           '-o', outfile,
                           '-epsg', '2927',  # input is in NAD83(HARN) / Washington South (ftUS)
                           '-elevation_surveyfeet',  # input has elevation in surveyfeet
                           '-target_epsg', str(target_epsg)],
                          stderr=subprocess.PIPE,
                          stdout=subprocess.PIPE)
    return proc
USFS_HA = glob.glob('../data/raw/lidar/dnr_portal/dnr_portal_tiles/usfs_plots/hectare_clips_epsg2927/*.laz')
USFS_PLOT = glob.glob('../data/raw/lidar/dnr_portal/dnr_portal_tiles/usfs_plots/plot_clips_epsg2927/*.laz')
DNR_HA = glob.glob('../data/raw/lidar/dnr_portal/dnrplot_clips/hectare_clips_epsg2927/*.laz')
DNR_PLOT = glob.glob('../data/raw/lidar/dnr_portal/dnrplot_clips/plot_clips_epsg2927/*.laz')
cluster=LocalCluster()
c = Client(cluster)
jobs_to_do = []
PLOT_SETS = [USFS_HA, USFS_PLOT, DNR_HA, DNR_PLOT]
ODIRS = ['../data/interim/lidar/plot_clips/dnr_portal/usfs_plots/hectare_clips/',
         '../data/interim/lidar/plot_clips/dnr_portal/usfs_plots/plot_clips/',
         '../data/interim/lidar/plot_clips/dnr_portal/dnr_plots/hectare_clips/',
         '../data/interim/lidar/plot_clips/dnr_portal/dnr_plots/plot_clips/']

for plot_set, odir in zip(PLOT_SETS, ODIRS):
    jobs_to_do += [reproject(p, odir=odir) for p in plot_set]
res = c.persist(jobs_to_do)
progress(res)
@dask.delayed
def get_boundary(infile, odir):
    proc = subprocess.run(['wine', '/storage/lidar/LAStools/bin/lasboundary.exe',
                           '-i', infile,
                           '-odir', odir,
                           '-oshp',
                           '-labels'], 
                          stderr=subprocess.PIPE,
                          stdout=subprocess.PIPE)
#     print(outfile)
    return proc
USFS_HA = glob.glob('../data/interim/lidar/plot_clips/dnr_portal/usfs_plots/hectare_clips/*.laz')
USFS_PLOT = glob.glob('../data/interim/lidar/plot_clips/dnr_portal/usfs_plots/plot_clips/*.laz')
DNR_HA = glob.glob('../data/interim/lidar/plot_clips/dnr_portal/dnr_plots/hectare_clips/*.laz')
DNR_PLOT = glob.glob('../data/interim/lidar/plot_clips/dnr_portal/dnr_plots/plot_clips/*.laz')

jobs_to_do = []
PLOT_SETS = [USFS_HA, USFS_PLOT, DNR_HA, DNR_PLOT]
ODIRS = ['../data/interim/lidar/plot_clips/dnr_portal/usfs_plots/hectare_clips/',
         '../data/interim/lidar/plot_clips/dnr_portal/usfs_plots/plot_clips/',
         '../data/interim/lidar/plot_clips/dnr_portal/dnr_plots/hectare_clips/',
         '../data/interim/lidar/plot_clips/dnr_portal/dnr_plots/plot_clips/']

for plot_set, odir in zip(PLOT_SETS, ODIRS):
    jobs_to_do += [get_boundary(p, odir=odir) for p in plot_set]
res = c.persist(jobs_to_do)
progress(res)
# c.close()
# cluster.close()<jupyter_output><empty_output>
StarcoderdataJupyterStructuredCleanDedup