building detection adds on

.gitignore

@@ -0,0 +1,3 @@
+*/MODEL*
+*/.ipynb_checkpoints
+*/__pycache__

Unet_buildings_detection/Dummy_buildings_dataset_generator.py

@@ -0,0 +1,70 @@
+import torch
+from torch.utils.data import Dataset
+import numpy as np
+import cv2
+import matplotlib.pyplot as plt
+
+
+
+COLORS = ['#F44336',"#E91E63",'#9C27B0','#673AB7','#3F51B5','#2196F3','#03A9F4','#00BCD4','#4CAF50',
+ '#8BC34A','#CDDC39','#FFEB3B','#FFC107','#FF9800','#FF5722']
+
+COLORS = ['#F44336',"#E91E63",'#9C27B0','#673AB7','#3F51B5','#2196F3','#03A9F4','#00BCD4','#4CAF50',
+ '#8BC34A','#CDDC39','#FFEB3B','#FFC107','#FF9800','#FF5722']
+
+
+r_min = 12
+r_max = 36
+line_width_min = 2
+line_width_max = 4
+background_intensity = 30.0 / 255.0
+def hex2rgb(h):
+    h = h.lstrip('#')
+    return tuple(int(h[i:i+2], 16) for i in (0, 2 ,4))
+
+def DrawRandomSquare(img,segments,r_min,r_max,alpha):
+    color = hex2rgb( np.random.choice(COLORS) )
+    t = np.random.random()
+    r = int(t * r_min + (1-t) * r_max)
+    i = int(np.random.random()*img.shape[0])
+    j = int(np.random.random()*img.shape[1])
+    theta = np.pi * np.random.random()
+    ri = r*np.cos(theta)
+    rj = r*np.sin(theta) 
+    pts = [(ri,rj),(-rj,ri),(-ri,-rj),(rj,-ri) ]
+    pts = [(i+y,j+x) for (y,x) in pts]
+    pts = np.array(pts, np.int32)
+    pts = pts.reshape((-1,1,2))
+    canvas = img.copy()
+    cv2.fillPoly(canvas,[pts],color)
+    cv2.fillPoly(segments,[pts],(0,1,0))
+    img = cv2.addWeighted(img, 1.0 - alpha, canvas, alpha, 0, img )
+    box = [min(pts[:,:,0])[0],min(pts[:,:,1])[0], max(pts[:,:,0])[0]-min(pts[:,:,0])[0], max(pts[:,:,1])[0] - min(pts[:,:,1])[0] ]
+    return img,segments,box
+
+def generateSegmentation(canvas_size, n_max, alpha = 0.5):
+    canvas = background_intensity * np.ones((canvas_size,canvas_size,3))
+    segments = np.zeros((canvas_size,canvas_size,2))
+
+    for _ in range(np.random.choice(range(n_max))):
+        canvas,segments,b = DrawRandomSquare(canvas,segments,r_min,r_max,alpha)
+    return canvas,segments
+
+
+
+class SimpleSegmentationDataset(Dataset):
+    """A simple dataset for image segmentation purpose"""
+    def __init__(self, patch_size, n_max, alpha =1.0,virtual_size=1000):
+        self.virtual_size = virtual_size
+        self.patch_size = patch_size
+        self.n_max = n_max
+        self.alpha = alpha
+
+
+    def __len__(self):
+        return self.virtual_size
+
+    def __getitem__(self,idx):
+        x,y= generateSegmentation(self.patch_size, self.n_max, self.alpha)
+        sample = {'input': x, 'groundtruth': y}
+        return sample

Unet_buildings_detection/IOU_computations.py

@@ -0,0 +1,95 @@
+from osgeo import gdal, osr,ogr
+import numpy as np
+import json
+
+def vectorize_raster(geoJsonFileName,array2d,layerName="BuildingID",fieldName="BuildingID"):
+
+
+    memdrv = gdal.GetDriverByName('MEM')
+    src_ds = memdrv.Create('', array2d.shape[1], array2d.shape[0], 1)
+    band = src_ds.GetRasterBand(1)
+    band.WriteArray(array2d)
+
+    dst_layername = "BuildingID"
+    drv = ogr.GetDriverByName("geojson")
+    dst_ds = drv.CreateDataSource(geoJsonFileName)
+    dst_layer = dst_ds.CreateLayer(layerName, srs=None)
+
+    fd = ogr.FieldDefn(fieldName, ogr.OFTInteger)
+    dst_layer.CreateField(fd)
+    dst_field = 0
+
+    gdal.Polygonize(band, None, dst_layer, dst_field, [], callback=None)
+
+
+
+    return
+def predict_score_batch(temporary_fold,batch_y,prediction):
+
+    """
+    Predicts scores considering IoU metrics following https://github.com/SpaceNetChallenge/utilities
+    temporary_fold: path of the folder where to store temporary used geojson files extraced from raster input batch_y
+    batch_y (np array): One hot pixel wise labels batch groundtruth(size_batch x WIDTH X HEIGHT x nb_channels) where nb_channels=2
+               batch_y[:,:,:,0] should represent the background mask
+               batch_y[:,:,:,1] should represent the buildings mask
+    prediction (np array):  One hot pixel wise labels batch predidction (size_batch x WIDTH X HEIGHT x nb_channels) where nb_channels=2
+    return averaged accuracy, f1 score, IoU of the batch
+
+    """
+    tot_score_batch=0
+    tot_f1_score_batch=0
+    tot_ious_batch=0
+
+    for i in range(len(batch_y)):
+        vectorize_raster(temporary_fold+str(i)+'_test_gt.geojson',batch_y[i])
+        vectorize_raster(temporary_fold+str(i)+'_test_pred.geojson',prediction[i])
+        with open(temporary_fold+str(i)+'_test_gt.geojson') as f:
+            geojson_groundtruth = json.load(f)
+        with open(temporary_fold+str(i)+'_test_pred.geojson') as f:
+            geojson_prediction = json.load(f)
+
+
+        M=len(geojson_prediction['features'])
+        N=len(geojson_groundtruth['features'])
+#         print('Image %d: %d predictions proposed and %d groundtruth'%(i,M,N))
+        score=0
+        IOUs_sum=0
+        for feature_pred in geojson_prediction['features']:   
+            IoUs=[]
+            IoUs_accu=[]
+#             print(ogr.CreateGeometryFromJson(json.dumps(feature_pred['geometry'])).GetArea())
+#             print('Polygone')
+
+            for feature_gt in geojson_groundtruth['features']:
+#                 print(ogr.CreateGeometryFromJson(json.dumps(feature_gt['geometry'])).GetArea())
+                poly1=ogr.CreateGeometryFromJson(json.dumps(feature_gt['geometry']))
+                poly2=ogr.CreateGeometryFromJson(json.dumps(feature_pred['geometry']))
+                intersection = poly1.Intersection(poly2)
+                union = poly1.Union(poly2)
+                if intersection is None:
+                    IoUs.append(0.0)
+                else:
+                    IoUs.append(intersection.GetArea()/union.GetArea())
+
+            IoUs=np.asarray(IoUs)
+#             print(IoUs)
+            IoUs_accu=(IoUs>0.5).astype(int)*IoUs
+#             print(IoUs_accu)
+            if (IoUs_accu.size and np.amax(IoUs_accu)>0):
+                index=np.argmax(IoUs_accu)
+#                 print('index %d'%index)
+                geojson_groundtruth['features'].remove(geojson_groundtruth['features'][index])
+#                 print('new size groundtruth %d'%len(geojson_groundtruth['features']))
+                score+=1
+                IOUs_sum+=IoUs[index]
+#         print('score: %f: '%score)
+#         print('IOUs_sum: %f: '%IOUs_sum)
+        tot_ious_batch+=IOUs_sum/N
+        tot_score_batch+=score/N
+        tot_f1_score_batch+=2*score/(M+N)
+    tot_ious_batch/=len(batch_y)
+    tot_score_batch/=len(batch_y)
+    tot_f1_score_batch/=len(batch_y)
+    return tot_score_batch*100,tot_f1_score_batch*100,tot_ious_batch
+
+

Unet_buildings_detection/image_utils.py

@@ -0,0 +1,60 @@
+import numpy as np
+from osgeo import gdal,osr,ogr
+import torch
+from scipy.ndimage.morphology import distance_transform_bf
+
+
+
+def standardize(data):
+    '''
+    Standardize the input data of the network
+    :param data to be standardized (size size_batch x WIDTH x HEIGHT x number of channels) 
+
+    returns data standardized size size_batch x WIDTH x HEIGHT x number of channels 
+
+    '''
+
+    WIDTH=data.shape[1]
+    HEIGHT=data.shape[2]
+    channels=data.shape[3]
+
+
+    mean_t=torch.mean(data.view(len(data)*WIDTH*HEIGHT,channels),0)
+    std_t=torch.std(data.view(len(data)*WIDTH*HEIGHT,channels), 0)
+    data=(data-mean_t)/std_t
+
+    #For normalization 
+    min_t=torch.min(data.view(len(data)*WIDTH*HEIGHT,channels), 0)
+    max_t=torch.max(data.view(len(data)*WIDTH*HEIGHT,channels), 0)
+    data=(data-min_t[0])/((max_t[0]-min_t[0]))
+
+    return data
+
+
+def distance_map_batch(Y,threshold=20,bins=15):
+
+    """
+    Compute the distance map following https://arxiv.org/pdf/1709.05932.pdf 
+    Y: One hot pixel wise labels ( size_batch x WIDTH x HEIGHT x nb_classes) with nb_classes=2 
+               Y[:,:,:,0] should represent the background mask
+               Y[:,:,:,1] should represent the buildings mask
+    threshold: distance threshold
+    bins: number of bins considered for the distance map
+    Default values are computed for resolution of input image of 50cm per pixel
+    return torch distance map for buildings (size_batch x WIDTH x HEIGHT x bins)
+    """
+    Y_dist=[]
+    for i in range(len(Y)):
+        distance=distance_transform_bf(np.asarray(Y)[i,:,:,1],sampling=2)
+        distance=np.minimum(distance,threshold*(distance>0))*(bins-1)/threshold
+        inp=torch.LongTensor(distance)
+        inp_ = torch.unsqueeze(inp, len(distance.shape))
+        one_hot = torch.FloatTensor(distance.shape[0],distance.shape[1], bins).zero_()
+        one_hot.scatter_(len(distance.shape), inp_, 1)
+        one_hot=np.asarray(one_hot)
+        Y_dist.append(one_hot)
+
+    return torch.FloatTensor(np.asarray(Y_dist))
+
+
+

Unet_buildings_detection/readme_buildings.md

@@ -0,0 +1,19 @@
+# Unet buildings detection
+
+## Required libraries installation
+On Linux
+
+```sh
+
+$pip install scipy
+$pip install matplotlib
+$conda install gdal
+$conda install -c ioos rtree 
+$pip install centerline
+$pip install osmnx
+$pip install torch torchvision
+$conda install -c menpo opencv
+```
+
+
+

Unet_buildings_detection/requirements.txt

@@ -0,0 +1,5 @@
+python 3.6
+pytorch
+scipy
+gdal
+cv2