# General structure of all figure codes:

# Code blocks are usually named after neuron types
# neuron skeleton identifier (skid) is read from .csv files that are generated by "conn_table.R"
# With skid, one need CATMAID access (https://fafb.catmaid.virtualflybrain.org/) to obtain the neuron skeleton and connectivity info.
# Usually, each neuron type has a top view, a side view, and some codes for generating a synapse depth table.
# The synapse depth table is used to generate a histogram
# Saving to .png files is commented out by default. 


# Here we use Dm9 as an example. 
# Run "startup.R" and "layer_centraleye.R" first

# Dm9 ---------------------------------------------------------------------

# scale bar
sbar <- matrix(c(40000,0,5000, 40000,0,15000), ncol=3,byrow=T) %*% t(xRot)

type <- "Dm9"
# read in the neuron info from csv, which is generated by "conn_table.R"
tb <- read.csv(paste("table_by_type_out/", "R7R8_outgoing_", type, ".csv", sep = ''))

skid <- tb$skid %>% na.omit() #neuron skeleton IDs
neu <-  read.neurons.catmaid(skid, .progress='text' ) #read neuron from CATMAID server
neu <- fix_jump(neu) # fix simple jumps in the data, mostly for cosmetics
neu_xform <- xEucl_neu(neu, me_pca$rotation, me_pca$center) # transform to the desired ref. frame
neu_xform <- nlapply(neu_xform, function(x) smooth_neuron(x, sigma = 5e2))


ii <- which(tb$name %in% eg_name) # look for the complete neuron(s) to plot

# - top view
nopen3d()
par3d('windowRect' = c(10, 10, 2110, 2110)) #window size
plot3d(neu_xform[[ii]], col='black', lwd=2, soma = T, WithNodes = F, lit=F) #draw neurons
pch3d(xyz_M5_avg_xform_yz[Mi1_ind_Nnb_p,], radius=5000,col="#9ecae1",pch=16,alpha=0.3) #seed column shading
pch3d(xyz_M5_avg_xform_yz[Mi1_ind_Nnb_y,], radius=5000,col="#fec44f",pch=16,alpha=0.3)
pch3d(xyz_M5_avg_xform_yz[Mi1_ind_p,], radius=5000,col='#9ecae1',pch=16,alpha=0.6)
pch3d(xyz_M5_avg_xform_yz[Mi1_ind_y,], radius=5000,col="#fec44f",pch=16,alpha=0.6)
segments3d(sweep(sbar,2,c(0,12000,17000), '+'), lwd=2) #scale bar
rgl.viewpoint(fov=0,zoom=0.55, userMatrix= rotationMatrix(-90/180*pi,0,1,0) %*% rotationMatrix((90+xAng)/180*pi,1,0,0)) #adjusting viewing angle

# load connectivity info (no. and [x,y,z]) from CATMAID
conn_in_R7 <- catmaid_get_connectors_between(pre_skids = c(skid_pR7,skid_yR7), post_skids = skid[ii])
if (!is.null(conn_in_R7)) {
  conn_in_R7 <- conn_in_R7[, c("post_node_x", "post_node_y", "post_node_z")]
  conn_in_R7_xform <- sweep(as.matrix(conn_in_R7), 2, me_pca$center) %*% me_pca$rotation
  pch3d(conn_in_R7_xform, pch=16, radius=1000, alpha=0.7, col=pal_syn['R7'])
}
conn_in_R8 <- catmaid_get_connectors_between(pre_skids = c(skid_pR8, skid_yR8), post_skids = skid[ii])
if (!is.null(conn_in_R8)) {
  conn_in_R8 <- conn_in_R8[, c("post_node_x", "post_node_y", "post_node_z")]
  conn_in_R8_xform <- sweep(as.matrix(conn_in_R8), 2, me_pca$center) %*% me_pca$rotation
  pch3d(conn_in_R8_xform, pch=16, radius=1000, alpha=0.8, col=pal_syn['R8'])
}

# saving to a png file
# rgl.snapshot(filename = paste("F2", type, skid[ii], "top_zoom.png", sep = '_'))


# - side view
nopen3d()
par3d('windowRect' = c(10,10, 1310, 2110))
plot3d(neu_xform[[ii]], col='black', lwd=2, soma = T, lit=F)
plot3d(pR8_xform[[1]], lwd=3, col=pal_syn['R8'])
shade3d(cp1, col='#9ecae1', alpha = .1, lit=F)
shade3d(cp2, col='#9ecae1', alpha = .1, lit=F)
shade3d(cy1, col='#fec44f', alpha = .1, lit=F)
shade3d(cy2, col='#fec44f', alpha = .1, lit=F)
segments3d(layers_side, lwd=1)
# text3d(layer_M_anno, texts = paste("M", seq(1,10), sep = ''), adj = 0)
segments3d(sweep(sbar,2,c(0,0,-10000), '+'), lwd=2)
# text3d(colMeans(sbar)+c(2000,0,0), texts = "10 um", cex = 1.5)
rgl.viewpoint(fov=0,zoom=0.4, userMatrix= rotationMatrix(90/180*pi,0,0,1) %*% rotationMatrix(xAng/180*pi,1,0,0))
# legend3d(x=0.7, y=0.7, c("R7", "R8"), pch = c(16, 16), col = pal_syn)

# separate R7 and R8 input
conn_in_R7 <- catmaid_get_connectors_between(pre_skids = c(skid_pR7,skid_yR7), post_skids = skid[ii])
if (!is.null(conn_in_R7)) {
  conn_in_R7 <- conn_in_R7[, c("post_node_x", "post_node_y", "post_node_z")]
  conn_in_R7_xform <- sweep(as.matrix(conn_in_R7), 2, me_pca$center) %*% me_pca$rotation
  pch3d(conn_in_R7_xform, pch=16, radius=1000, alpha=0.7, col=pal_syn['R7'])
}
conn_in_R8 <- catmaid_get_connectors_between(pre_skids = c(skid_pR8, skid_yR8), post_skids = skid[ii])
if (!is.null(conn_in_R8)) {
  conn_in_R8 <- conn_in_R8[, c("post_node_x", "post_node_y", "post_node_z")]
  conn_in_R8_xform <- sweep(as.matrix(conn_in_R8), 2, me_pca$center) %*% me_pca$rotation
  pch3d(conn_in_R8_xform, pch=16, radius=1000, alpha=0.8, col=pal_syn['R8'])
}

# rgl.snapshot(filename = paste("F2", type, skid[ii],  "side_incoming.png", sep = '_') )


# separate R7 and R8 output
conn_out_R7 <- catmaid_get_connectors_between(pre_skids = skid[ii], post_skids = c(skid_pR7,skid_yR7))
if (!is.null(conn_out_R7)) {
  conn_out_R7 <- conn_out_R7[, c("post_node_x", "post_node_y", "post_node_z")]
  conn_out_R7_xform <- sweep(as.matrix(conn_out_R7), 2, me_pca$center) %*% me_pca$rotation
  pch3d(conn_out_R7_xform, pch=16, radius=1000, alpha=0.7, col=pal_syn['R7'])
}
conn_out_R8 <- catmaid_get_connectors_between(pre_skids = skid[ii], post_skids = c(skid_pR8, skid_yR8))
if (!is.null(conn_out_R8)) {
  conn_out_R8 <- conn_out_R8[, c("post_node_x", "post_node_y", "post_node_z")]
  conn_out_R8_xform <- sweep(as.matrix(conn_out_R8), 2, me_pca$center) %*% me_pca$rotation
  pch3d(conn_out_R8_xform, pch=16, radius=1000, alpha=0.7, col=pal_syn['R8'])
}

# SAVE
# rgl.snapshot(filename = paste("F2", type, skid[ii],  "side_outgoing.png", sep = '_') )
 

# # - connection depth table
# tb <- mktb_in_py78(skid)
# write.csv(tb, "F2_Dm9_all_incoming.csv", row.names=F)
# tb <- mktb_out_py78(skid)
# write.csv(tb, "F2_Dm9_all_outgoing.csv", row.names=F)
# 
# 
# for (ii in 1:length(skid)) {
#   tb <- mktb_in_py78(skid[ii])
#   write.csv(tb, paste("F2_Dm9_", skid[ii], "_incoming.csv",sep = ''), row.names=F)
#   tb <- mktb_out_py78(skid[ii])
#   write.csv(tb, paste("F2_Dm9_", skid[ii], "_outgoing.csv",sep = ''), row.names=F)
# }