Passed
Push — master ( a2a58c...34f3b6 )
by Marcin
06:20 queued 14s
created

has_high_rna_content()   A

Complexity

Conditions 4

Size

Total Lines 20
Code Lines 12

Duplication

Lines 0
Ratio 0 %

Importance

Changes 0
Metric Value
cc 4
eloc 12
nop 2
dl 0
loc 20
rs 9.8
c 0
b 0
f 0
1
#!/usr/bin/env python
2
# -*- coding: utf-8 -*-
3
4
"""rna_pdb_tools - a swiss army knife to manipulation of RNA pdb structures
5
6
Usage::
7
8
   $ rna_pdb_tools.py --delete A:46-56 --inplace *.pdb
9
10
    $ rna_pdb_tools.py --get-seq *
11
    # BujnickiLab_RNApuzzle14_n01bound
12
    > A:1-61
13
    # BujnickiLab_RNApuzzle14_n02bound
14
    > A:1-61
15
    CGUUAGCCCAGGAAACUGGGCGGAAGUAAGGCCCAUUGCACUCCGGGCCUGAAGCAACGCG
16
    [...]
17
18
See `rna_pdb_merge_into_one.py` to merge PDB files in the order as you like into one NMR-like (multimodel) file
19
20
Examples::
21
22
    rna_pdb_tools.py --backbone-only   --get-rnapuzzle-ready  --inplace --suffix=bo examples/4GXY_min.pdb
23
24
To extract specific atoms for each residue and write them to separate PDB file (next to the input files, following syntax "<input>_<resid>.pdb")::
25
26
    rna_pdb_tools.py --rpr input/4GXY_min.pdb --save-single-res --ref-frame-only
27
28
Atoms presets::
29
30
  --backbone-only       used only with --get-rnapuzzle-ready, keep only backbone (= remove bases)
31
  --ref-frame-only      used only with --get-rnapuzzle-ready, keep only reference frames, OP1 OP2 P
32
  --no-backbone         used only with --get-rnapuzzle-ready, remove atoms of backbone (define as P OP1 OP2 O5')
33
  --bases-only          used only with --get-rnapuzzle-ready, keep only atoms of bases
34
      
35
.. image:: ../pngs/276411138-236435ff-2944-4bec-ab75-dca0d1e3aacf.jpg
36
37
38
-v is for verbose, --version for version ;)
39
"""
40
import argparse
41
import textwrap
42
import os
43
import shutil
44
import sys
45
import tempfile
46
import glob
47
import os
48
from rna_tools.rna_tools_lib import edit_pdb, add_header, get_version, \
49
                          collapsed_view, fetch, fetch_ba, replace_chain, RNAStructure, \
50
                          select_pdb_fragment, sort_strings, set_chain_for_struc
51
from rna_tools.tools.rna_x3dna.rna_x3dna import x3DNA
52
53
54
def get_parser():
55
    version = os.path.basename(os.path.dirname(os.path.abspath(__file__))), get_version(__file__)
56
    version = version[1].strip()
57
    parser = argparse.ArgumentParser(
58
        description=__doc__, formatter_class=argparse.RawTextHelpFormatter)
59
60
    parser.add_argument('--version', help='', action='version', version=version)
61
62
    parser.add_argument('-r', '--report', help='get report',
63
                        action='store_true')
64
65
    parser.add_argument('--no-progress-bar', help='for --no-progress-bar for --rpr', action='store_true')
66
67
    parser.add_argument('--renum-atoms', help='renumber atoms, tested with --get-seq',
68
                        action='store_true')
69
70
    parser.add_argument('--renum-nmr', help='',
71
                        action='store_true')
72
73
    parser.add_argument('--renum-residues-dirty', help='',  action='store_true')
74
75
    parser.add_argument('--undo', help='undo operation of action done --inplace, , rename "backup files" .pdb~ to pdb, ALL files in the folder, not only ~ related to the last action (that you might want to revert, so be careful)',  action='store_true')
76
77
    parser.add_argument('--delete-anisou', help='remove files with ANISOU records, works with --inplace',
78
                        action='store_true')
79
80
    parser.add_argument('--remove0', help='remove atoms of X=0 position',
81
                        action='store_true')
82
83
    parser.add_argument('--fix', help='fix a PDB file, ! external program, pdbfixer used to fix missing atoms',
84
                        action='store_true')
85
86
    parser.add_argument('--to-mol2', help='fix a PDB file, ! external program, pdbfixer used to fix missing atoms',
87
                        action='store_true')
88
89
    parser.add_argument('--split-alt-locations', help='splits atoms, e.g. for alt locs A and B, it splits atoms two MODELS (all localizations A goes into MODEL1 and all localizations B goes into MODEL2',
90
                        action='store_true')
91
92
    parser.add_argument('-c', '--clean', help='get clean structure',
93
                        action='store_true')
94
95
    parser.add_argument('--is-pdb', help='check if a file is in the pdb format',
96
                        action='store_true')
97
98
    parser.add_argument('--is-nmr', help='check if a file is NMR-style multiple model pdb',
99
                        action='store_true')
100
101
    parser.add_argument('--nmr-dir', help='make NMR-style multiple model pdb file from a set of files \n\n' +
102
                        "  rna_pdb_tools.py --nmr-dir . 'cwc15_u5_fragments*.pdb' > ~/Desktop/cwc15-u5.pdb\n\n" +
103
                        "please use '' for pattern file recognition, this is a hack to deal with folders with\n"
104
                        "thousands of models, if you used only *.pdb then the terminal will complain that you\n"
105
                        "selected to many files.")
106
107
    parser.add_argument('--un-nmr', help="""split NMR-style multiple model pdb files into individual models [biopython],
108
109
   rna_pdb_tools.py  --un-nmr  split.pdb
110
   2
111
   /Users/magnus/Desktop/3hl2/split_1.pdb
112
   /Users/magnus/Desktop/3hl2/split_2.pdb
113
114
""",
115
                        action='store_true')
116
117
    parser.add_argument('--orgmode', help='get a structure in org-mode format <sick!>',
118
                        action='store_true')
119
120
    parser.add_argument('--get-chain', help='get chain, one or many, e.g, A, but now also ABC works')
121
122
    parser.add_argument('--fetch', action='store_true', help='fetch file from the PDB db, e.g., 1xjr,\nuse \'rp\' to fetch, fetch a given join, 4w90:C or 4w90_C'  +
123
                        'the RNA-Puzzles standardized_dataset [around 100 MB]')
124
125
    parser.add_argument('--fetch-ba', action='store_true',
126
                        help='fetch biological assembly from the PDB db')
127
128
    parser.add_argument('--fetch-chain', action='store_true', help='fetch a structure in extract chain, e.g. 6bk8 H')
129
130
    parser.add_argument('--fetch-fasta', action='store_true', help='fetch a fasta/sequence for given PDB ID, e.g. 6bk8')
131
132
    parser.add_argument('--get-seq', help='get seq', action='store_true')
133
134
    parser.add_argument('--color-seq', help='color seq, works with --get-seq', action='store_true')
135
136
    parser.add_argument('--ignore-files', help='files to be ingored, .e.g, \'solution\'')
137
    
138
    parser.add_argument('--compact',
139
                        help=textwrap.dedent("""with --get-seq, get it in compact view'
140
$ rna_pdb_tools.py --get-seq --compact *.pdb
141
# 20_Bujnicki_1
142
ACCCGCAAGGCCGACGGCGCCGCCGCUGGUGCAAGUCCAGCCACGCUUCGGCGUGGGCGCUCAUGGGU # A:1-68
143
# 20_Bujnicki_2
144
ACCCGCAAGGCCGACGGCGCCGCCGCUGGUGCAAGUCCAGCCACGCUUCGGCGUGGGCGCUCAUGGGU # A:1-68
145
# 20_Bujnicki_3
146
ACCCGCAAGGCCGACGGCGCCGCCGCUGGUGCAAGUCCAGCCACGCUUCGGCGUGGGCGCUCAUGGGU # A:1-68
147
# 20_Bujnicki_4
148
149
"""), action='store_true')
150
151
    parser.add_argument('--hide-warnings', help='hide warnings, works with --get-chain, it hides warnings that given changes are not detected in a PDB file', action='store_true')
152
153
    parser.add_argument('--get-ss', help='get secondary structure', action='store_true')
154
155
    parser.add_argument('--rosetta2generic', help='convert ROSETTA-like format to a generic pdb',
156
                        action='store_true')
157
158
    parser.add_argument('--no-hr', help='do not insert the header into files',
159
                        action='store_true')
160
161
    parser.add_argument('--renumber-residues', help='by defult is false',
162
                        action='store_true')
163
164
    parser.add_argument('--dont-rename-chains',
165
                        help=textwrap.dedent("""used only with --get-rnapuzzle-ready.
166
By default:
167
   --get-rnapuzzle-ready rename chains from ABC.. to stop behavior switch on this option
168
"""),
169
                        action='store_true')
170
171
    parser.add_argument('--dont-fix-missing-atoms',
172
                        help="""used only with --get-rnapuzzle-ready""",
173
                        action='store_true')
174
175
    parser.add_argument('--inspect',
176
                        help="inspect missing atoms (technically decorator to --get-rnapuzzle-ready without actually doing anything but giving a report on problems)",
177
                        action='store_true')
178
    
179
    parser.add_argument('--collapsed-view', help='',
180
                        action='store_true')
181
182
    parser.add_argument('--cv', help='alias to collapsed_view',
183
                        action='store_true')
184
185
    parser.add_argument('-v', '--verbose', help='tell me more what you\'re doing, please!',
186
                        action='store_true')
187
188
    parser.add_argument('--mutate', help=textwrap.dedent("""mutate residues,
189
e.g.,
190
      --mutate "A:1a+2a+3a+4a,B:1a"
191
to mutate to adenines the first four nucleotides of the chain A
192
and the first nucleotide of the chain B"""))
193
194
    parser.add_argument('--edit',
195
                        dest="edit",
196
                        default='',
197
                        help="edit 'A:6>B:200', 'A:2-7>B:2-7'")
198
199
    parser.add_argument('--rename-chain',
200
                        help="edit 'A>B' to rename chain A to chain B")
201
202
    parser.add_argument('--swap-chains', help='B>A, rename A to _, then B to A, then _ to B')
203
204
    parser.add_argument('--set-chain', help='set chain for all ATOM lines and TER (quite brutal function)')
205
206
    parser.add_argument('--replace-chain',
207
                        default='',
208
                        help=textwrap.dedent("""a file PDB name with one chain that will be used to
209
replace the chain in the original PDB file,
210
the chain id in this file has to be the same with the chain id of the original chain"""))
211
212
    parser.add_argument('--delete',  # type="string",
213
                        dest="delete",
214
                        default='',
215
                        help="delete the selected fragment, e.g. A:10-16, or for more than one fragment --delete 'A:1-25+30-57'")
216
217
    parser.add_argument('--extract',  # type="string",
218
                        dest="extract",
219
                        default='',
220
                        help="extract the selected fragment, e.g. A:10-16, or for more than one fragment --extract 'A:1-25+30-57', or even 'A:1-25+B:30-57'")
221
222
    parser.add_argument('--extract-chain',
223
                         help="extract chain, e.g. A")
224
225
    parser.add_argument('--uniq', help=textwrap.dedent("""
226
rna_pdb_tools.py --get-seq --uniq '[:5]' --compact --chain-first * | sort
227
A:1-121        ACCUUGCGCAACUGGCGAAUCCUGGGGCUGCCGCCGGCAGUACCC...CA # rp13nc3295_min.out.1
228
A:1-123        ACCUUGCGCGACUGGCGAAUCCUGAAGCUGCUUUGAGCGGCUUCG...AG # rp13cp0016_min.out.1
229
A:1-123        ACCUUGCGCGACUGGCGAAUCCUGAAGCUGCUUUGAGCGGCUUCG...AG # zcp_6537608a_ALL-000001_AA
230
A:1-45 57-71   GGGUCGUGACUGGCGAACAGGUGGGAAACCACCGGGGAGCGACCCGCCGCCCGCCUGGGC # solution
231
"""))
232
233
    parser.add_argument('--chain-first', help="", action='store_true')
234
    parser.add_argument('--oneline', help="", action='store_true')
235
236
    parser.add_argument('--replace-htm', help="", action='store_true')
237
238
    parser.add_argument('--fasta',
239
                        help= textwrap.dedent("""with --get-seq, show sequences in fasta format,
240
can be combined with --compact (mind, chains will be separated with ' ' in one line)
241
242
$ rna_pdb_tools.py --get-seq --fasta --compact input/20_Bujnicki_1.pdb
243
> 20_Bujnicki_1
244
ACCCGCAAGGCCGACGGC GCCGCCGCUGGUGCAAGUCCAGCCACGCUUCGGCGUGGGCGCUCAUGGGU
245
246
"""), action='store_true')
247
248
    parser.add_argument('--cif2pdb', help="[PyMOL Python package required]", action='store_true')
249
    parser.add_argument('--pdb2cif', help="[PyMOL Python package required]", action='store_true')
250
251
    parser.add_argument('--mdr', help='get structures ready for MD (like rpr but without first)',
252
                        action='store_true')
253
254
    x = parser.add_argument_group('RNAPUZZLE-READY')
255
    x.add_argument('--get-rnapuzzle-ready',
256
                        help=textwrap.dedent("""get RNApuzzle ready (keep only standard atoms).'
257
Be default it does not renumber residues, use --renumber-residues
258
[requires BioPython]"""), action='store_true')
259
260
    x.add_argument('--rpr', help='alias to get_rnapuzzle ready)',
261
                        action='store_true')
262
263
264
    rpr = parser.add_argument_group('CAN BE COMBINED WITH')# --get-rnapuzzle-ready (or --rpr) can be combined with')
265
266
    rpr.add_argument('--keep-hetatm', help='keep hetatoms', action='store_true')
267
268
    rpr.add_argument('--inplace', help=textwrap.dedent("""in place edit the file! [experimental,
269
only for get_rnapuzzle_ready, --delete, --get-ss, --get-seq, --edit-pdb]"""),
270
                        action='store_true')
271
272
    rpr.add_argument('--here', help=textwrap.dedent("""save a file next to the original file with auto suffix
273
for --extract it's .extr.pdb"""),
274
                        action='store_true')
275
276
    rpr.add_argument('--suffix', help=textwrap.dedent("""when used with --inplace allows you to change a name of a new file, --suffix del will give <file>_del.pdb (mind added _)"""))
277
278
    rpr.add_argument('--replace-hetatm', help="replace 'HETATM' with 'ATOM' [tested only with --get-rnapuzzle-ready]",
279
                        action="store_true")
280
281
    rpr.add_argument('--dont-report-missing-atoms',
282
                        help="""used only with --get-rnapuzzle-ready""",
283
                        action='store_true')
284
285
    rpr.add_argument('--backbone-only',
286
                        help="""used only with --get-rnapuzzle-ready, keep only backbone (= remove bases)""",
287
                        action='store_true')
288
289
    rpr.add_argument('--p-only',
290
                        help="""used only with --get-rnapuzzle-ready, keep p backbone (= remove bases)""",
291
                        action='store_true')
292
293
    rpr.add_argument('--ref-frame-only',
294
                        help="""used only with --get-rnapuzzle-ready, keep only reference frames, OP1 OP2 P""",
295
                        action='store_true')
296
297
    rpr.add_argument('--no-backbone',
298
                        help="""used only with --get-rnapuzzle-ready, remove atoms of backbone (define as P OP1 OP2 O5')""",
299
                        action='store_true')
300
301
    rpr.add_argument('--bases-only',
302
                        help="""used only with --get-rnapuzzle-ready, keep only atoms of bases""",
303
                        action='store_true')
304
305
    rpr.add_argument('--save-single-res',
306
                        help="""used only with --get-rnapuzzle-ready, for each residue create a new pdb output file, you can combine it with --bases-only etc.""",
307
                        action='store_true')
308
309
    parser.add_argument('file', help='file', nargs='+')
310
    #parser.add_argument('outfile', help='outfile')
311
    return parser
312
313
314
# main
315
if __name__ == '__main__':
316
    # get version
317
    version = os.path.basename(os.path.dirname(os.path.abspath(__file__))), get_version(__file__)
318
    version = version[1].strip()
319
320
    # get parser and arguments
321
    parser = get_parser()
322
323
    args = parser.parse_args()
324
325
    # quick fix for one files vs file-s
326
    if list == type(args.file) and len(args.file) == 1:
327
        args.file = args.file[0]
328
329
    if args.report:
330
        s = RNAStructure(args.file)
331
        print(s.get_report)
332
        print(s.get_preview())
333
        print(s.get_info_chains())
334
335
    if args.clean:
336
        s = RNAStructure(args.file)
337
        s.decap_gtp()
338
        s.std_resn()
339
        s.remove_hydrogen()
340
        s.remove_ion()
341
        s.remove_water()
342
        s.renum_atoms()
343
        s.fix_O_in_UC()
344
        s.fix_op_atoms()
345
        # print s.get_preview()
346
        # s.write(args.outfile)
347
        if not args.no_hr:
348
            print(add_header(version))
349
        print(s.get_text())
350
351
    if args.get_seq:
352
        # quick fix - make a list on the spot
353
        if list != type(args.file):
354
            args.file = [args.file]
355
        ##################################
356
        analyzed = []
357
        for f in args.file:
358
            #####################################
359
            if args.uniq:
360
                subname = eval('f' + args.uniq)
361
                if subname in analyzed:
362
                    continue
363
                else:
364
                    analyzed.append(subname)
365
            ########
366
            s = RNAStructure(f)
367
            if not s.is_pdb():
368
                print('Error: Not a PDB file %s' % f)
369
                sys.exit(1)
370
            s.decap_gtp()
371
            s.std_resn()
372
            s.remove_hydrogen()
373
            s.remove_ion()
374
            s.remove_water()
375
            if args.renum_atoms:
376
                s.renum_atoms()
377
            s.fix_O_in_UC()
378
            s.fix_op_atoms()
379
380
            output = ''
381
382
            # with # is easier to grep this out
383
            if args.fasta:
384
                # s.fn vs s.name
385
                output += s.get_seq(compact=args.compact, chainfirst=args.chain_first, fasta=args.fasta, addfn=s.name, color=args.color_seq) + '\n'
386
            elif args.oneline:
387
                output += s.get_seq(compact=args.compact, chainfirst=args.chain_first, color=args.color_seq).strip() + ' # '+ os.path.basename(f.replace('.pdb', '')) + '\n'
388
            else:
389
                output += '# ' + os.path.basename(f.replace('.pdb', '')) + '\n'
390
                output += s.get_seq(compact=args.compact, chainfirst=args.chain_first, color=args.color_seq) + '\n'
391
392
            try:
393
                sys.stdout.write(output)
394
                sys.stdout.flush()
395
            except IOError:
396
                pass
397
398
    if args.get_ss:
399
        # quick fix - make a list on the spot
400
        if list != type(args.file):
401
            args.file = [args.file]
402
        ##################################
403
        for f in args.file:
404
            output = f + '\n'
405
            p = x3DNA(f)
406
            output += p.get_secstruc() + '\n'
407
            try:
408
                sys.stdout.write(output)
409
                sys.stdout.flush()
410
            except IOError:
411
                pass
412
413
    # getchain
414
    if args.get_chain:
415
        s = RNAStructure(args.file)
416
        ## s.std_resn()
417
        ## s.remove_hydrogen()
418
        ## s.remove_ion()
419
        ## s.remove_water()
420
        ## s.renum_atoms()
421
        ## s.fix_O_in_UC()
422
        ## s.fix_op_atoms()
423
        # print s.get_preview()
424
        warningmsg = ''
425
        for chain in list(args.get_chain):
426
            chain_txt = s.get_chain(chain)
427
            if not chain_txt.strip():
428
                warningmsg += 'Warning: Chain %s not detected!' % chain
429
            else:
430
                print(chain_txt)
431
        if not args.hide_warnings:
432
            if warningmsg:
433
                print(warningmsg)
434
435
    if args.rosetta2generic:
436
        s = RNAStructure(args.file)
437
        s.std_resn()
438
        s.remove_hydrogen()
439
        s.remove_ion()
440
        s.remove_water()
441
        s.fix_op_atoms()
442
        s.renum_atoms()
443
        # print s.get_preview()
444
        # s.write(args.outfile)
445
        if not args.no_hr:
446
            print(add_header(version))
447
        print(s.get_text())
448
449
    remarks_only = False
450
    if args.inspect:
451
        args.get_rnapuzzle_ready = True
452
        fix_missing_atom = False
453
        remarks_only = True
454
        args.dont_rename_chains = True
455
        args.renumber_residues = False
456
457
    if args.get_rnapuzzle_ready or args.rpr or args.mdr:  # rnapuzzle
458
        # quick fix - make a list on the spot
459
        if list != type(args.file):
460
            args.file = [args.file]
461
        ##################################
462
        # progress bar only in --inplace mode!
463
        if args.inplace:
464
            if not args.no_progress_bar:
465
                import progressbar
466
                bar = progressbar.ProgressBar(max_value=len(args.file))
467
                bar.update(0)
468
469
        for c, f in enumerate(args.file):
470
            if args.verbose:
471
                print(f)
472
            if args.inplace:
473
                shutil.copy(f, f + '~')
474
475
            # keep previous edits
476
            previous_edits = []
477
            with open(f) as fx:
478
                for l in fx:
479
                    if l.startswith('HEADER --'):
480
                        previous_edits.append(l.strip())
481
            ######################
482
483
            s = RNAStructure(f)
484
            if args.replace_hetatm:
485
                s.replace_hetatms()
486
487
            s.remove_hydrogen()
488
            s.decap_gtp()
489
            s.std_resn()
490
            s.fix_op_atoms()
491
492
            s.remove_ion()
493
            s.remove_water()
494
            # s.renum_atoms()
495
496
            s.shift_atom_names()
497
            s.prune_elements()
498
499
            # s.write('tmp.pdb')
500
            
501
            rename_chains = False if args.dont_rename_chains else True
502
503
            report_missing_atoms = not args.dont_report_missing_atoms
504
            fix_missing_atom = not args.dont_fix_missing_atoms
505
506
            ignore_op3 = False
507
            if args.mdr:
508
                ignore_op3 = True
509
                
510
            remarks = s.get_rnapuzzle_ready(args.renumber_residues, fix_missing_atoms=fix_missing_atom,
511
                                            rename_chains=rename_chains,
512
                                            report_missing_atoms=report_missing_atoms,
513
                                            backbone_only=args.backbone_only,
514
                                            no_backbone=args.no_backbone,
515
                                            bases_only=args.bases_only,
516
                                            p_only=args.p_only,
517
                                            keep_hetatm=args.keep_hetatm,
518
                                            ignore_op3=ignore_op3,
519
                                            save_single_res=args.save_single_res,
520
                                            ref_frame_only = args.ref_frame_only,
521
                                            verbose=args.verbose)
522
523
            if args.inplace:
524
                if args.suffix:
525
                    f = f.replace('.pdb', '_' + args.suffix + '.pdb')
526
                    if args.verbose: print(f)
527
                        
528
                with open(f, 'w') as f:
529
                    if not args.no_hr:
530
                        f.write(add_header(version) + '\n')
531
                    if previous_edits:
532
                        f.write('\n'.join(previous_edits) + '\n')
533
                    if remarks:
534
                        f.write('\n'.join(remarks) + '\n')
535
                    f.write(s.get_text())
536
                # progress bar only in --inplace mode!
537
                if not args.no_progress_bar:
538
                   bar.update(c)
0 ignored issues
show
introduced by
The variable bar does not seem to be defined for all execution paths.
Loading history...
539
540
            else:
541
                output = ''
542
                if not args.no_hr:
543
                    output += add_header(version) + '\n'
544
                if remarks:
545
                    output += '\n'.join(remarks) + '\n'
546
                output += s.get_text() + '\n'
547
548
                if remarks_only:
549
                    sys.stdout.write('\n'.join(remarks))
550
                    sys.stdout.flush()
551
                else:
552
                    if args.here:
553
                        if '_rpr' not in f:  # good idea?
554
                            nf = f.replace('.pdb', '_rpr.pdb')
555
                            with open(nf, 'w') as fio:
556
                                print(nf)
557
                                fio.write(output)
558
                    else:
559
                        try:
560
                            sys.stdout.write(output)
561
                            sys.stdout.flush()
562
                        except IOError:
563
                            pass
564
        # hmm... fix for problem with renumbering, i do renumbering
565
        # and i stop here
566
        # i'm not sure that this is perfect
567
        sys.exit(0)
568
569
570
    if args.renumber_residues:
571
        s = RNAStructure(args.file)
572
        s.remove_hydrogen()
573
        s.remove_ion()
574
        s.remove_water()
575
        s.get_rnapuzzle_ready(args.renumber_residues)
576
        s.renum_atoms()
577
        if not args.no_hr:
578
            print(add_header(version))
579
        print(s.get_text())
580
581
# args.renumber_resides_dirty
582
    if args.renum_residues_dirty:
583
        # quick fix - make a list on the spot
584
        if list != type(args.file):
585
            args.file = [args.file]
586
        ##################################
587
        for f in args.file:
588
            if args.inplace:
589
                shutil.copy(f, f + '~')
590
591
            s = RNAStructure(f)
592
593
            output = ''
594
            #if not args.no_hr:
595
            #    output += add_header(version) + '\n'
596
            #    output += 'HEADER --delete ' + args.delete + '\n'  # ' '.join(str(selection))
597
            c = 1
598
            old_resi = -1
599
            for l in s.lines:
600
                if l.startswith('ATOM') or l.startswith("HETATOM"):
601
                    resi = int(l[23:26].strip())
602
                    if resi != old_resi:
603
                        old_resi = resi
604
                        c += 1
605
                    # print(resi, c)
606
                    #resi = c
607
                    #if chain in selection:
608
                    #    if resi in selection[chain]:
609
                    #        continue  # print chain, resi
610
                    output += l[:23] + str(c).rjust(3) + l[26:] + '\n'
611
            # write: inplace
612
            if args.inplace:
613
                with open(f, 'w') as f:
614
                    f.write(output)
615
            else:  # write: to stdout
616
                try:
617
                    sys.stdout.write(output)
618
                    sys.stdout.flush()
619
                except IOError:
620
                    pass
621
622
623
    if args.undo:
624
        # quick fix - make a list on the spot
625
        dir = args.file #os.path.abso(os.path.dirname(args.file))
626
        for f in glob.glob(dir + '/*.pdb~'):
627
            if args.verbose:
628
                print(f, '->',f.replace('.pdb~', '.pdb'))
629
            os.rename(f, f.replace('.pdb~', '.pdb'))
630
631
    if args.delete:
632
        # quick fix - make a list on the spot
633
        if list != type(args.file):
634
            args.file = [args.file]
635
        ##################################
636
        for f in args.file:
637
            if args.inplace:
638
                shutil.copy(f, f + '~')
639
640
            selection = select_pdb_fragment(args.delete)
641
            s = RNAStructure(f)
642
643
            output = ''
644
            if not args.no_hr:
645
                output += add_header(version) + '\n'
646
                output += 'HEADER --delete ' + args.delete + '\n'  # ' '.join(str(selection))
647
            for l in s.lines:
648
                if l.startswith('ATOM'):
649
                    chain = l[21]
650
                    resi = int(l[23:26].strip())
651
                    if chain in selection:
652
                        if resi in selection[chain]:
653
                            continue  # print chain, resi
654
                    output += l + '\n'
655
656
            # write: inplace
657
            if args.inplace:
658
                os.rename(f, f.replace('.pdb', '.pdb~'))                
659
                if args.suffix:
660
                    f = f.replace('.pdb', '_' + args.suffix + '.pdb')
661
                with open(f, 'w') as f:
662
                    f.write(output)
663
            else:  # write: to stdout
664
                try:
665
                    sys.stdout.write(output)
666
                    sys.stdout.flush()
667
                except IOError:
668
                    pass
669
670
    if args.replace_chain:
671
        # quick fix - make a list on the spot
672
        if list != type(args.file):
673
            args.file = [args.file]
674
        ##################################
675
        for f in args.file:
676
            if args.inplace:
677
                shutil.copy(f, f + '~')
678
679
            # --replace_chain <file> without A:<file> it will be easier than --x "A:<file>"
680
            s = RNAStructure(args.replace_chain)
681
            chain_ids = (s.get_all_chain_ids())
682
            if len(chain_ids) > 1:
683
                raise Exception('There is more than one chain in the inserted PDB file. There should be only one chain, the one you want to insert to the PDB.')
684
            out = replace_chain(f, args.replace_chain, list(chain_ids)[0])
685
            print(out)
686
687
    if args.mutate:
688
        # quick fix - make a list on the spot
689
        if list != type(args.file):
690
            args.file = [args.file]
691
        ##################################
692
        from rna_tools.tools.mini_moderna3.moderna import *
693
694
        for f in args.file:
695
            if args.ignore_files:
696
                if args.ignore_files in f:
697
                    continue
698
                
699
            if args.inplace:
700
                shutil.copy(f, f + '~')  # create a backup copy if inplace
701
702
            # create a working copy of the main file
703
            ftf = tempfile.NamedTemporaryFile(delete=False).name  # f temp file
704
            shutil.copy(f, ftf)  # create a backup copy if inplace
705
706
            # go over each chain
707
            # rna_pdb_tools.py --mutate 'A:1CB:1G,A:1U+B:1A' CG_AB.pdb > ~/Desktop/a.pdb
708
            args.mutate = args.mutate.upper()
709
            for m in args.mutate.split(','):  # A:1A, B:1A
710
                chain, resi_mutate_to = m.strip().split(':')  # A:1A+2C
711
                resi_mutate_to_list = resi_mutate_to.split('+')  # A:1A+2C
712
713
                model = load_model(f, chain)
0 ignored issues
show
Comprehensibility Best Practice introduced by
The variable load_model does not seem to be defined.
Loading history...
714
                # go over each mutation in this chain
715
                for resi_mutate_to in resi_mutate_to_list:
716
                    resi = resi_mutate_to[:-1]
717
                    mutate_to = resi_mutate_to[-1]
718
                    if args.verbose:
719
                        print('Mutate', model[resi], 'to', mutate_to)
720
                    exchange_single_base(model[resi], mutate_to)
0 ignored issues
show
Comprehensibility Best Practice introduced by
The variable exchange_single_base does not seem to be defined.
Loading history...
721
722
                # multi mutation in one chain
723
                tf = tempfile.NamedTemporaryFile(delete=False)
724
                model.write_pdb_file(tf.name)
725
726
                # work on the copy of f, ftf
727
                output = replace_chain(ftf, tf.name, chain)
728
                with open(ftf, 'w') as tmp:
729
                    tmp.write(output)
730
731
            # write: inplace
732
            if args.inplace:
733
                # ftf now is f, get ready for the final output
734
                if args.suffix:
735
                    f = f.replace('.pdb', '_' + args.suffix + '.pdb')
736
                # rpr on the file
737
                shutil.copy(ftf, f)
738
                os.system('rna_pdb_tools.py --rpr --no-progress-bar --inplace ' + f)
739
            else:  # write: to stdout
740
                try:
741
                    sys.stdout.write(output)
0 ignored issues
show
introduced by
The variable output does not seem to be defined for all execution paths.
Loading history...
742
                    sys.stdout.flush()
743
                except IOError:
744
                    pass
745
746
747
# extract
748
    if args.extract:
749
        # quick fix - make a list on the spot
750
        if list != type(args.file):
751
            args.file = [args.file]
752
        ##################################
753
        for f in args.file:
754
            if args.inplace:
755
                shutil.copy(f, f + '~')
756
757
            selection = select_pdb_fragment(args.extract)
758
            s = RNAStructure(f)
759
760
            output = ''
761
            if not args.no_hr:
762
                output += add_header(version) + '\n'
763
                output += 'HEADER extract ' + args.extract + '\n'  # ' '.join(str(selection))
764
            for l in s.lines:
765
                if l.startswith('ATOM'):
766
                    chain = l[21]
767
                    resi = int(l[23:26].strip())
768
                    if chain in selection:
769
                        if resi in selection[chain]:
770
                            # continue  # print chain, resi
771
                            output += l + '\n'
772
773
            # write: inplace
774
            if args.inplace:
775
                with open(f, 'w') as f:
776
                    f.write(output)
777
            elif args.here:
778
                if '_extr' not in f:  # good idea?
779
                    nf = f.replace('.pdb', '_extr.pdb')
780
                    with open(nf, 'w') as fio:
781
                        print(nf)
782
                        fio.write(output)
783
            else:  # write: to stdout
784
                try:
785
                    sys.stdout.write(output)
786
                    sys.stdout.flush()
787
                except IOError:
788
                    pass
789
790
791
    if args.extract_chain:
792
        # quick fix - make a list on the spot
793
        if list != type(args.file):
794
            args.file = [args.file]
795
        ##################################
796
        for f in args.file:
797
            if args.inplace:
798
                shutil.copy(f, f + '~')
799
800
            selection = select_pdb_fragment(args.extract)
801
            s = RNAStructure(f)
802
803
            output = ''
804
            if not args.no_hr:
805
                output += add_header(version) + '\n'
806
                output += 'HEADER extract ' + args.extract + '\n'  # ' '.join(str(selection))
807
            for l in s.lines:
808
                if l.startswith('ATOM') or l.startswith('TER') or l.startswith('HETATM'):
809
                    chain = l[21]
810
                    if chain in args.extract_chain:
811
                        output += l + '\n'
812
813
            # write: inplace
814
            if args.inplace:
815
                with open(f, 'w') as f:
816
                    f.write(output)
817
            else:  # write: to stdout
818
                try:
819
                    sys.stdout.write(output)
820
                    sys.stdout.flush()
821
                except IOError:
822
                    pass
823
824
    if args.is_pdb:
825
        s = RNAStructure(args.file)
826
        output = str(s.is_pdb())
827
        sys.stdout.write(output + '\n')
828
829
    if args.un_nmr:
830
        s = RNAStructure(args.file)
831
        str(s.un_nmr())
832
833
    if args.is_nmr:
834
        struc = RNAStructure(args.file)
835
        output = str(struc.is_nmr(args.verbose))
836
        sys.stdout.write(output + '\n')
837
838
    #edit
839
    if args.edit:
840
        if list != type(args.file):
841
            args.file = [args.file]
842
        ##################################
843
        for f in args.file:
844
            if args.verbose:
845
                print(f)
846
            if args.inplace:
847
                shutil.copy(f, f + '~')
848
849
            output = edit_pdb(f, args)
850
851
            if args.inplace:
852
                with open(f, 'w') as f:
853
                    f.write(output)
854
            else:  # write: to stdout
855
                try:
856
                    sys.stdout.write(output)
857
                    sys.stdout.flush()
858
                except IOError:
859
                    pass
860
861
862
    if args.fetch:
863
        fn = fetch(args.file)
864
            
865
    if args.fetch_fasta:
866
        pdb_id = args.file
867
        pdb_id = pdb_id.replace('.pdb', '')
868
869
        import urllib3
870
        http = urllib3.PoolManager()
871
872
        response = http.request('GET', 'https://www.rcsb.org/fasta/entry/' + pdb_id.upper())
873
        if not response.status == 200:
874
            raise PDBFetchError()
0 ignored issues
show
Comprehensibility Best Practice introduced by
The variable PDBFetchError does not seem to be defined.
Loading history...
875
        else:
876
            txt = response.data.decode("utf-8").strip()
877
            print(txt)            
878
            if 'No valid PDB IDs were submitted' in txt:
879
                pass
880
            else:
881
                with open(pdb_id + '.fa', 'w') as f:
882
                    f.write(txt)
883
884
        
885
    if args.fetch_chain:
886
        fn = fetch(args.file[0])
887
        chain = args.file[1]
888
        nfn = fn.replace('.pdb','') + '_' + chain + '.pdb' # 6bk8_H.pdb
889
        #cmd = 'rna_pdb_tools.py --get-chain %s %s' > %s' % (chain, fn, nfn)
890
        cmd = 'rna_pdb_tools.py --get-chain %s %s' % (chain, fn)
891
        os.system(cmd)
892
        # print(nfn)
893
894
    if args.fetch_ba:
895
        fetch_ba(args.file)
896
897
    if args.collapsed_view or args.cv:
898
        collapsed_view(args)
899
900
    if args.delete_anisou:
901
        # quick fix - make a list on the spot
902
        if list != type(args.file):
903
            args.file = [args.file]
904
        ##################################
905
        for f in args.file:
906
            if args.inplace:
907
                shutil.copy(f, f + '~')
908
909
            s = RNAStructure(f)
910
911
            output = ''
912
            if not args.no_hr:
913
                output += add_header(version) + '\n'
914
915
            for l in s.lines:
916
                if l.startswith('ANISOU'):
917
                    continue
918
                else:
919
                    output += l + '\n'
920
921
            if args.inplace:
922
                with open(f, 'w') as f:
923
                    f.write(output)
924
            else:  # write: to stdout
925
                try:
926
                    sys.stdout.write(output)
927
                    sys.stdout.flush()
928
                except IOError:
929
                    pass
930
931
932
    if args.swap_chains:
933
        # quick fix - make a list on the spot
934
        if list != type(args.file):
935
            args.file = [args.file]
936
        ##################################
937
        for f in args.file:
938
            if args.inplace:
939
                shutil.copy(f, f + '~')
940
            # rename_chain 'A>B'
941
            s = RNAStructure(f)
942
            output = ''
943
            if not args.no_hr:
944
                output += add_header(version) + '\n'
945
946
            chain_id_old, chain_id_new = args.swap_chains.split('>')
947
            if chain_id_old == chain_id_new:
948
                output = open(f).read() # return the file ;-) itself unchanged
949
            else:
950
                s.rename_chain(chain_id_new, '_')
951
                s.rename_chain(chain_id_old, chain_id_new)
952
                output += s.rename_chain('_', chain_id_old)
953
954
            if args.inplace:
955
                with open(f, 'w') as f:
956
                    f.write(output)
957
            else:  # write: to stdout
958
                try:
959
                    sys.stdout.write(output)
960
                    sys.stdout.flush()
961
                except IOError:
962
                    pass
963
964
    if args.rename_chain:
965
        # quick fix - make a list on the spot
966
        if list != type(args.file):
967
            args.file = [args.file]
968
        ##################################
969
        for f in args.file:
970
            if args.inplace:
971
                shutil.copy(f, f + '~')
972
            # rename_chain 'A>B'
973
            s = RNAStructure(f)
974
            chain_id_old, chain_id_new = args.rename_chain.split('>')
975
            output = ''
976
            if not args.no_hr:
977
                output += add_header(version) + '\n'
978
            output += s.rename_chain(chain_id_old, chain_id_new)
979
            if args.inplace:
980
                with open(f, 'w') as f:
981
                    f.write(output)
982
            else:  # write: to stdout
983
                try:
984
                    sys.stdout.write(output)
985
                    sys.stdout.flush()
986
                except IOError:
987
                    pass
988
989
990
    if args.split_alt_locations:
991
        # quick fix - make a list on the spot
992
        if list != type(args.file):
993
            args.file = [args.file]
994
        ##################################
995
        for f in args.file:
996
            if args.inplace:
997
                shutil.copy(f, f + '~')
998
999
            #s = RNAStructure(f)
1000
            s = open(f)
1001
            output = ''
1002
            if not args.no_hr:
1003
                output += add_header(version) + '\n'
1004
1005
            # First, collect all alternative locations.
1006
            alts = set()
1007
            for l in s:
1008
                if l.startswith('ATOM'):
1009
                    if l[16].strip():
1010
                        alts.add(l[16])
1011
            s.close()
1012
1013
            if args.verbose:
1014
                print('alts:', alts)
1015
1016
            for index, alt in enumerate(alts):
1017
                output += 'MODEL %i' % (index + 1)
1018
                s = open(f)
1019
                for l in s:
1020
                    if l.startswith('ATOM') or l.startswith('HETATM'):
1021
                        # if empty, then print this line
1022
                        if l[16] == ' ':
1023
                            output += l
1024
                        if l[16] == alt:
1025
                            output += l
1026
                    else:
1027
                        output += l
1028
                output += 'ENDMDL\n'
1029
                s.close()
1030
1031
            if args.inplace:
1032
                with open(f, 'w') as f:
1033
                    f.write(output)
1034
            else:  # write: to stdout
1035
                try:
1036
                    sys.stdout.write(output)
1037
                    sys.stdout.flush()
1038
                except IOError:
1039
                    pass
1040
1041
    if args.orgmode:
1042
        if args.inplace:
1043
            shutil.copy(args.file, args.file + '~')
1044
        s = RNAStructure(args.file)
1045
        s.decap_gtp()
1046
        s.std_resn()
1047
        s.remove_hydrogen()
1048
        s.remove_ion()
1049
        s.remove_water()
1050
        s.fix_op_atoms()
1051
        s.renum_atoms()
1052
        s.shift_atom_names()
1053
        s.prune_elements()
1054
        # print s.get_preview()
1055
        # s.write(args.outfile)
1056
        # for l in s.lines:
1057
        #    print l
1058
1059
        remarks = s.get_rnapuzzle_ready(
1060
            args.renumber_residues, fix_missing_atoms=True, rename_chains=True, verbose=args.verbose)
1061
1062
        with open(args.file + '~', 'w') as f:
1063
            if not args.no_hr:
1064
                f.write(add_header(version) + '\n')
1065
1066
            f.write('\n'.join(remarks) + '\n')
1067
            f.write(s.get_text())
1068
1069
        try:
1070
            from Bio import PDB
1071
            from Bio.PDB import PDBIO
1072
            import warnings
1073
            warnings.filterwarnings('ignore', '.*Invalid or missing.*',)
1074
            warnings.filterwarnings('ignore', '.*with given element *',)
1075
        except:
1076
            sys.exit('Error: Install biopython to use this function (pip biopython)')
1077
1078
        parser = PDB.PDBParser()
1079
        struct = parser.get_structure('', args.file + '~')
1080
        model = struct[0]
1081
1082
        # chains [['A', 'seq', [residues]]]
1083
        chains = []
1084
        for c in model.get_list():
1085
            seq = ''
1086
            chain = []
1087
            for r in c:
1088
                chain.append(str(r.get_resname().strip()) + str(r.id[1]))
1089
                seq += r.get_resname().strip()
1090
            chains.append([c.id, seq, chain])
1091
1092
        t = []
1093
        #[['A', 'CCGCCGCGCCAUGCCUGUGGCGG', ['C1', 'C2', 'G3', 'C4', 'C5', 'G6', 'C7', 'G8', 'C9', 'C10', 'A11', 'U12', 'G13', 'C14', 'C15', 'U16', 'G17', 'U18', 'G19', 'G20', 'C21', 'G22', 'G23']], ['B', 'CCGCCGCGCCAUGCCUGUGGCGG', ['C1', 'C2', 'G3', 'C4', 'C5', 'G6', 'C7', 'G8', 'C9', 'C10', 'A11', 'U12', 'G13', 'C14', 'C15', 'U16', 'G17', 'U18', 'G19', 'G20', 'C21', 'G22', 'G23']]]
1094
        for c in chains:
1095
            t.append('* ' + c[0] + ':' + c[2][0][1:] + '-' + c[2][-1][1:] + ' ' + c[1])
1096
            for r in c[2]:
1097
                t.append('** ' + c[0] + ':' + r)
1098
        print('\n'.join(t))
1099
1100
    if args.fix:
1101
        # quick fix - make a list on the spot
1102
        if list != type(args.file):
1103
            args.file = [args.file]
1104
        ##################################
1105
        for f in args.file:
1106
            cmd = 'pdbfixer ' + f + ' --add-atoms all --add-residues'
1107
            print(cmd)
1108
            os.system(cmd)
1109
            if args.inplace:
1110
                shutil.move("output.pdb", f)
1111
            else:
1112
                shutil.move("output.pdb", f.replace('.pdb', '_fx.pdb'))
1113
                            
1114
    if args.to_mol2:
1115
        # quick fix - make a list on the spot
1116
        if list != type(args.file):
1117
            args.file = [args.file]
1118
        ##################################
1119
        for f in args.file:
1120
            cmd = 'obabel -i pdb ' + f + ' -o mol2 -O ' + f.replace('.pdb', '.mol2')
1121
            print(cmd)
1122
            os.system(cmd)
1123
1124
1125
    if args.nmr_dir:
1126
        files = sort_strings(glob.glob(args.nmr_dir + '/' + args.file))
1127
1128
        c = 1
1129
1130
        for f in files:
1131
            #s = RNAStructure(f)
1132
            print("MODEL        " + str(c))
1133
            # at some point I could use RNAStructure for this
1134
            print(open(f).read())
1135
            #print(s.get_text(add_end=False))
1136
            print('ENDMDL')
1137
            c += 1
1138
        print('END')
1139
1140
    if args.set_chain:
1141
        if list != type(args.file):
1142
            args.file = [args.file]
1143
1144
        for f in args.file:
1145
            txt = set_chain_for_struc(f, args.set_chain)
1146
            if args.inplace:
1147
                shutil.move(f, f.replace('.pdb', '.pdb~'))
1148
                with open(f, 'w') as fi:
1149
                    fi.write(txt)
1150
            else:
1151
                print(txt)
1152
1153
    if args.remove0:
1154
        if list != type(args.file):
1155
            args.file = [args.file]
1156
1157
        for f in args.file:
1158
            s = RNAStructure(f)
1159
            txt = ''
1160
            for l in s.lines:
1161
                if s.get_atom_coords(l) == (0.0, 0.0, 0.0):
1162
                    continue
1163
                else:
1164
                    txt += l + '\n'
1165
1166
            if args.inplace:
1167
                shutil.move(f, f.replace('.pdb', '.pdb~'))
1168
                with open(f, 'w') as fi:
1169
                    fi.write(txt)
1170
            else:
1171
                print(txt)
1172
1173
    if args.replace_htm:
1174
        if list != type(args.file):
1175
            args.file = [args.file]
1176
1177
        for f in args.file:
1178
            with open(f) as fn:
1179
                txt = fn.read()
1180
            txt = txt.replace('HETATM', 'ATOM  ')
1181
            if args.inplace:
1182
                shutil.move(f, f.replace('.pdb', '.pdb~'))
1183
                with open(f, 'w') as fi:
1184
                    fi.write(txt)
1185
            else:
1186
                print(txt)
1187
1188
    if args.renum_nmr:
1189
        if list != type(args.file):
1190
            args.file = [args.file]
1191
        txt = ''
1192
        for f in args.file:
1193
            c = 1
1194
            for l in open(f):
1195
                if l.startswith('MODEL'):
1196
                    txt += "MODEL       " + str(c) + '\n'
1197
                    c += 1
1198
                elif l.strip() == 'END':
1199
                    pass
1200
                else:
1201
                    txt += l
1202
        if args.inplace:
1203
            shutil.move(f, f.replace('.pdb', '.pdb~'))
0 ignored issues
show
introduced by
The variable f does not seem to be defined for all execution paths.
Loading history...
1204
            with open(f) as fi:
1205
                fi.write(txt)
1206
        else:
1207
            print(txt)
1208
            
1209
    from rna_tools.rna_tools_config import PYMOL_PATH
1210
    sys.path.insert(0, PYMOL_PATH)
1211
    if args.cif2pdb:
1212
        # quick fix - make a list on the spot
1213
        if list != type(args.file):
1214
            args.file = [args.file]
1215
        ##################################
1216
        for cif_file in args.file:
1217
            from Bio.PDB import MMCIFParser, PDBIO
1218
            parser = MMCIFParser()
1219
            structure = parser.get_structure("structure_id", cif_file)
1220
            pdb_file = cif_file.replace('.cif', '_fCIF.pdb')
1221
1222
            remarks = []
1223
            try:
1224
                # Save to PDB format
1225
                io = PDBIO()
1226
                io.set_structure(structure)
1227
                io.save(pdb_file)
1228
            except TypeError as e:
1229
                print('Warning: some of the chains in this mmCIF file has chain names with more char than 1, e.g. AB, and the PDB format needs single-letter code, e.g. A.')
1230
1231
                def has_high_rna_content(chain, threshold=0.8):
1232
                    # RNA nucleotides: A, C, G, U, and X (you can modify as needed)
1233
                    rna_nucleotides = ['A', 'C', 'G', 'U', 'X']
1234
                    total_residues = 0
1235
                    rna_residues = 0
1236
1237
                    # Count the total number of residues and RNA-like residues
1238
                    for residue in chain:
1239
                        total_residues += 1
1240
                        if residue.get_resname().strip() in rna_nucleotides:
1241
                            rna_residues += 1
1242
1243
                    # Calculate the proportion of RNA residues
1244
                    if total_residues == 0:
1245
                        return False  # Avoid division by zero if chain has no residues
1246
1247
                    rna_percentage = rna_residues / total_residues
1248
1249
                    # Check if the percentage of RNA residues is greater than or equal to the threshold (80% by default)
1250
                    return rna_percentage >= threshold
1251
1252
                from Bio.PDB.MMCIFParser import MMCIFParser
1253
                from Bio.PDB import MMCIFParser, Structure, Model, Chain
1254
                
1255
                # Initialize the parser
1256
                parser = MMCIFParser()
1257
1258
                # Parse the structure
1259
                structure = parser.get_structure("structure", cif_file)
1260
1261
                # Create a list of single-letter chain identifiers
1262
                import string
1263
                letters = list(string.ascii_uppercase)
1264
1265
                # New structure
1266
                new_structure = Structure.Structure("new_structure")
1267
                new_model = Model.Model(0)  # Create a new model
1268
                new_structure.add(new_model)  # Add the new model to the new structure
1269
1270
                for model in structure:
1271
                    for chain in model:
1272
                        if has_high_rna_content(chain):
1273
                            chain_id_new = letters.pop(0)
1274
                            chain_id = chain.get_id()
1275
1276
                            atom_count = 0
1277
                            for residue in chain:
1278
                                  for atom in residue:
1279
                                       atom_count += 1
1280
                            print(f'rna chain {chain.id} -> {chain_id_new} # of atoms: {atom_count}')
1281
                            remarks.append(f'REMARK rna chain {chain.id} -> {chain_id_new}')
1282
                            chain.id = chain_id_new
1283
                            new_model.add(chain)
1284
1285
                io = PDBIO()
1286
                io.set_structure(new_structure)
1287
                io.save(pdb_file)
1288
1289
            print(f'saved: {pdb_file}')
1290
1291
            # open a file add remarks
1292
            new_file = ''
1293
            with open(pdb_file, 'r') as f:
1294
                if not args.no_hr:
1295
                    new_file += add_header(version) + '\n'
1296
                if remarks:
1297
                    new_file += '\n'.join(remarks) + '\n'
1298
                new_file += f.read()
1299
1300
            with open(pdb_file, 'w') as f:
1301
                f.write(new_file)
1302
1303
    if args.pdb2cif:
1304
        try:
1305
            from pymol import cmd
1306
        except ImportError:
1307
            print('This functionality needs PyMOL. Install it or if installed, check your setup')
1308
            sys.exit(1)
1309
1310
        # quick fix - make a list on the spot
1311
        if list != type(args.file):
1312
            args.file = [args.file]
1313
        ##################################
1314
        for f in args.file:
1315
            cmd.load(f)
1316
            fo = f.replace('.pdb', '.cif')
1317
            cmd.save(fo, '(all)')
1318
            cmd.delete('all')
1319
            print(fo, 'saved')
1320