comparison PLIDflow/scripts/pdbcenter_npts_finder.R @ 0:6fcfa4756040 draft

Uploaded
author bitlab
date Tue, 14 Jan 2020 06:09:42 -0500
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:6fcfa4756040
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