<jupyter_start><jupyter_text>Functions for data processing<jupyter_code>def raster2array(rasterfn):
    '''Converts raster to 2d array'''
    #print('converting raster to array...')
    raster = gdal.Open(rasterfn)
    band = raster.GetRasterBand(1)
    array = band.ReadAsArray()
    return array



    
class RasterProp:
    def __init__(self,
                 rasterFile,
                 sliceClass=None, slicing = False):
        self.raster = gdal.Open(rasterFile)
        self.geotransform = self.raster.GetGeoTransform()
        self.projRef      = self.raster.GetProjectionRef()
        self.originX = self.geotransform[0]
        self.originY = self.geotransform[3] 
        self.pixelWidth = self.geotransform[1] 
        self.pixelHeight = self.geotransform[5]
        
        
def coord2pixelOffset(rasProp,x,y):
    print('coordinate to pixel offsetting...')
    originX = rasProp.originX
    originY = rasProp.originY
    
    pixelWidth = rasProp.pixelWidth 
    pixelHeight = rasProp.pixelHeight

    xOffset = int((x - originX)/pixelWidth)
    yOffset = int((y - originY)/pixelHeight)
    return xOffset,yOffset


def pixel2coord(geoTrans, x, y):
    xoff, a, b, yoff, d, e = geoTrans

    xp = a * x + b * y + a * 0.5 + b * 0.5 + xoff
    yp = d * x + e * y + d * 0.5 + e * 0.5 + yoff
    return(int(xp), int(yp))
def array2raster(array, rasProp,newRasterfn):
    print('converting array to raster...')
    cols = array.shape[1]
    rows = array.shape[0]
    driver = gdal.GetDriverByName('GTiff')
    outRaster = driver.Create(newRasterfn, 
                              cols, rows,
                              bands=1, 
                              eType= gdal.GDT_Float32)

    outRaster.SetGeoTransform((rasProp.originX, 
                               rasProp.pixelWidth, 
                               0, rasProp.originY, 
                               0, rasProp.pixelHeight))
    
    outband = outRaster.GetRasterBand(1)
    outband.WriteArray(array)
    outRasterSRS = osr.SpatialReference()
    outRasterSRS.ImportFromWkt(rasProp.projRef)
    outRaster.SetProjection(outRasterSRS.ExportToWkt())
    outband.FlushCache()
    
def createTotalCostRaster(factorPathList, 
                          weightList, 
                          rasProp,
                          rasterName, 
                          slicing=False, 
                          strPoint=None, 
                          endPoint=None):

        costArray = maxMinScale(raster2array(factorPathList[0]))*weightList[0]
        costArray[np.isnan(costArray)]=0
        
        for fpos in range(1,len(factorPathList)):
            #print(factorNames[fpos])
            factorArray = maxMinScale(raster2array(factorPathList[fpos]))*weightList[fpos]
            factorArray[np.isnan(factorArray)]=0
            #plt.imshow(factorArray)
            costArray = np.add(costArray, factorArray)
        costArray[np.isnan(costArray)]=0
        plt.imshow(costArray)
        plt.colorbar()
        array2raster(costArray, rasProp, rasterName)
        return costArray, rasProp
            
        np.place(costArray, costArray==nan,0)
        array2raster(costArray, sliceRasProp, rasterName)
        plt.imshow(costArray)
        return costArray, sliceRasProp
    
    
def maxMinScale(array):
    return (array/abs(array.max()-array.min()))
def createPath(rasProp, costSurfaceArray,
               startCoord,stopCoord):
    '''returns an array of the same shape as costSurfaceArray with
    1 for path and 0 for other cells'''
    print('creating path...')
    # coordinates to array index
    startCoordX = startCoord[0]
    startCoordY = startCoord[1]
    startIndexX,startIndexY = coord2pixelOffset(rasProp,
                                                startCoordX,
                                                startCoordY)

    stopCoordX = stopCoord[0]
    stopCoordY = stopCoord[1]
    stopIndexX,stopIndexY = coord2pixelOffset(rasProp,
                                              stopCoordX,stopCoordY)
    # create path
    indices, weight = route_through_array(costSurfaceArray, 
                                          (startIndexY,startIndexX), 
                                          (stopIndexY,stopIndexX),
                                          geometric=True,
                                          fully_connected=True)
    
    
    indices = np.array(indices).T
    path = np.zeros_like(costSurfaceArray)
    path[indices[0], indices[1]] = 1
    print('path created...')
    return path
def getStartEndCord(file):
    '''For reading 'start' and 'end' coordindates from shape files - 
    used specifically for DC connection files'''
    
    shape = fiona.open(file)
    first = shape.next()
    strX, strY =first.get('properties').get('CoordX'), first.get('properties').get('CoordY')
    second = shape.next()
    endX, endY =second.get('properties').get('CoordX'), second.get('properties').get('CoordY')
    #return first
    return ((strX,strY) ,(endX,endY))<jupyter_output><empty_output><jupyter_text>reading Dimension rasters and converting to arraysexogenously including city and sea area to give high weights to them<jupyter_code>ecoPath = os.path.abspath('01_Data500/fac_eco_onlySlope.tif')
envPath = os.path.abspath('01_Data500/fac_env.tif')
pubPath = os.path.abspath('01_Data500/fac_pub.tif')
infPath = os.path.abspath('01_Data500/fac_inf.tif')
citPath = os.path.abspath('01_Data500/city.tif')
seaPath = os.path.abspath('01_Data500/sea_coast.tif')

ecoFac = raster2array(ecoPath)
envFac = raster2array(envPath)
pubFac = raster2array(pubPath)
infFac = raster2array(infPath)
citAre = raster2array(citPath)
seaAre = raster2array(seaPath)*100

dcProjects = os.path.abspath('02_DC_Project_DE/')
dc5Path = str(dcProjects+'/DC_5.shp')
dcNorthProjects = os.path.abspath('03_DC_Project_North/')
dcNorPath = str(dcNorthProjects+'/DC_N.shp')
dcNorPath<jupyter_output><empty_output><jupyter_text>Path Creation  Individual Paths Saving full raster for visualization<jupyter_code># environment
envEdgeCost = (8/10)*envFac + \
              (1/10)*pubFac + (1/10)*infFac + citAre + seaAre

envLine = createPath(RasterProp(envPath),
                    costSurfaceArray = envEdgeCost,
                    startCoord= getStartEndCord(dcNorPath)[0],
                    stopCoord= getStartEndCord(dcNorPath)[1])

array2raster(envLine, rasProp=RasterProp(envPath),
            newRasterfn = os.path.abspath('03_DC_Project_North/02_paths/envLine.tif'))
# Acceptance line
pubEdgeCost = (1/10)*envFac + \
              (8/10)*pubFac + (1/10)*infFac + citAre + seaAre

pubLine = createPath(RasterProp(pubPath),
                    costSurfaceArray = pubEdgeCost,
                    startCoord= getStartEndCord(dcNorPath)[0],
                    stopCoord= getStartEndCord(dcNorPath)[1])

array2raster(pubLine, rasProp=RasterProp(pubPath),
            newRasterfn = os.path.abspath('03_DC_Project_North/02_paths/pubLine.tif'))
# Infrastructure line
infEdgeCost = (1/10)*envFac + \
              (1/10)*pubFac + (8/10)*infFac + citAre + seaAre

infLine = createPath(RasterProp(infPath),
                    costSurfaceArray = infEdgeCost,
                    startCoord= getStartEndCord(dcNorPath)[0],
                    stopCoord= getStartEndCord(dcNorPath)[1])

array2raster(infLine, rasProp=RasterProp(infPath),
            newRasterfn = os.path.abspath('03_DC_Project_North/02_paths/infLine.tif'))
# All Equal line
equEdgeCost = (1/3)*envFac + \
              (1/3)*pubFac + (1/3)*infFac + citAre + seaAre

equLine = createPath(RasterProp(infPath),
                     costSurfaceArray = equEdgeCost,
                     startCoord = getStartEndCord(dc5Path)[0],
                     stopCoord  = getStartEndCord(dc5Path)[1])

array2raster(equLine, rasProp=RasterProp(infPath),
            newRasterfn = os.path.abspath('02_DC_Project_DE/02_dc5_paths/equLine.tif'))<jupyter_output>/Users/samarth/miniconda3/envs/gis/lib/python3.7/site-packages/ipykernel_launcher.py:6: DeprecationWarning: Collection.__next__() is buggy and will be removed in Fiona 2.0. Switch to `next(iter(collection))`.
  
/Users/samarth/miniconda3/envs/gis/lib/python3.7/site-packages/ipykernel_launcher.py:8: DeprecationWarning: Collection.__next__() is buggy and will be removed in Fiona 2.0. Switch to `next(iter(collection))`.<jupyter_text>Lines A, B, C from graph  Visualization of the lines unsed in description in the paper<jupyter_code>pointACost = 0.25*envFac + \
             0.75*pubFac + 0*infFac + citAre + seaAre

pointBCost = 0.50*envFac + \
             0.25*pubFac + 0.25*infFac + citAre + seaAre

pointCCost = 0.0*envFac + \
             0.25*pubFac + 0.75*infFac + citAre + seaAre



A_line = createPath(RasterProp(envPath),
                    costSurfaceArray = pointACost,
                    startCoord= getStartEndCord(dc5Path)[0],
                    stopCoord= getStartEndCord(dc5Path)[1])

B_line = createPath(RasterProp(envPath),
                    costSurfaceArray = pointBCost,
                    startCoord= getStartEndCord(dc5Path)[0],
                    stopCoord= getStartEndCord(dc5Path)[1])

C_line = createPath(RasterProp(envPath),
                    costSurfaceArray = pointCCost,
                    startCoord= getStartEndCord(dc5Path)[0],
                    stopCoord= getStartEndCord(dc5Path)[1])


array2raster(A_line, rasProp=RasterProp(envPath),
            newRasterfn = os.path.abspath('02_DC_Project_DE/02_dc5_paths/line_A.tif'))

array2raster(B_line, rasProp=RasterProp(envPath),
            newRasterfn = os.path.abspath('02_DC_Project_DE/02_dc5_paths/line_B.tif'))

array2raster(C_line, rasProp=RasterProp(envPath),
            newRasterfn = os.path.abspath('02_DC_Project_DE/02_dc5_paths/line_C.tif'))<jupyter_output>/Users/samarth/miniconda3/envs/gis/lib/python3.7/site-packages/ipykernel_launcher.py:6: DeprecationWarning: Collection.__next__() is buggy and will be removed in Fiona 2.0. Switch to `next(iter(collection))`.
  
/Users/samarth/miniconda3/envs/gis/lib/python3.7/site-packages/ipykernel_launcher.py:8: DeprecationWarning: Collection.__next__() is buggy and will be removed in Fiona 2.0. Switch to `next(iter(collection))`.<jupyter_text>Sweep Sweep over different weighting factors for the different dimensionsPath created and indexes of the paths saved as npy<jupyter_code>## %%time
allPaths = []
strCoord = getStartEndCord(dcNorPath)[0]
endCoord = getStartEndCord(dcNorPath)[1]
for env in range(0,11,1):
    for inf in range(0,11,1):
        for pub in range(0,11,1):
            if (env!=0 or inf!=0 or pub!=0):
                if (env == inf==pub==0):
                    continue;
                c_env = env/(env+inf+pub)
                c_inf = inf/(env+inf+pub)
                c_pub = pub/(env+inf+pub)
                print([env,inf,pub])
                totalCost = c_env*envFac + c_pub*pubFac + \
                c_inf*infFac + citAre + seaAre
                path = createPath(RasterProp(ecoPath), totalCost,
                                  strCoord, 
                                  endCoord)
                pathidx =np.nonzero(path)
                fileName = str('eip_')+str(env)+str(inf)+str(pub)
                comName = os.path.abspath('03_DC_Project_North/02_paths_2/'+fileName+'.npy')
                np.save(comName,pathidx)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.10979045606620884, 'fraction_numerical': 0.02271526677231907, 'mean_word_length': 2.5134549953603464, 'pattern_counts': {'":': 0, '<': 19, '<?xml version=': 0, '>': 19, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 20, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '471146', 'license': 'MIT', 'n_tokens_mistral': 3839, 'n_tokens_neox': 3634, 'n_words': 611, 'path': '02_PathCreation.ipynb', 'repo_name': 'samarthiith/PowerLine-GIS', 'size': '11358'}