Nothing
"mktrj.enma" <- function(enma=NULL, # enma data structure
pdbs=NULL, # pdbs object
s.inds=NULL, # structure ids
m.inds=NULL, # modes ids
mag=10, # magnification factor
step=1.25, # step size
file=NULL, # output pdb file
rock=TRUE, ncore=NULL,
... ) { # args for write.pdb
## make a trjactory of atomic displacments along a given mode
if(!inherits(enma, "enma"))
stop("mktrj.enma: must supply 'enma' object, i.e. from 'nma.pdbs'")
if(!is.null(pdbs)) {
if(!inherits(pdbs, "pdbs"))
stop("input 'pdbs' should be a list object as obtained from 'read.fasta.pdb'")
}
## Parallelized by parallel package
ncore <- setup.ncore(ncore, bigmem = FALSE)
if(ncore>1)
mylapply <- mclapply
else
mylapply <- lapply
if(is.null(s.inds))
s.inds <- 1:nrow(enma$fluctuations)
if(is.null(m.inds))
m.inds <- 1:5
if(is.null(file) & length(s.inds)==1 & length(m.inds)==1)
file <- paste("mode_", m.inds+6, "-s", s.inds, ".pdb", sep="")
## check dimensions to determine rm.gaps
if(any(is.gap(enma$xyz))
& (dim(enma$U.subspace)[1] == dim(enma$xyz)[2])) {
rm.gaps <- FALSE
}
else {
rm.gaps <- TRUE
}
if(!rm.gaps & length(s.inds)>1)
stop(paste("enma object must be calculated with argument rm.gaps=TRUE", "\n",
"for trajectory generation of multiple structures and modes"))
if(any(enma$L[s.inds, m.inds]<=0))
warning("Mode with eigenvalue <=0 detected. Check 'mode' index.")
nstep <- c(seq(step, to=mag, by=step))
zcoor <- cbind(1) %*% nstep
scor <- function(x,u,m) { return(x*u+m) }
myMktrj <- function(i) {
coor <- NULL
ind <- s.inds[i]
for(j in 1:length(m.inds)) {
mode <- m.inds[j]
u.inds <- which(!is.na(enma$U.subspace[,mode,ind]))
if(rm.gaps)
xyz.inds <- gap.inspect(enma$xyz)$f.inds
else
xyz.inds <- u.inds
plus <- sapply(c(zcoor), scor, u=enma$U.subspace[u.inds,mode,ind], m=enma$xyz[ind,xyz.inds])
minus <- sapply(c(-zcoor), scor, u=enma$U.subspace[u.inds,mode,ind], m=enma$xyz[ind,xyz.inds])
if(rock) {
tmp <- cbind(enma$xyz[ind,xyz.inds],
plus, plus[,rev(1:ncol(plus))],
enma$xyz[ind,xyz.inds],
minus, minus[,rev(1:ncol(minus))])
}
else {
tmp <- cbind(plus[,rev(1:ncol(plus))],
enma$xyz[ind,xyz.inds],
minus)
}
coor <- rbind(coor, t(tmp))
}
return(coor)
}
## do the calc
coor <- mylapply(1:length(s.inds), myMktrj)
coor <- do.call(rbind, coor)
class(coor) <- "xyz"
if(!is.null(file)) {
if(rm.gaps)
xyz.inds <- gap.inspect(enma$xyz)$f.inds
else
xyz.inds <- which(!is.na(enma$U.subspace[,m.inds[1],s.inds[1]]))
if(is.null(pdbs))
write.pdb(xyz=coor, file=file, ...)
else {
write.pdb(xyz=coor, file=file,
chain=pdbs$chain[s.inds[1], xyz2atom(xyz.inds)],
resno=pdbs$resno[s.inds[1], xyz2atom(xyz.inds)],
resid=pdbs$resid[s.inds[1], xyz2atom(xyz.inds)],
b=enma$fluctuations[s.inds[1], !is.gap(enma$fluctuations[s.inds[1],])],
...)
}
invisible(coor)
}
else {
return(coor)
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.