import sys
from osgeo import ogr
from shapely.wkb import loads
from shapely.geometry import *
from random import uniform
import fipscodes
import sqlite3

# Import the module that converts spatial data between formats

sys.path.append("filepath")
from globalmaptiles import GlobalMercator

# Main function that reads the shapefile, obtains the population counts,
# creates a point object for each person by race, and exports to a SQL database.

def load_geo_race(state):
    geo = {}
    for line in open("%spop.txt" % state):
        l = line.strip() 
        l = l.split(",")
        geo[l[0]] = map(int, l[2:])
    return geo

def main(input_filename, state):
    geo = load_geo_race(state)
    # Create a GlobalMercator object for later conversions
    
    merc = GlobalMercator()

    # Open the shapefile
    #print input_filename
    ds = ogr.Open(input_filename)
    
    if ds is None:
        print "Open failed.\n"
        sys.exit( 1 )

    # Obtain the first (and only) layer in the shapefile
    
    lyr = ds.GetLayerByIndex(0)

    lyr.ResetReading()

    # Obtain the field definitions in the shapefile layer

    feat_defn = lyr.GetLayerDefn()
    field_defns = [feat_defn.GetFieldDefn(i) for i in range(feat_defn.GetFieldCount())]

    # Obtain the index of the field for the count for whites, blacks, Asians, 
    # Others, and Hispanics.
    
    for i, defn in enumerate(field_defns):
        
        if defn.GetName() == "GEOID20":
            geo_field = i

    # Set-up the output file
    
    # Obtain the number of features (Census Blocks) in the layer
    
    n_features = len(lyr)

    # Iterate through every feature (Census Block Ploygon) in the layer,
    # obtain the population counts, and create a point for each person within
    # that feature.

    for j, feat in enumerate( lyr ):
        
        # Print a progress read-out for every 1000 features and export to hard disk
        
        #if j % 1000 == 0:
        #    conn.commit()
        #    print "%s/%s (%0.2f%%)"%(j+1,n_features,100*((j+1)/float(n_features)))
            
        # Obtain total population, racial counts, and state fips code of the individual census block

        g = feat.GetField(geo_field)
        hispanic, white, black, asian, other = geo[g]
        statefips = fipscodes.state[state.upper()]
        # Obtain the OGR polygon object from the feature

        geom = feat.GetGeometryRef()
        
        if geom is None:
            continue
        
        # Convert the OGR Polygon into a Shapely Polygon
        
        poly = loads(geom.ExportToWkb())
        
        if poly is None:
            continue        
            
        # Obtain the "boundary box" of extreme points of the polygon

        bbox = poly.bounds
        
        if not bbox:
            continue
     
        leftmost,bottommost,rightmost,topmost = bbox
    
        # Generate a point object within the census block for every person by race
       
        race_list = [(white, 'w'), (black, 'b'), (asian, 'a'), (hispanic, 'h'), (other, 'o')]
        for r in race_list:
            
            for i in range(r[0]):
                    
                # Choose a random longitude and latitude within the boundary box
                # and within the orginial ploygon of the census block
                    
                while True:
                        
                    samplepoint = Point(uniform(leftmost, rightmost),uniform(bottommost, topmost))
                        
                    if samplepoint is None:
                        break
                    
                    if poly.contains(samplepoint):
                        break
        
                # Convert the longitude and latitude coordinates to meters and
                # a tile reference
        
                x, y = merc.LatLonToMeters(samplepoint.y,samplepoint.x)
                tx,ty = merc.MetersToTile(x, y, 21)
                    
                # Create a unique quadkey for each point object
                    
                quadkey = merc.QuadTree(tx, ty, 21)
                    
                # Create categorical variable for the race category
                       
                race_type = r[1]         
        
                # Export data to the database file
                print ",".join(map(str, (statefips, int(x),int(y),quadkey,race_type) ))

# Execution code...

if __name__=='__main__':
    state = sys.argv[1].lower()
    fips = fipscodes.state[state.upper()]    
    main("/vsizip/shapes/tl_2020_%s_tabblock20.zip" % fips, state)    
#    for state in ['01','02','04','05','06','08','09','10','11','12','13','15',
#        '16','17','18','19','20','21','22','23','24','25','26','27','28','29','30',
#        '31','32','33','34','35','36','37','38','39','40','41','42','44','45','46',
#        '47','48','49','50','51','53','54','55','56']:
#        
#        print "state:%s"%state
#        
#        main( ".../Census Data/statefile_"+state+".shp", 
#        ".../Map Data/people_by_race5.db")