#=======================================================#
# Copyright (c) 2018. Yu-Cheng Ku. All Rights Reserved. #
# Melody of a River Town                                #
# Written by Yu-Cheng Ku                                #
#=======================================================#

rm(list=ls())

library(mvtnorm)

## Hill 
track = 
function(x)
{
	-0.05 * x^2 + 0.2
}

## Function for the track of Lady Chang-E 
track2 = 
function(x)
{
	-0.03 * x^2 + 0.2
}



## Fractal Tree here
Tree = 
function(A)
{
	B = cbind(A, A, A, A, 0*A, 0*A)
	C = cbind(A, A, 0*A, 0*A, 0*A, 0*A)
	D = cbind(A, A, A, A, A, A)
	E = rbind(C, B, D, B, C)
	return(E)
}

iter = 3 # Number of interations
for(i in 1:iter)
{
	T = matrix(1,1,1)

	# Recursion here
	for(i in 1:i)
	{
		T = Tree(T)
	}
}


## Start the graphing
win.graph(width = 6.75, height = 5.5)

par(mai = rep(0, 4), bg = "grey55")
plot(c(0,0), type="n", xlim=c(-0.5, 3.25), ylim=c(-0.48, 1.72)) 


n_scale = 25
for(i in 1:n_scale)
{
	range_u = 1.1 + i * (2.2-(1.1))/n_scale

	lines(rbind(c(-1.5,range_u),c(4.5,range_u)), col = paste("grey", 55-i*1,sep=""), lwd = 25)
}


#########
# stars #
#########
st.x.s = runif(60,-1,4)
st.y.s = runif(60,-0.5,2)

st.x.l = runif(8,-1,4)
st.y.l = runif(8,-0.5,2)

points(cbind(st.x.s, st.y.s), pch=".", col="white")
points(cbind(st.x.s, st.y.s), pch=".", col="white", cex=2)
points(cbind(st.x.l, st.y.l), pch=-as.hexmode(9733), cex=0.5, col="white")


########
# Moon #
########
points(2.4, 1.3, pch = 19, col="yellow2", cex=12)


##########
# Ground #
##########
x1 = seq(-5, 1, l=30)
x2 = seq(-9, 5, l=50)
lines(3.5+x1, track(x1)-.5,    col="#44571E", lwd=200)
lines(-1-x2,  track2(-x2)-.65, col="#4F620F", lwd=180)



###############
# Tree shadow #
###############
p = NULL

for (j in c(rep(8,8),c(9:1)))
{
      p = c(p, rnorm(1400, mean = 0, sd = sqrt(j)/2))
}

p = c(p, rnorm(400, mean = 0, sd = 0.5/2), rnorm(400, mean = 0, sd = 0.3/2), rnorm(400, mean = 0, sd = 0.1/2))

p = p/8
t = 1.9 * (1/length(p)) * c(1:length(p))

shadow = cbind(t, p) 
angle2 = -30 * pi/180

rotate = matrix(c(cos(angle2), sin(angle2), -sin(angle2), cos(angle2)), 2, 2) 

shadow = shadow %*% rotate

shadow[,1] = shadow[,1] - 1.095 
shadow[,2] = shadow[,2]/2 - 0.664

points(shadow, pch = ".", col = "#2C4B09")


############
# The Tree #
############

# tree trunk
lines(rbind(c(0.5,-0.13),  c(0.5,-0.01)), lwd = 12, col = "#722809")
points(0.495, -0.15, pch=17, cex=2.6, col="#722809")

# tree trunk shadow
lines(rbind(c(0.49,-0.14), c(0.49,-0.01)), lwd = 5, col = "#523305")
lines(rbind(c(0.49,-0.14), c(0.46,-0.17)), lwd = 7, col = "#523305")


# The Tree
image(T, col=c("grey55","#05440B"), axes=FALSE, xlim=c(-1, 4), ylim=c(-0.5, 1.9), add=TRUE, cex = 0.8)


#########
# Cloud #
#########
# Bigger
cl.x = rnorm(5000, mean = 0, sd = 2.2/9)
cl.x = cl.x + 1.8

cl.y = rnorm(5000, mean = 0, sd = 0.5/9)
cl.y = cl.y + 0.93

cloud = cbind(cl.x, cl.y)
points(cloud, pch = ".", col = "gray85", cex = 1)

# Smaller
cl.x.s = rnorm(2000, mean = 0, sd = 1/9)
cl.x.s = cl.x.s + 2.8

cl.y.s = rnorm(2000, mean = 0, sd = 0.2/9)
cl.y.s = cl.y.s + 1.3

cloud.s = cbind(cl.x.s, cl.y.s)
points(cloud.s, pch = ".", col = "gray85", cex = 1)



# Copyright
text(x = 2.9, y = -.45, label="by Yu-Cheng Ku", col = "grey60", cex = 1)