ee1ed09f9f5287d3ee64da8f69bf65afdfb4ab9d
[django_unres.git] / django_simple / todo / jobfiles.py
1 import logging
2
3 from django.db.models.signals import post_save, post_delete
4 from django.dispatch import receiver
5 from .models import Task
6 import os
7 import shutil
8 import subprocess
9 import json
10 from string import Template
11 import textwrap
12
13 logging.basicConfig()
14 logger = logging.getLogger(__name__)
15
16 def seq_2d_xx(line,seq):
17    line=line.replace(" ","")
18    lista=list(line)
19    for ii in range(0,len(line)):
20     if seq[ii]=='X':
21       if ii>0:
22             lista[ii-1]='-'
23       if ii+1<len(line):
24             lista[ii+1]='-'
25    return("".join(lista)) 
26
27 def seq_2d_write(line):
28    line=line.replace(" ", "")
29    out='{:d}\n{:6.1f}\n'.format(len(line),50.0)
30    for i in range(1,4):
31      out=out+'{:5d}{:8.3f}{:8.3f}\n'.format(i,0.0,180.0)
32    for ii in range(1,len(line)-2):
33      i+=1
34      if line[ii:ii+2]=='HH':
35        out=out+'{:5d}{:8.3f}{:8.3f}\n'.format(i,50.0,20.0)
36      elif line[ii:ii+2]=='EE':
37        out=out+'{:5d}{:8.3f}{:8.3f}\n'.format(i,180.0,40.0)
38      else:
39        out=out+'{:5d}{:8.3f}{:8.3f}\n'.format(i,0.0,180.0)
40    return(out)
41
42 def seq_2d_write_50(line):
43    line=line.replace(" ", "")
44    out='{:d}\n'.format(len(line))
45    for i in range(1,4):
46      out=out+'{:5d}{:8.3f}{:8.3f}{:8.3f}\n'.format(i,0.0,180.0,50.0)
47    for ii in range(1,len(line)-2):
48      i+=1
49      if line[ii:ii+2]=='HH':
50        out=out+'{:5d}{:8.3f}{:8.3f}{:8.3f}\n'.format(i,50.0,20.0,50.0)
51      elif line[ii:ii+2]=='EE':
52        out=out+'{:5d}{:8.3f}{:8.3f}{:8.3f}\n'.format(i,180.0,40.0,50.0)
53      else:
54        out=out+'{:5d}{:8.3f}{:8.3f}{:8.3f}\n'.format(i,0.0,180.0,50.0)     
55    return(out)
56
57 def write_ssbond(ssbond):
58    list=[]
59    if ssbond!='':
60      ssbond=json.loads(ssbond)
61      for e in ssbond:
62            list.append(e[0])
63            list.append(e[1])
64      list_sorted=sorted(list)
65      
66    if len(list)>0:
67        return str(len(list))+' '+' '.join(map(str,list_sorted))+'\n'+\
68             str(len(ssbond))+' '+' '.join(map(str,list))+'\n'
69    else: 
70        return '0\n'
71
72
73 @receiver(post_save, sender=Task)
74 def write_on_task_save(sender, instance, **kwargs):
75    if instance.ready and not instance.done and instance.running==0:
76      if (not os.path.isdir(instance.jobdirname)):
77         os.mkdir(instance.jobdirname)
78 #     os.chdir(instance.jobdirname)
79  
80
81
82 #md=Template('nstep=$nstep ntwe=100  ntwx=100   dt=0.20 damax=10.0 lang=1'+
83 #   't_bath=300 reset_vel=0 respa ntime_split=1 maxtime_split=512 mdpdb')
84 #wrapper = textwrap.TextWrapper(width=70)
85 #word_list = wrapper.wrap(text=md.substitute(nstep=str(nstep)))
86
87      if  instance.unres_ff =='FF2':
88       w="""WSC=1.00000 WSCP=1.23315 WELEC=0.84476 WBOND=1.00000 WANG=0.62954              &
89 WSCLOC=0.10554 WTOR=1.34316 WTORD=1.26571 WCORRH=0.19212 WCORR5=0.00000        &
90 WCORR6=0.00000 WEL_LOC=0.37357 WTURN3=1.40323 WTURN4=0.64673 WTURN6=0.00000    &
91 WVDWPP=0.23173 WHPB=1.00000 WSCCOR=0.25                                        &
92 CUTOFF=7.00000 WCORR4=0.00000"""
93      else:
94       w="""WSC=0.82686 WSCP=0.96947 WELEC=0.79373 WBOND=1.00000 WANG=0.46542              &
95 WSCLOC=0.07969 WTOR=0.81684 WTORD=0.67806 WCORRH=0.00000 WCORR5=0.00000        &
96 WCORR6=0.00000 WEL_LOC=0.71100 WTURN3=2.30298 WTURN4=0.86517 WTURN6=0.00000    &
97 WSCCOR=0.14577 WVDWPP=0.17781 WHPB=1.00000 WSCP14=0.00000                      &
98 CUTOFF=7.00000 WCORR4=0.00000"""
99      
100
101      with open(instance.jobdirname+'/file.inp','w') as f:
102        f.write(instance.name + ' UNRES server job'+ '\n')
103
104        if instance.type == 'min' :
105
106           control_line = 'SEED='+str(instance.md_seed)+' ' \
107              + instance.min_choice \
108              + ' pdbstart pdbref ' \
109              + ' BOXX='+str(instance.boxx) \
110              + ' BOXY='+str(instance.boxy) \
111              + ' BOXZ='+str(instance.boxz) +' '
112           if not instance.min_overlap:
113             control_line = control_line+'overlap '
114           if not instance.min_searchsc:
115             control_line = control_line+'nosearchsc '
116           if instance.min_pdbout:
117             control_line = control_line+'pdbout '
118           if instance.min_unres_pdb:
119             control_line = control_line+'unres_pdb '
120             
121           type_line = 'print_min_ini print_min_res print_min_stat'\
122             +' MAXMIN='+str(instance.min_maxmin)\
123             +' MAXFUN='+str(instance.min_maxfun)
124            
125
126        seq=instance.md_seq.replace(" ","")
127        dimen3=len(seq)
128        if dimen3>1:
129           if seq[0]=='X':
130             dimen3-=1
131           if seq[-1]=='X':
132             dimen3-=1
133   
134        nogly ='CDSQKIPTFNHLRWAVEYMcdsqkiptfnhlrwaveym'
135        seq2 = ''.join([c for c in seq if c in nogly])
136   
137        dimen3=(dimen3+len(seq2))*3
138          
139             
140        if instance.type == 'md' :
141           control_line = 'SEED='+str(instance.md_seed)+' '\
142            +instance.md_start +' md one_letter'\
143            + ' BOXX='+str(instance.boxx) \
144            + ' BOXY='+str(instance.boxy) \
145            + ' BOXZ='+str(instance.boxz) +' '
146
147           if instance.md_pdbref:
148              control_line = control_line+'pdbref '
149             
150           type_line = '        reset_vel=0'\
151              +' t_bath='+str(instance.md_temp)\
152              +' nstep='+str(instance.md_nstep)\
153              +' ntwe='+str(instance.md_ntwe)\
154              +' ntwx='+str(instance.md_ntwx)\
155              +' dt='+str(instance.md_dt)
156           if instance.md_lang == 'langevin':
157             type_line = type_line + ' lang=1 '
158             type_line = type_line + 'scal_fric='+str(instance.md_scal_fric)
159           if instance.md_lang == 'berendsen':
160             type_line = type_line + ' lang=0 tbf '
161             type_line = type_line + 'tau_bath='+str(instance.md_tau)
162           if instance.md_respa:
163             type_line = type_line + ' respa '
164           if instance.md_mdpdb:
165             type_line = type_line + ' mdpdb '
166
167        if instance.type == 'remd' or instance.type == 'dock':
168           control_line = 'SEED='+str(instance.md_seed)+' '\
169            +instance.md_start +' re one_letter'\
170            + ' BOXX='+str(instance.boxx) \
171            + ' BOXY='+str(instance.boxy) \
172            + ' BOXZ='+str(instance.boxz) +' '
173           
174           if instance.type == 'dock':
175              control_line = control_line+'CONSTR_HOMOL=20 '             
176            
177           if instance.md_pdbref:
178              control_line = control_line+'pdbref '
179           
180           cntrl_saxs=''
181           if instance.saxs_data != '':
182              lines=instance.saxs_data.split('\n')
183              nsaxs=0
184              for line in lines:
185               if len(line.split())==2:
186                nsaxs+=1 
187              cntrl_saxs=' nsaxs='+str(nsaxs)+' scal_rad='+str(instance.scal_rad) 
188              control_line += cntrl_saxs             
189             
190           type_line = '        reset_vel='+ str(instance.remd_nstex)\
191              +' nstep='+str(instance.md_nstep)\
192              +' ntwe='+str(instance.md_ntwe)\
193              +' ntwx='+str(instance.md_ntwx)\
194              +' dt='+str(instance.md_dt)
195           if instance.md_lang == 'langevin':
196             type_line = type_line + ' lang=1 '
197             type_line = type_line + 'scal_fric='+str(instance.md_scal_fric)
198           if instance.md_lang == 'berendsen':
199             type_line = type_line + ' lang=0 tbf '
200             type_line = type_line + 'tau_bath='+str(instance.md_tau)
201           if instance.md_respa:
202             type_line = type_line + ' respa '
203 #          if instance.md_mdpdb:
204 #            type_line = type_line + ' mdpdb '
205           
206           type_line2 = 'nrep='+str(instance.remd_nrep) \
207              +' nstex='+str(instance.remd_nstex) \
208              +' tlist mlist sync nsyn='+str(instance.remd_nstex)\
209              +' traj1file  rest1file TRAJCACHE=200'
210
211        if instance.md_start == 'pdbstart':
212              type_line = type_line + ' preminim cart'
213
214        wrapper = textwrap.TextWrapper(width=70)
215
216        word_list = wrapper.wrap(text=control_line)
217        for element in word_list[:-1]:
218            f.write('{:79}'.format(element)+'&\n')
219        f.write(word_list[-1]+'\n')
220
221        word_list = wrapper.wrap(text=type_line)
222        for element in word_list[:-1]:
223            f.write('{:79}'.format(element)+'&\n')
224        f.write(word_list[-1]+'\n')
225           
226        if instance.md_start == 'pdbstart':
227         if instance.type == 'remd' or instance.type == 'dock':
228          f.write('print_min_ini print_min_res print_min_stat maxfun='+
229            str(instance.min_maxfun)+'\n')
230         else:
231          f.write('print_min_ini print_min_res print_min_stat\n')
232
233        
234        if instance.type == 'remd' or instance.type == 'dock':
235           word_list = wrapper.wrap(text=type_line2)
236           for element in word_list[:-1]:
237               f.write('{:79}'.format(element)+'&\n')
238           f.write(word_list[-1]+'\n')
239
240           tmp1=json.loads(instance.remd_multi_t)
241           string1=" "
242           word_list = wrapper.wrap(text=string1.join(tmp1))
243           for element in word_list[:-1]:
244               f.write('{:79}'.format(element)+'&\n')
245           f.write(word_list[-1]+'\n')
246
247           tmp1=json.loads(instance.remd_multi_m)
248           string1=" "
249           word_list = wrapper.wrap(text=string1.join(tmp1))
250           for element in word_list[:-1]:
251               f.write('{:79}'.format(element)+'&\n')
252           f.write(word_list[-1]+'\n')
253        
254        if instance.saxs_data != '':
255          w+=' wsaxs='+str(instance.wsaxs)
256        f.write(w+'\n')
257
258
259        if instance.type == 'min' or instance.md_start == 'pdbstart':
260           if instance.type == 'dock':
261             f.write('model01.pdb\n')
262           else:  
263             f.write('plik.pdb\n')
264           f.write(write_ssbond(instance.ssbond))
265           if instance.md_2d == '':
266             f.write('0\n')
267           else:
268             f.write(seq_2d_write(seq_2d_xx(instance.md_2d,seq)))
269         
270        if instance.type == 'dock':
271           f.write('HOMOL_DIST=0.1    HOMOL_ANGLE=1.0 HOMOL_THETA=1.0 HOMOL_SCD=1.0 DIST_CUT=15.0  &\n')
272           f.write('READ2SIGMA START_FROM_MODELS\n')
273           for i in range(1,10):
274            f.write('model0'+str(i)+'.pdb\n')
275           for i in range(10,21):
276            f.write('model'+str(i)+'.pdb\n')
277
278        if instance.type == 'md' and instance.md_start != 'pdbstart':
279
280           if instance.md_pdbref:
281              f.write('plik.pdb\n')
282           f.write(str(len(seq))+'\n')
283           for i in range(0,len(seq),80):
284              f.write(seq[i:i+80]+'\n')
285           f.write('0\n')
286           if instance.md_2d == '':
287             f.write('0\n')
288           else:
289             f.write(seq_2d_write(seq_2d_xx(instance.md_2d,seq)))
290
291
292        if instance.type == 'remd' and instance.md_start != 'pdbstart':
293
294           if instance.md_pdbref:
295              f.write('plik.pdb\n')             
296           f.write(str(len(seq))+'\n')
297           for i in range(0,len(seq),80):
298              f.write(seq[i:i+80]+'\n')
299           f.write('0\n')
300           if instance.md_2d == '':
301             f.write('0\n')
302           else:
303             f.write(seq_2d_write(seq_2d_xx(instance.md_2d,seq)))
304
305
306        if instance.saxs_data != '':
307              lines=instance.saxs_data.split('\n')
308              for line in lines:
309               tmp=line.split()
310               if len(tmp)==2:
311                f.write(tmp[0]+' '+tmp[1]+'\n') 
312        
313         
314      if instance.type == 'min':
315       os.chdir(instance.jobdirname)
316       if  instance.unres_ff =='FF2':
317         shutil.copy('../files/pbs.csh','.')
318       else:
319         shutil.copy('../files/pbs_new.csh','pbs.csh')
320       ret_code = subprocess.Popen(' /opt/torque/bin/qsub pbs.csh', shell=True, stdout=subprocess.PIPE, stderr=subprocess.STDOUT)
321       os.chdir('..')
322      elif instance.type == 'md':
323       os.chdir(instance.jobdirname)
324       if  instance.unres_ff =='FF2':
325         f1 = open('../files/pbs_md.csh', 'r')
326       else:
327         f1 = open('../files/pbs_md_new.csh', 'r') 
328       f2 = open('pbs.csh', 'w')
329       for line in f1:
330            f2.write(line.replace('$dimen3',str(dimen3)+' '+str(instance.md_temp) ))
331       f1.close()
332       f2.close()  
333
334       f1 = open('file.seq', 'w')
335       for i in range(0,len(seq),80):
336              f1.write(seq[i:i+80]+'\n')
337       f1.write('\n')       
338       f1.close()
339
340       ret_code = subprocess.Popen(' /opt/torque/bin/qsub pbs.csh', shell=True, stdout=subprocess.PIPE, stderr=subprocess.STDOUT)      
341       os.chdir('..')
342      else:
343 # instance.type == 'remd'
344       tmp1=json.loads(instance.remd_multi_m)
345       nreplicas=sum(map(int, tmp1))
346 # hardcoded for 4 core nodes with FGPROC=2
347       nodes=int(nreplicas*2/4)
348       reszta = nreplicas*2-nodes*4
349       if reszta == 0 :
350         pbs = '#PBS -l nodes='+str(nodes)+':ppn=4'
351       else:
352         pbs = '#PBS -l nodes='+str(nodes)+':ppn=4+1:ppn='+str(reszta)
353         
354       tmp_ss=''
355       if instance.ssbond!='':
356         ix=[]
357         ii=0
358         for i in range(0,len(seq)):
359            if seq[i]=='X':
360               ii=ii+1
361            ix.append(ii)
362         l=[]
363         ss=json.loads(instance.ssbond)
364         for e in ss:
365            l.append(e[0]-ix[e[0]])
366            l.append(e[1]-ix[e[1]])
367         tmp_ss=' '.join(map(str,l))
368
369       os.chdir(instance.jobdirname)
370       if  instance.unres_ff =='FF2':
371         f1 = open('../files/pbs8.csh', 'r')
372       else:
373         f1 = open('../files/pbs8_new.csh', 'r')
374       f2 = open('pbs8.csh', 'w')
375       for line in f1:
376           if '#PBS -l nodes=4:ppn=4' in line:
377            f2.write(line.replace('#PBS -l nodes=4:ppn=4', pbs))
378           elif '$ssbond' in line:
379            f2.write(line.replace('$ssbond',tmp_ss))
380           elif 'cd $PBS_O_WORKDIR' in line and instance.type == 'dock':
381             f2.write(line)
382             f2.write('../files/generator_v7 plik.pdb plik2.pdb\n')
383           else:
384            tmp1=json.loads(instance.remd_multi_t)
385            f2.write(line.replace('$temperatures','"'+" ".join(tmp1)+'" '+str(nreplicas)))
386       f1.close()
387       f2.close()
388
389       f1 = open('file.seq', 'w')
390       for i in range(0,len(seq),80):
391              f1.write(seq[i:i+80]+'\n')
392       f1.write('\n')       
393       f1.close()
394       
395       os.chdir('..')  
396
397 # write wham & cluster input
398       with open(instance.jobdirname+'/file_wham.inp','w') as f:
399           if instance.md_nstep/instance.md_ntwx*nreplicas<=8000:
400               isampl=1
401           else:
402               isampl=int(instance.md_nstep/instance.md_ntwx*nreplicas/8000)
403           f.write('{:79}'.format(' isampl='+str(isampl))+'&\n')
404           if instance.md_pdbref:
405             f.write('{:79}'.format('SEED='+str(instance.md_seed)+
406                ' einicheck=1 rescale=2 delta=0.02 cxfile classify')+'&\n')
407           else:    
408             f.write('{:79}'.format('SEED='+str(instance.md_seed)+
409               ' einicheck=1 rescale=2 delta=0.02 cxfile')+'&\n')
410           if instance.md_2d != '':
411                f.write('{:79}'.format(' WITH_DIHED_CONSTR')+'&\n')
412           if instance.type == 'dock':
413                f.write('{:79}'.format(' n_ene=19 CONSTR_HOMOL=20')+'&\n')
414           f.write('BOXX='+str(instance.boxx)+' BOXY='+str(instance.boxy)+
415                     ' BOXZ='+str(instance.boxz)+cntrl_saxs +'\n')
416                  
417           f.write('nres='+str(len(seq))+' one_letter\n')
418           for i in range(0,len(seq),80):
419              f.write(seq[i:i+80]+'\n')
420
421           f.write(write_ssbond(instance.ssbond))
422
423           if instance.type == 'dock':
424             f.write('HOMOL_DIST=0.1    HOMOL_ANGLE=1.0 HOMOL_THETA=1.0 HOMOL_SCD=1.0 DIST_CUT=15.0  &\n')
425             f.write('READ2SIGMA \n')
426             for i in range(1,10):
427               f.write('model0'+str(i)+'.pdb\n')
428             for i in range(10,21):
429               f.write('model'+str(i)+'.pdb\n')
430
431           if instance.md_2d != '':
432              f.write(seq_2d_write_50(seq_2d_xx(instance.md_2d,seq)))
433
434           if instance.saxs_data != '':
435              saxs_r=[]
436              saxs_P=[]
437              saxs_norm=0
438              r_last=None
439              fsaxs = open(instance.jobdirname+'/saxs.data', 'w')
440              lines=instance.saxs_data.split('\n')
441              for line in lines:
442               tmp=line.split()
443               if len(tmp)==2:
444                f.write(tmp[0]+' '+tmp[1]+'\n')
445                saxs_r.append(tmp[0])
446                saxs_P.append(float(tmp[1]))
447                if r_last==None:
448                   r_last=float(tmp[0])
449                else:
450                   saxs_norm += float(tmp[1])*(float(tmp[0])-r_last)
451                   r_last=float(tmp[0])
452              for r,p in zip(saxs_r,saxs_P):
453                fsaxs.write('{} {}\n'.format(r,p/saxs_norm)) 
454              fsaxs.close()
455
456           f.write(w+'\n')
457           f.write('\n')
458           f.write('nt='+str(instance.remd_nrep)+' replica read_iset\n')
459           tmp1=json.loads(instance.remd_multi_t)
460           for element in tmp1:
461               f.write('nr=1 temp='+element+' fi=0.0 0.0 0.0 0.0 0.0\n')
462               f.write('kh= 0.0 Q0=0.0\n')
463           rec=instance.md_nstep/instance.md_ntwx
464           f.write('nfile_cx=1 rec_start='+str(rec/10)
465                   +' rec_end='+str(rec)+' totraj='+str(nreplicas)+'\n')
466           f.write('file_MD000\n')
467           if instance.md_pdbref:
468              f.write('pdbref\n')
469              f.write('plik.pdb\n')
470
471            
472       with open(instance.jobdirname+'/file_cluster.inp','w') as f:
473           f.write(instance.name + ' UNRES server job'+ '\n')
474           f.write('{:79}'.format('nres='+str(len(seq))
475              +' one_letter rescale=2 PRINT_CART PDBOUT=1 iopt=1'
476              +' temper='+str(instance.remd_cluter_temp))+'&\n')
477           if instance.md_2d != '':
478                f.write('{:79}'.format(' WITH_DIHED_CONSTR')+'&\n')
479           if instance.type == 'dock':
480                f.write('{:79}'.format(' CONSTR_HOMOL=20')+'&\n')
481           f.write('BOXX='+str(instance.boxx)+' BOXY='+str(instance.boxy)+
482                   ' BOXZ='+str(instance.boxz)+cntrl_saxs +'\n')
483           f.write(w+'\n')
484           for i in range(0,len(seq),80):
485              f.write(seq[i:i+80]+'\n')
486           f.write(write_ssbond(instance.ssbond))
487           if instance.type == 'dock':
488             f.write('HOMOL_DIST=0.1    HOMOL_ANGLE=1.0 HOMOL_THETA=1.0 HOMOL_SCD=1.0 DIST_CUT=15.0  &\n')
489             f.write('READ2SIGMA \n')
490             for i in range(1,10):
491               f.write('model0'+str(i)+'.pdb\n')
492             for i in range(10,21):
493               f.write('model'+str(i)+'.pdb\n')
494           if instance.md_2d != '':
495              f.write(seq_2d_write_50(seq_2d_xx(instance.md_2d,seq)))
496           if instance.saxs_data != '':
497              lines=instance.saxs_data.split('\n')
498              for line in lines:
499               tmp=line.split()
500               if len(tmp)==2:
501                f.write(tmp[0]+' '+tmp[1]+'\n') 
502       
503       os.chdir(instance.jobdirname)  
504       ret_code = subprocess.Popen(' /opt/torque/bin/qsub pbs8.csh', shell=True, stdout=subprocess.PIPE, stderr=subprocess.STDOUT)
505       os.chdir('..')
506      
507      print ret_code.stdout.readlines()
508      logger.warning("Dir \"%s\" has been saved." % instance.jobdirname)
509      
510 @receiver(post_delete, sender=Task)
511 def delete_on_task_del(sender, instance, **kwargs):     
512    if instance.jobdirname:
513       if os.path.isdir(instance.jobdirname):
514          shutil.rmtree(instance.jobdirname)
515          logger.warning("Dir \"%s\" has been rm." % instance.jobdirname)         
516