############################
# ggplot method  center
############################
gmap = get_map(location=c(-95.4142, 29.7907), source="google",zoom=13)
ggmap(gmap)
master = readRDS("/home/ajackson/mirrors/ajackson/crime/data/masterdataset.rds")

# join lat longs to master table and then delete extra columns
mapdata = left_join(master, geotable, by="Address")
mapdata$Add1 <- NULL
mapdata$Add2 <- NULL
mapdata$Street.y <- NULL
mapdata$lon1 <- NULL
mapdata$lat1 <- NULL
mapdata$lon2 <- NULL
mapdata$lat2 <- NULL
mapdata$ShortPremise = str_replace_all(mapdata$ShortPremise,"^ ","")
mapdata$Offense_Type = str_replace_all(mapdata$Offense_Type," ","_")
# filter data to what is needed for a plot
plotpath="/home/ajackson/Dropbox/CrimeStats/CrimeMaps"
premiselist = c("Parking","Business","Residence","Street")
offenselist = c("Theft","Burglary","Robbery","Auto_Theft")

# Filter out stuff not needed
offensefilter = "Auto_Theft"
premisefilter = "Parking"

for (offensefilter in offenselist)
  {
  for (premisefilter in premiselist)
  {
  pointdata = mapdata %>%
  filter(!is.na(avglat)) %>%  # drop rows without a lat long
  filter(Offense_Type == offensefilter) %>% # pick out a particular offense
  filter(as.character(ShortPremise) == premisefilter) %>% #  pick out a premise
  filter(avglat>29.78 & avglat<29.82 & avglon> -95.415 & avglon< -95.3880) # eliminate wrong neighborhood records

if (nrow(pointdata) == 0) {
  print(paste("###### break", premisefilter, offensefilter)) 
  break
  }

#   summarize by location and find peak value
pointdatasum <- pointdata %>%
  group_by(Address) %>%
  summarise(total=n())
peak = as.numeric(max(pointdatasum$total))

ggmap(gmap, extent='normal', maprange=FALSE) %+% pointdata + aes(x = avglon, y = avglat) +
  geom_density2d() + aes_(peak=peak) +
  stat_density2d(aes(fill = ..level../max(..level..)*peak, alpha = ..level..), size=.01, geom='polygon', n=50, contour = TRUE , show.legend = FALSE) +
  scale_fill_gradient(low = "green", high = "red",name=offensefilter ) +
  scale_alpha(range = c(0.00, 0.35), guide = FALSE) +
  labs(title=paste(offensefilter,"in",premisefilter)) +
  coord_map(projection="mercator", xlim=c(-95.42, -95.38), ylim=c(29.775, 29.815)) +
  theme(legend.position = "right", axis.title = element_blank(), text = element_text(size = 12))

ggsave(paste(offensefilter,"_",premisefilter,".jpeg", sep=""), device="jpeg", path = plotpath)
  }
}
for (offensefilter in offenselist)
  {
  for (premisefilter in premiselist)
  {
  pointdata = mapdata %>%
  filter(!is.na(avglat)) %>%  # drop rows without a lat long
  filter(Offense_Type == offensefilter) %>% # pick out a particular offense
  filter(as.character(ShortPremise) == premisefilter) %>% #  pick out a premise
  filter(avglat>29.78 & avglat<29.82 & avglon> -95.415 & avglon< -95.3880) # eliminate wrong neighborhood records

if (nrow(pointdata) == 0) {
  print(paste("###### break", premisefilter, offensefilter)) 
  break
}
  
pointdatasum <- pointdata %>%
  group_by(Address) %>%
  summarise(total=n())
pointdatasum = left_join(pointdatasum,geotable, by="Address")
pointdatasum <- pointdatasum[complete.cases(pointdatasum),]
pointdatasum = distinct(pointdatasum)

ggmap(gmap, extent='normal', maprange=FALSE, show.legend=FALSE) %+% pointdatasum + aes(x = avglon, y = avglat) +
  geom_point(data=pointdatasum, aes(size=total,color=total)) +
  scale_color_continuous(guide=FALSE) +
  scale_alpha(range = c(0.00, 0.35), guide = FALSE) +
  labs(title=paste(offensefilter,"in",premisefilter)) +
  coord_map(projection="mercator", xlim=c(-95.42, -95.38), ylim=c(29.775, 29.815)) +
  theme(legend.position = "right", axis.title = element_blank(), text = element_text(size = 12))

ggsave(paste("bubble_",offensefilter,"_",premisefilter,".jpeg", sep=""), device="jpeg", path = plotpath)
  }
}