#
#---------------------------------------------------
#
# Geo_Profiling v.1.4.5 - 23 / 01 / 2013
#this version DOES NOT check if points are too close to the border (it is for simulated set of data)
# Authors:
# Alessio Papini - Department of Plant Biology University of Florence Italy, Via La Pira, 4 Firenze, mail alpapiniATunifi.it
# Ugo Santosuosso - Department of Anatomny, Istology and Forensic Medicine, Largo Brambilla, 1 Firenze, mail ugoATunifi.it
# 
# Implements:
# - euclidean distance
# - cvs support
# - test input crimes data are not to close to image border
# - test the correct arguments numbers
# - auto generate output file name
# - preserve elaboration data (usefull for blending images with clustering)
# - Euclidean and Manhattan distance
# - defaults parameters
# Geoprof1_4_4OKcsv.py no weight!!!
#---------------------------------------------------
#
import os, sys
import math
#
# imports libraries for using Matrix, images and csv files
#
from numpy import *
import csv
import Image
import ImageDraw
import ImageFont
# 
# Default parameters
#
NUM_ARGS=3
BLEND_PERCENTAGE=0.7
#
B=6
f=0.4
g=0.4
Dist="E"
#
# Retrives the name of the program from the command line
#
nome=sys.argv[0]
if (len(sys.argv) < NUM_ARGS ):
#
 print" ------------------------------- "
 print" "
 print" ------ "+nome+" --------"
 print" " 
 print" The coordinates must be provided in a simple text file in csv format."
 print" Values on a row separated by a comma"
 print"        First value in row is the sequential number of event"
 print"        Second value in row is X value"
 print"        Third value in row is Y value"
 print"        Fourth value in row is weight value ( 0<weight<=100 )"
 print" Origin of the data (0,0) is upper left corner of the image"
 print" How the program works: insert parameters following the program name: that is "
 print" "
 print" python "+nome+" nameimagemap.bmp coordinates.csv name_prefix_of_results_NO_EXT Dist B f g "
 print" "
 print" REQUIRED:"
 print""
 print" * nameimagemap.bmp coordinates.csv name_prefix_of_results_NO_EXT *"
 print" "
 print" Dist = type of distance E (euclidean-default) M (Manhattan)"
 print""  
 print" Suggested defaults are B = ",B
 print"                        f = ",f
 print"                        g = ",g
 print" "
 print" No path specifications = current directory"
 print" "
 print" ------------------------------- "
 sys.exit (0)
#
# Retrieves the parameters for the reconstruction from the command line
#
# i2 - name of the image containing the map where the coordinates have been calculated, 
#  with its path to directory if different from the current working directory
#
i2=sys.argv[1]
#
# c1 = raw_input('Write the name of the simple text file containing the coordinates of the crimes:
# the first line containing the x coordinates and the second line containing the y coordinates
# (starting from left, up, as in GIMP images; the values must be separated by spaces)
#
c1=sys.argv[2]
#
# RES_PREFIX - the prefix name of the results with probability of finding the origin of the events
#
RES_PREFIX=sys.argv[3]
final_image=RES_PREFIX+".bmp"
#
# controllo che i dati siano stati inseriti altrimenti lascio il default
#
# parameters required from the Rossmo Formula
#
if (len(sys.argv) > (NUM_ARGS+1)):
    if (sys.argv[4]<>"E"):
     Dist = sys.argv[4]
if (len(sys.argv) > (NUM_ARGS+2)):
    B=float(sys.argv[5])
if (len(sys.argv) > (NUM_ARGS+3)):
    f=float(sys.argv[6])
if (len(sys.argv) > (NUM_ARGS+4)):
    g=float(sys.argv[7])
#
print " Elaboration with parameters B=",B," f=",f," g=",g," Distance=", Dist
#
#take as input the files containing the coordinates of the crimes and open it in read mode
#
c = open(c1, 'r')
#
# Create the lists of the coordinates of crimes from a csv file
#
crimesx=[]
crimesy=[]
weight=[]
#
# read the csv file
#
with open(c1, 'rb') as fcsv:
  reader = csv.reader(fcsv)
  for row in reader:
#   num, a, b, c= row
   a, b= row
   crimesx.append(a)
   crimesy.append(b)
#   weight.append(c)
#
# ora le devo trasfromare in integer
# inoltre posso contare il numero di elementi presenti nel csv: questo metodo scambia inoltre righe con colonne
#
# test that X-coords and Y-coord are the same numbers
#
NB_CRIMES=len(crimesx)
NB_CRIMESbis = len(crimesy)
if NB_CRIMES != NB_CRIMESbis:
 print('The number of x coordinates is different from the number of y coordinates!! Please check the data')
#
#
#
crime_x=[]
crime_y=[]
pesi=[]
#
# fill the lists of crimes coordinates
#
for it in range (NB_CRIMES):
  crime_x.append(float(crimesx[it]))
  crime_y.append(float(crimesy[it]))
#  pesi.append((float(weight[it])/100))  
#
# end for
#
im2 = Image.open(i2)
WIDTH, HEIGHT = im2.size
#
# tests that X-coords and Y-coords are NOT out of the image
#
uscita = 0
i=0
#for i in range (NB_CRIMES):	
#
# tests points are not near 3 pixel to the border of the image
#
#  if (crime_x[i] > (WIDTH-3)):
#    uscita = 1
#  if (crime_y[i] > (HEIGHT-3)):
#    uscita = 1
#  if (crime_x[i] < 3):
#    uscita = 1
#  if (crime_y[i] < 3):
#    uscita = 1
#
# end for
#
#if (uscita > 0):
#
# if exists a input data error exit
#
#  print ""
#  print " one or more points are OUT of the image OR are too near the image border (less 3 pixels)- check input data - exiting!"
#  print ""
#  sys.exit (0)
#
# end if
#
# if input data are correct execute te reconstruction 
#
# qualcosa non mi torna qui - *ugo* - create an empty white image
#
i1=Image.new("RGB",(WIDTH,HEIGHT), "white")
i1.save("Imgwhite.bmp")
i1 = "Imgwhite.bmp"
#
# create the final image containing exclusively the results of elaboration
#
i4 = RES_PREFIX+"def.bmp"
im1 = Image.open(i1)
im1.save(i4)
#
# defines a zeros filled matrix with the same size of the image 
#
immagine= zeros( (WIDTH,HEIGHT) )
print ' Image dimensions in pixel:',WIDTH,HEIGHT
#------
# create a file with the same size of im1 
# it will be bodified durig the reconstruction
# and it will contains the final elaboration
#
im4 = Image.open(i4)
#
# here begins the cycles of the rossmo function
if Dist == "M" :
#---versione manhattan distance --------------------------------------------------------------------------------------
 for i in range (WIDTH):
  for j in range (HEIGHT):
   result = 0
   for n in range (NB_CRIMES):
#    manhattan distance 
     distance = (math.fabs(i - crime_x[n])+math.fabs(j-crime_y[n]))
     if distance > B:
       term1 = 1 / math.pow(distance,f)
       result = result + term1
     else:
       term2=(math.pow(B,(g-f)))/math.pow((2*B-distance),g)
       result = result + term2
#
# memorizzo il risultato in una matrice
#
   immagine[i][j]=result
# 
# it finds the max an the min value of the image
#
   if(i==0 and j==0):
     biggest_p = result;
     smallest_p = result;
   if(result > biggest_p):
     biggest_p = result; 
   if(smallest_p > result):
     smallest_p = result;
#-----End Manhattan---------------------------------------------------------------------------------------------------------
else:
#------ euclidean distance -------------------------------------------------------------------------------------------
#    teorema pitagora: (math.pow((math.fabs(i - crime_x[n])+ math.pow((math.fabs(j-crime_y[n])), 2)
#    the first term  is the x catete; the other the y catete
#---------------------------------------------------------------------------------------------------------------------
#
 for i in range (WIDTH):
  for j in range (HEIGHT):
   result = 0
   for n in range (NB_CRIMES):
#    euclidean distance 
     distance = math.sqrt(math.pow((math.fabs(i - crime_x[n])), 2) + math.pow((math.fabs(j-crime_y[n])),2))
     if distance > B:
       term1 = 1 / math.pow(distance,f)
       result = result + term1
     else:
       term2=(math.pow(B,(g-f)))/math.pow((2*B-distance),g)
       result = result + term2
    
# memorizzo il risultato in una matrice
#
   immagine[i][j]=result
# 
# it finds the max an the min value of the image
#
   if(i==0 and j==0):
     biggest_p = result;
     smallest_p = result;
   if(result > biggest_p):
     biggest_p = result; 
   if(smallest_p > result):
     smallest_p = result;
# 
#---------------------------------------------------------------------------------------------------------------------
# hence i must use biggest_p and smallest_p to standardize the value of result: smallest_p=0; biggest_P=>255 
#
# retrieves the result from the matrix and normalizes it
print " End of Calculus.. ",
print " Normalizzation and graphing.. ",
for i in range (WIDTH):
  for j in range (HEIGHT):
   result = immagine[i][j]
   resultrel = ((result-smallest_p)*255/(biggest_p-smallest_p))
#---- end for normalizzazione ------
#
# second part to insert colors that is yellow for low probability and red for high probability
#
# DA MODIFICARE in questo modo *ugo* : 
# massima probabilita' (>95%)= rosso (255,0,0)
# fra il 90 ed il 95        = arancio o rosa o qualche sfumatura del genere
# fra 85 e 90               = giallo
# fra 75 ed 80              = verde
# sotto 80                  = azzurro (sfumato)
#
   fracres, intres = math.modf(resultrel) 
#
# the above is the integer part of result, that is intres that i will use underneath
# i must force the double variable to become integer: let's see how it is
   intintres = int(intres)
#
# nuova versione di put pixel - si comincia dal valore di pixel piu basso perche piu frequente.
#
   if (intintres <= 219):
    # disegna il pixel in toni di azzurro
    im4.putpixel((i,j),(0, intintres, 255-intintres))
   elif ((intintres > 219)and (intintres <= 231)):
    # disegna il pixel in verde
    im4.putpixel((i,j),(0, intintres, 0))
   elif ((intintres>231)and(intintres <=243)):
    # disegna il pixel in giallo
    im4.putpixel((i,j),(255, 255, 0)) 
   else:
    # if (intintres > 243):
    # disegna il pixel in rosso
    im4.putpixel((i,j),(255, 0, 0))
#
#
# End FOR (doppio) -----
#
# it puts a small cross on the crimes coordinates
# rifare quest aparte di codice mettendo if crime_x[cross]< (WIDTH-3): per primo
# poi elif per passare al caso successivo: so that i controlli successivi a -2 e -1 li fa only if the first condition is not satisfied (more common
# option) may be i have to insert a return before the successive elif but probably not
# successive step to produce a function that draws a cross in that position x,y 
#
# Plots a small white circle in the point of interest
#
print" Drawing points.. ",
draw = ImageDraw.Draw(im4)
for cross in range (NB_CRIMES):
  draw.ellipse((crime_x[cross]-1, crime_y[cross]-1, crime_x[cross]+1, crime_y[cross]+1), fill=(255,255,255))
#  testo="n."+str(int(cross+1))+" ("+str(weight[cross])+"%)"
#  lung_testo=len(testo)
#  draw.text((crime_x[cross]-int(5*lung_testo/2), crime_y[cross]-13),testo,fill=(255,255,0))
#
# with the above instruction it says that 255-result is not an integer but a double and i agree
# may be i should normalize tha value but in the c instruction ai don't see anything like that
# im4.putpixel((i,j),(255, 255, 0)) if i use this instruction i obtain an image that is completely yellow
#
print" Saving image and deleting temporary files"
im4.save(i4)
im1=Image.open(i2)
im2=Image.open(i4)
im3=Image.blend(im1, im2, BLEND_PERCENTAGE)
im3.save(final_image)
#
# delete temporary file
#
os.remove(i1)
#
# non va cancellata la IM4 *ugo*
#
#uso blend the value of the mask of 0.8 apparently works well,, It can be changed
#
print" "
print" --------------------------------- "
print" "
print(' Resulting image saved with the indicated name - end of execution')
print" "
print" --------------------------------- "
#
# END PROGRAM
#
print" "