Mercurial > repos > bitlab > plidflow
comparison PLIDflow/scripts/pdbcenter_npts_finder.R @ 2:afd5b5ffc38f draft
Uploaded
author | bitlab |
---|---|
date | Tue, 14 Jan 2020 07:52:48 -0500 |
parents | 6fcfa4756040 |
children |
comparison
equal
deleted
inserted
replaced
1:eda62adfc858 | 2:afd5b5ffc38f |
---|---|
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 |