-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathstep.4.pairing.GSPs.R
61 lines (52 loc) · 1.74 KB
/
step.4.pairing.GSPs.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
###########################
## pair GSP1-GPS2
## by pair-panelty
###########################
currdir = getwd()
system('rm -rf pairing', ignore.stderr =T)
system('mkdir pairing')
setwd('pairing')
## generate list of targets from primer.psl
system("sort -k1,1 ../primer.psl > primer.psl.srt")
system("cut -d ' ' -f 1 primer.psl.srt > _targets")
system("sort -u _targets > target.list")
system("rm _targets")
## For multi-threading
pair.targets = readLines('target.list')
n.targets = length(pair.targets)
n.tar.cpu = ceiling(n.targets/ncpu)
seq1 = seq(1, n.targets, n.tar.cpu)
seq2 = seq(n.tar.cpu, n.targets, n.tar.cpu)
if (length(seq2) != length(seq1)){
seq2 = c(seq2, n.targets)
}
n.seq = length(seq1)
# prep sub targets
for (cpu in 1:n.seq){
system(paste('mkdir ', cpu, sep=''))
target.n1 = seq1[cpu]
target.n2 = seq2[cpu]
targets.cpu = pair.targets[target.n1 : target.n2]
writeLines(targets.cpu, paste(cpu, '/targets', sep=''))
}
# multithread running
for (cpu in 1:n.seq){
system(paste('touch _running_', cpu, sep=''))
system(paste('Rscript --vanilla ', ampdir, '/scripts/pairing.GSPs.sub.R '
, cpu, ' ', ampdir, ' &', sep=''))
}
#=== wait
Sys.sleep(5)
nrun=length(list.files('./', '_running_'))
while (nrun >0){
Sys.sleep(5)
nrun=length(list.files('./', '_running_'))
}
##=== all cpus paired
system('cat */missed.targets.cpu > ../unpaired.targets')
system('cat */paired.ranked.cpu > paired.ranked.allcpus')
system('head -n1 paired.ranked.allcpus > ../ranked.all.candidate.pairs')
system('grep -v penalty paired.ranked.allcpus >> ../ranked.all.candidate.pairs')
setwd(currdir)
# remove intermediate data
system('rm -rf pairing')