## Makefile ---
P=Road name numbers of Louth area of St. Cathatines and Lincoln, Niagara District, Ontario, Canada
## These are from an underlying Louth Township Lot ... Concession ...
## Author: Dan Jacobson https://www.jidanni.org/
## Copyright: https://www.gnu.org/licenses/gpl.html
## Created: 2024-01-12T12:35:03+0000
## Last-Updated: 2024-01-14T23:30:07+0000
##     Update #: 103

# I made a great KMZ map. I'm real proud of it. That the underlying
# grid is a little crooked... you can't blame me. However what I'm not
# proud of is the methods I used. The axes and positive negative stuff
# and offsets got my head spinning. In the end I just tried many
# combinations until I got it right!

louth.vik:
#Intersections. (Avenues run east-west, streets north-south.)
S01A01=43.183015 -79.277227 #Um, projected.
S01A08=43.11805  -79.274971 #Underlying this is Louth Township Lot 1 Concession 8.
S03A03=43.164099 -79.28721
S21A08=43.114974 -79.37716
G=geod +ellps=WGS84 -I -p -f %.6f
#central_line_azimuth=$$(echo $(S21A08) $(S01A08)| $G | perl -anwle 'print $$F[0]') #too crooked, use instead:
central_line_azimuth=$$(echo $(S01A01) $(S01A08)| $G | perl -anwle 'print $$F[0] - 90')
streets_per_meter=$$(echo $(S01A08) $(S21A08) | $G | perl -anwle 'print 0 + ( 21 - 1 ) / $$F[-1]')
avenues_per_meter=$$(echo $(S01A01) $(S01A08) | $G | perl -anwle 'print 0 + ( 8 - 1 ) / $$F[-1]')
s0=+proj=omerc +alpha=$(central_line_azimuth) +gamma=0 +lat_0=$(firstword $(S01A01)) +lonc=$(lastword $(S01A01))
s1=+proj=affine +s11=$(avenues_per_meter) +s22=$(streets_per_meter)
s2=+proj=affine +xoff=1 +yoff=-1
p=+proj=pipeline +step $(s0) +step $(s1) +step $(s2)
louth.csv:louth0.csv #check.csv
	proj -r -f %.6f -I $p < $< | perl -pwle 's/\s+/,/g;' > $@
#	cat check.csv >> $@
#OK, what if we want to not use the pipeline, and just use one +proj string with no +steps?
#Well, yes, but then we would have to do each axis separately, as they would have different +k's.
louth0.csv: Makefile
	perl -wle 'BEGIN { $$, = " " }; for ( 1 .. 21 ) { print -$$_, 8, $$_; } $(\
	)for ( 1 .. 8 ) { print -1, $$_, $$_; }' > $@
N=_POSSIBLE_NAMES=field_
KO = 	-f LIBKML -dsco NAME="$P" -dsco DESCRIPTION="See $(subst $(HOME)/,https://www.,$(PWD))"
CO=	-oo HEADERS=NO -oo KEEP_GEOM_COLUMNS=NO -oo X$N1 -oo Y$N2 \
	-sql 'SELECT field_3 AS Name FROM "$(basename $<)"'
#Yes, we would like to make a single declaration, but as mentioned above, the two axes have different +k's.
info:; @gdalsrsinfo -V -o all "$$(echo $(s0))" #so we just explore our basic transformation.
check.csv:Makefile
	perl -nwle 'BEGIN{$$,=","}next unless /(^S..A..)=(\S+)\s+(\S+)/; print reverse @{^CAPTURE};' $< > $@
v=viking
%.kmz %.kml:%.csv; ogr2ogr $@ $? $(KO) $(CO)
%.vik:%.kmz
	if pidof -q  $v; then killall $v; else :; fi
	cd /tmp && nohup $v $(PWD)/$? &
	sleep 19
.PRECIOUS: %.kmz