Mercurial > repos > bitlab > plidflow

comparison PLIDflow/scripts/pdbcenter_npts_finder.R @ 6:795e11fac81b draft default tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Included new tools for standardization
author bitlab
date Wed, 22 Apr 2020 06:12:00 -0400
parents 6fcfa4756040
children
comparison
equal deleted inserted replaced
5:97f12f7cc852 6:795e11fac81b
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