genpointsbygeo.py

#!/usr/bin/env python

from faker import Factory
fake = Factory.create('ja_JP')

import math
import csv

# shape
import fiona
from shapely.geometry import Point,asShape

## args
import argparse
parser = argparse.ArgumentParser(description='Generater of coordinate points.')
parser.add_argument('--samples', type=int, nargs='?', default=10, help='number of points')
parser.add_argument('--csv', type=str, nargs='?', default="data.csv", help='csv file name')
parser.add_argument('--radius', type=int, nargs='?', default=100, help='radius(m)')
# Tokyo Station(35.681382, 139.766084)
parser.add_argument('--latitude', type=float, nargs='?', default=35.681382, help='center latitude of the circle')
parser.add_argument('--longitude', type=float, nargs='?', default=139.766084, help='center longitude of the circle')
parser.add_argument('--shapefile', type=str, nargs='?', default=None, help='shapefile')
args = parser.parse_args()

##
# 
samples       = args.samples
limit_samples = samples * samples
csvfile       = args.csv
centerlat     = args.latitude
centerlong    = args.longitude
radius_m      = args.radius
shapefile     = args.shapefile

#print args, limit_samples

#########
earthradius = 6378137
lat1radm   = ((2*math.pi*earthradius)/360)
latradius  = radius_m/lat1radm
long1radm  = ((earthradius*math.cos(centerlat/180*math.pi)*2*math.pi)/360)
longradius = radius_m/long1radm

with open(csvfile, "w+") as f:
    csv_writer = csv.writer(f)
    
    counter = 0
    
    if shapefile is not None:
        with fiona.open(shapefile) as fiona_collection:
            shapefile_record = fiona_collection.next()
            shape = asShape(shapefile_record['geometry'])
            
            for _ in range(0,limit_samples):
                geolat  = fake.geo_coordinate(center=centerlat , radius=latradius)
                geolong = fake.geo_coordinate(center=centerlong, radius=longradius)
                
                point = Point(geolong, geolat)
                if shape.contains(point):
                    d = [geolat, geolong]
                    csv_writer.writerow(d)
                    counter = counter + 1
                                
                if counter >= samples:
                    break
    
    # shape file is None
    else:
        for _ in range(0,limit_samples):
            geolat  = fake.geo_coordinate(center=centerlat , radius=latradius)
            geolong = fake.geo_coordinate(center=centerlong, radius=longradius)
            
            r = math.sqrt(math.pow((float(geolat)-centerlat)*lat1radm,2)+math.pow((float(geolong)-centerlong)*long1radm,2))
            if r < radius_m:
                d = [geolat, geolong]
                csv_writer.writerow(d)
                counter = counter + 1

            if counter >= samples:
                break