## process AFM data obtained from experiment and calibration to apply to CafeMol calculation
## Usage: Rscript process_afm_data.R afm.txt dimer_template.txt
## afm.txt is the AFM data file in Angstroms unit, dimer_template.txt is a template header file
## the output will be an input file (.inp) for cafemol

#!/usr/bin/env Rscript
args = commandArgs(trailingOnly=TRUE)
#print(args[1])
#length(args)
## test if there is at least one argument: if not, return an error
if (length(args) != 2) {
  stop("Please provide afm data file and template header file", call.=FALSE)
}

#install.packages("dplyr")
library(dplyr)

afm.file = args[1]
d = read.table(afm.file, 
               sep="\t", 
               fill=FALSE, 
               strip.white=TRUE)

# min(d$V3) 
d$V3 = d$V3 - min(d$V3)
# min(d$V3)
# dim(d)
d$id = rep("AFM_IMAGE", nrow(d))
#move last column to first
df.new <- d %>%   select(id, everything())
# dim(df.new)
# get every 5th row
#df.new = df[seq(1, nrow(df), 5), ]
df.new$V1 = round(df.new$V1*10, digits = 1)
df.new$V2 = round(df.new$V2*10, digits = 1)
## offset x and y pixel position
df.new$V1 = df.new$V1 - min(df.new$V1) + 5
df.new$V2 = df.new$V2 - min(df.new$V2) + 5
# min(df.new$V1)
# min(df.new$V2)

## set z values less than 1A to zero
V3.sparse = df.new$V3*10
V3.sparse[V3.sparse < 1.0] = 0.000

df.new$V1 = format(df.new$V1, digits = 1, nsmall = 1)
df.new$V2 = format(df.new$V2, digits = 1, nsmall = 1)
df.new$V3 = format(round(V3.sparse, digits = 3))


# head(df.new)
# dim(df.new)
## get unique combination of x y pixel position
df.new.unique = df.new[!duplicated(df.new[2:3]),]
# dim(df.new.unique)

#### Generate the input file
header.file = args[2]
top_part = readLines(header.file)

out.file = paste0(tools::file_path_sans_ext(basename(afm.file)),"_cafemol.inp")
writeLines(top_part, out.file)
write.table(df.new.unique, out.file, 
            sep = " ", quote = FALSE, row.names = FALSE, col.names = FALSE, append = TRUE)
write(">>>>", out.file, append = TRUE)

#### copy to biowulf
# scp -pr /Users/bhandariyr/Documents/projects/afm/cafemol/ring/cbl_afm_fit.inp bhandariyr@biowulf.nih.gov:/gpfs/gsfs12/users/bhandariyr/cafemol/cafemol/example/cbl_afm_fit/

#### Now make x y pixel position column divisible by n
# div_4 <- function(x) {
#   for (i in  1:length(x)) {
#     if (x[i] %% 4 == 0){  }
#     else if((x[i]+1) %% 4 == 0){ x[i] = x[i]+1 }
#     else if((x[i]-1) %% 4 == 0){ x[i] = x[i]-1 }
#     else if((x[i]+2) %% 4 == 0){ x[i] = x[i]+2 }
#     else { print(x[i])
#       stop("division error occured") }
#   }
#   return(x)
# }
# div_2 <- function(x) {
#   for (i in  1:length(x)) {
#     if (x[i] %% 2 == 0){  }
#     else if((x[i]+1) %% 2 == 0){ x[i] = x[i]+1 }
#     else { print(x[i])
#       stop("division error occured") }
#   }
#   return(x)
# }
# 
# # df.new$V1 = div_2(df.new$V1)
# # df.new$V2 = div_2(df.new$V2)
# 
