|
0
|
1 #pdbcenter_npts_finder.R calculates the geometric center of a PDB file
|
|
|
2
|
|
|
3 #!/usr/bin/env Rscript
|
|
|
4 args = commandArgs(trailingOnly=TRUE)
|
|
|
5
|
|
|
6 if(length(args) < 1){
|
|
|
7 stop("USE: Rscript pdbcenter_npts_finder.R <receptor_name> <session_dir>")
|
|
|
8 }
|
|
|
9
|
|
|
10 #Arguments definition
|
|
|
11 receptor_name <- args[1]
|
|
|
12 session_dir <- args[2]
|
|
|
13
|
|
|
14 #Read PDB file
|
|
|
15 archivo <- scan(args[1], what = character(), quiet = TRUE)
|
|
|
16
|
|
|
17 #Make a table where positions of carbon (C) atoms and their x, y, z coordenates will be placed
|
|
|
18 salto <- 0
|
|
|
19 for(i in 1:(length(archivo)-2)){
|
|
|
20 if(archivo[i] == "ATOM" && archivo[i+2] == "C"){
|
|
|
21 salto <- salto + 1
|
|
|
22 }
|
|
|
23 }
|
|
|
24 tablaC <- matrix(1, nrow = salto, ncol = 4)
|
|
|
25
|
|
|
26 #Header. Note: first column: C atom position (will be filled in when the number of C will be known)
|
|
|
27 #Row labels
|
|
|
28 colnames(tablaC) <- c("C Position","X", "Y", "z")
|
|
|
29 rownames(tablaC) <- 1:salto
|
|
|
30
|
|
|
31 #Fill in tablaC with x, y,z coordenates of the C atoms as numeric
|
|
|
32 salto <- 0
|
|
|
33 posicionC <- c()
|
|
|
34 for(i in 1:(length(archivo)-2)){
|
|
|
35 if(archivo[i] == "ATOM" && archivo[i+2] == "C"){
|
|
|
36 salto <- salto + 1
|
|
|
37 posicionC <- c(posicionC, archivo[i+1])
|
|
|
38
|
|
|
39 tablaC[salto,2] <- as.numeric(archivo[i+6])
|
|
|
40 tablaC[salto,3] <- as.numeric(archivo[i+7])
|
|
|
41 tablaC[salto,4] <- as.numeric(archivo[i+8])
|
|
|
42 }
|
|
|
43 }
|
|
|
44
|
|
|
45 #Add C atom positions in the column 1
|
|
|
46 for(i in 1:salto){
|
|
|
47 tablaC[i,1] <- as.numeric(posicionC[i])
|
|
|
48 }
|
|
|
49
|
|
|
50 #Geometric center calculation as grid center x, y, z to be introduced in the GPF file
|
|
|
51 pto_x_medio <- (max(tablaC[,2]) + min(tablaC[,2]))/2
|
|
|
52 pto_y_medio <- (max(tablaC[,3]) + min(tablaC[,3]))/2
|
|
|
53 pto_z_medio <- (max(tablaC[,4]) + min(tablaC[,4]))/2
|
|
|
54
|
|
|
55 pto_x_medio_r <- round(pto_x_medio,digits=3)
|
|
|
56 pto_y_medio_r <- round(pto_y_medio,digits=3)
|
|
|
57 pto_z_medio_r <- round(pto_z_medio,digits=3)
|
|
|
58
|
|
|
59 #print(paste(pto_x_medio_r, pto_y_medio_r, pto_z_medio_r, sep = ";"))
|
|
|
60
|
|
|
61 #Calculate distances in each coordenate to have de number of points (npts) x, y, z to be introduced in the GPF file
|
|
|
62 eu_x <- max(tablaC[,2])-min(tablaC[,2])
|
|
|
63 eu_y <- max(tablaC[,3])-min(tablaC[,3])
|
|
|
64 eu_z <- max(tablaC[,4])-min(tablaC[,4])
|
|
|
65
|
|
|
66 eu_x_r <- round(eu_x,digits=0)
|
|
|
67 eu_y_r <- round(eu_y,digits=0)
|
|
|
68 eu_z_r <- round(eu_z,digits=0)
|
|
|
69
|
|
|
70 #print(paste(eu_x, eu_y, eu_z, sep = ";"))
|
|
|
71 #print(paste(eu_x_r, eu_y_r, eu_z_r, sep = ";"))
|
|
|
72
|
|
|
73 if(eu_x > 126 || eu_y > 126 || eu_x > 126){
|
|
|
74 print("At least one dimension in number of points is out of range (>126) according AutoDock software")
|
|
|
75 break
|
|
|
76 }
|
|
|
77
|
|
|
78 #Write PDB geometric center and eu_x_r, eu_y_r, eu_z_r in a file
|
|
|
79 #write(paste(pto_x_medio_r, pto_y_medio_r, pto_z_medio_r, eu_x_r, eu_y_r, eu_z_r,sep= "\n"), file = "/home/galaxy/galaxy/tools/proteindocking/scripts/pdbcenter_npts.txt")
|
|
|
80 write(paste(pto_x_medio_r, pto_y_medio_r, pto_z_medio_r, eu_x_r, eu_y_r, eu_z_r,sep= "\n"), file = "pdbcenter_npts.txt")
|
|
|
81
|
|
|
82
|
|
|
83
|
|
|
84
|
