Rosetta
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
pepspec.cc
Go to the documentation of this file.
1 // -*- mode:c++;tab-width:2;indent-tabs-mode:t;show-trailing-whitespace:t;rm-trailing-spaces:t -*-
2 // vi: set ts=2 noet:
3 //
4 // (c) Copyright Rosetta Commons Member Institutions.
5 // (c) This file is part of the Rosetta software suite and is made available under license.
6 // (c) The Rosetta software is developed by the contributing members of the Rosetta Commons.
7 // (c) For more information, see http://www.rosettacommons.org. Questions about this can be
8 // (c) addressed to University of Washington UW TechTransfer, email: license@u.washington.edu.
9 
10 /// @file
11 /// @brief
12 
13 
14 // libRosetta headers
15 #include <protocols/frags/VallData.hh>
16 #include <protocols/frags/TorsionFragment.hh>
17 
18 #include <core/fragment/ConstantLengthFragSet.hh>
19 #include <core/fragment/FragSet.hh>
20 #include <core/fragment/Frame.hh>
21 #include <core/fragment/picking_old/vall/util.hh>
22 #include <core/fragment/picking_old/FragmentLibraryManager.hh>
23 
24 #include <core/scoring/LREnergyContainer.hh>
25 #include <core/scoring/methods/Methods.hh>
26 #include <core/scoring/ScoreType.hh>
27 
28 #include <protocols/simple_moves/BackboneMover.hh>
29 #include <protocols/simple_moves/MinMover.hh>
30 #include <protocols/moves/MonteCarlo.hh>
31 #include <protocols/moves/Mover.hh>
32 #include <protocols/moves/MoverContainer.hh>
33 #include <protocols/moves/OutputMovers.hh>
34 #include <protocols/rigid/RigidBodyMover.hh>
35 // #include <protocols/moves/rigid_body_moves.hh>
36 #include <protocols/moves/TrialMover.hh>
37 #include <protocols/simple_moves/PackRotamersMover.hh>
38 #include <protocols/simple_moves/RotamerTrialsMover.hh>
39 #include <protocols/moves/RepeatMover.hh>
40 #include <protocols/simple_moves/FragmentMover.hh>
41 
42 #include <protocols/viewer/viewers.hh>
43 
44 #include <core/types.hh>
45 
46 #include <core/scoring/sasa.hh>
47 
48 // #include <core/util/prof.hh> // profiling
49 // #include <core/util/CacheableData.hh> // profiling
50 
51 // #include <core/sequence/SequenceMapping.hh>
52 
53 #include <core/chemical/AtomTypeSet.hh>
54 
55 #include <core/chemical/AA.hh>
56 #include <core/conformation/Residue.hh>
57 #include <core/conformation/ResidueMatcher.hh>
58 #include <core/pack/rotamer_set/RotamerCouplings.hh>
59 #include <core/chemical/ResidueTypeSet.hh>
60 #include <core/conformation/ResidueFactory.hh>
61 #include <core/chemical/VariantType.hh>
62 #include <core/chemical/util.hh>
63 #include <core/chemical/ChemicalManager.hh>
64 
65 #include <core/scoring/rms_util.hh>
66 #include <core/scoring/etable/Etable.hh>
67 #include <core/scoring/ScoringManager.hh>
68 #include <core/scoring/ScoreFunction.hh>
69 #include <core/scoring/ScoreFunctionFactory.hh>
70 #include <core/scoring/Ramachandran.hh>
71 #include <core/scoring/hbonds/HBondSet.hh>
72 #include <core/scoring/hbonds/hbonds.hh>
73 #include <core/scoring/etable/count_pair/CountPairFactory.hh>
74 #include <core/scoring/etable/count_pair/CountPairAll.hh>
75 #include <core/scoring/etable/count_pair/CountPairFunction.hh>
76 #include <core/scoring/etable/EtableEnergy.hh>
77 #include <core/scoring/etable/Etable.hh>
78 #include <core/scoring/etable/BaseEtableEnergy.tmpl.hh>
79 
80 #include <core/pack/rotamer_trials.hh>
81 #include <core/pack/pack_rotamers.hh>
82 #include <core/pack/task/PackerTask.hh>
83 #include <core/pack/task/TaskFactory.hh>
84 #include <core/pack/task/operation/TaskOperations.hh>
85 #include <protocols/toolbox/IGEdgeReweighters.hh>
86 #include <core/pack/task/IGEdgeReweightContainer.hh>
87 
88 #include <core/kinematics/FoldTree.hh>
89 #include <core/kinematics/MoveMap.hh>
90 #include <core/kinematics/util.hh>
91 
92 #include <core/pose/Pose.hh>
93 #include <core/pose/util.hh>
94 #include <core/pose/PDBPoseMap.hh>
95 #include <core/pose/PDBInfo.hh>
96 #include <core/pose/metrics/CalculatorFactory.hh>
97 
98 #include <basic/options/util.hh>//option.hh>
99 // #include <basic/options/after_opts.hh>
100 
101 #include <basic/basic.hh>
102 #include <basic/MetricValue.hh>
103 
104 #include <basic/database/open.hh>
105 
106 #include <protocols/init/init.hh>
107 
108 #include <utility/vector1.hh>
109 #include <utility/file/file_sys_util.hh>
110 
111 #include <numeric/xyzVector.hh>
112 #include <numeric/random/random.hh>
113 
114 
115 #include <core/scoring/constraints/ConstraintSet.hh>
116 #include <core/scoring/func/FlatHarmonicFunc.hh>
117 #include <core/scoring/constraints/AtomPairConstraint.hh>
118 #include <core/scoring/constraints/CoordinateConstraint.hh>
119 #include <core/scoring/constraints/ConstraintIO.hh>
120 
121 // // C++ headers
122 #include <cstdlib>
123 #include <fstream>
124 #include <iostream>
125 #include <string>
126 
127 #include <basic/Tracer.hh>
128 
129 #include <basic/options/keys/in.OptionKeys.gen.hh>
130 #include <basic/options/keys/out.OptionKeys.gen.hh>
131 #include <basic/options/keys/score.OptionKeys.gen.hh>
132 #include <basic/options/keys/pepspec.OptionKeys.gen.hh>
133 #include <basic/options/keys/constraints.OptionKeys.gen.hh>
134 
135 #include <protocols/toolbox/pose_metric_calculators/DecomposeAndReweightEnergiesCalculator.hh>
136 #include <protocols/toolbox/pose_metric_calculators/ResidueDecompositionByChainCalculator.hh>
137 #include <core/pose/metrics/simple_calculators/InterfaceNeighborDefinitionCalculator.hh>
138 #include <core/pose/metrics/simple_calculators/InterfaceSasaDefinitionCalculator.hh>
139 #include <core/pose/metrics/simple_calculators/InterfaceDeltaEnergeticsCalculator.hh>
140 
141 //Auto Headers
142 #include <core/import_pose/import_pose.hh>
143 #include <core/io/pdb/pose_io.hh>
144 #include <core/util/SwitchResidueTypeSet.hh>
145 #include <utility/vector1.hh>
146 
147 #include <utility/excn/Exceptions.hh>
148 
149 
150 using namespace core;
151 using namespace protocols;
152 using namespace ObjexxFCL;
153 //using namespace protocols;
154 
155 using utility::vector1;
156 using std::string;
157 using io::pdb::dump_pdb;
158 
159 using namespace basic::options;
160 using namespace basic::options::OptionKeys;
161 
162 using namespace pose;
163 using namespace chemical;
164 using namespace conformation;
165 using namespace scoring;
166 using namespace optimization;
167 using namespace kinematics;
168 using namespace protocols::moves;
169 using namespace id;
170 using namespace protocols::frags;
171 
172 static THREAD_LOCAL basic::Tracer TR( "apps.public.pepspec" );
173 
174 
175 Size prot_chain( 0 );
176 Size prot_begin( 0 );
177 Size prot_anchor( 0 );
178 Size prot_end( 0 );
179 
180 Size pep_jump( 2 );
181 Size pep_chain( 0 );
182 Size pep_begin( 0 );
183 Size pep_anchor( 0 );
184 Size pep_end( 0 );
185 
186 std::string input_seq;
187 
188 //load homol coord cst data
189 //parse cst line as vector1[ 10 ] name1 res1 x y z x0 sd tol
190 struct pep_coord_cst {
191  std::string atom_name;
192  int pep_pos;
193  Real x;
194  Real y;
195  Real z;
196  Real x0;
197  Real sd;
198  Real tol;
199 };
200 vector1< pep_coord_cst > pep_coord_csts;
201 
202 class myMC{
203 public:
204  myMC(
205  Pose pose,
206  Real score,
207  Real temp_init
208  );
209  Pose last_pose_;
210  Real last_score_;
211  Pose low_pose_;
212  Real low_score_;
213  Real temp_;
214  bool accept_;
215  bool accept();
216  void roll(
217  Pose & pose,
218  Real & score
219  );
220  void recover_low(
221  Pose & pose
222  );
223 };
224 
225 myMC::myMC(
226  Pose pose,
227  Real score,
228  Real temp_init
229 )
230 {
231  last_pose_ = pose;
232  last_score_ = score;
233  low_pose_ = pose;
234  low_score_ = score;
235  temp_ = temp_init;
236  accept_ = true;
237 }
238 
239 void
240 myMC::roll(
241  Pose & pose,
242  Real & score
243 )
244 {
245  accept_ = false;
246  if ( score <= last_score_ ) {
247  accept_ = true;
248  if ( score <= low_score_ ) {
249  low_pose_ = pose;
250  low_score_ = score;
251  }
252  } else {
253  Real dE = score - last_score_;
254  Real p = std::exp( -1 * dE / temp_ );
255  if ( numeric::random::rg().uniform() < p ) accept_ = true;
256  }
257  if ( accept_ ) {
258  last_pose_ = pose;
259  last_score_ = score;
260  } else {
261  pose = last_pose_;
262  score = last_score_;
263  }
264 }
265 
266 void
267 myMC::recover_low(
268  Pose & pose
269 )
270 {
271  pose = low_pose_;
272 }
273 
274 bool
275 myMC::accept()
276 {
277  return accept_;
278 }
279 
280 
281 Size
282 aa2index(
283  chemical::AA aa
284 )
285 {
286  Size index( 0 );
287  switch( aa ){
288  case chemical::aa_ala :
289  index = 1;
290  break;
291  case chemical::aa_cys :
292  index = 2;
293  break;
294  case chemical::aa_asp :
295  index = 3;
296  break;
297  case chemical::aa_glu :
298  index = 4;
299  break;
300  case chemical::aa_phe :
301  index = 5;
302  break;
303  case chemical::aa_gly :
304  index = 6;
305  break;
306  case chemical::aa_his :
307  index = 7;
308  break;
309  case chemical::aa_ile :
310  index = 8;
311  break;
312  case chemical::aa_lys :
313  index = 9;
314  break;
315  case chemical::aa_leu :
316  index = 10;
317  break;
318  case chemical::aa_met :
319  index = 11;
320  break;
321  case chemical::aa_asn :
322  index = 12;
323  break;
324  case chemical::aa_pro :
325  index = 13;
326  break;
327  case chemical::aa_gln :
328  index = 14;
329  break;
330  case chemical::aa_arg :
331  index = 15;
332  break;
333  case chemical::aa_ser :
334  index = 16;
335  break;
336  case chemical::aa_thr :
337  index = 17;
338  break;
339  case chemical::aa_val :
340  index = 18;
341  break;
342  case chemical::aa_trp :
343  index = 19;
344  break;
345  case chemical::aa_tyr :
346  index = 20;
347  break;
348  default :
349  break;
350  }
351  return index;
352 }
353 
354 Size
355 get_n_pep_nbrs(
356  pose::Pose const & pose,
357  vector1< bool > const is_pep,
358  Real const cutoff_cg
359 )
360 {
361  Size n_pep_nbrs( 0 );
362  for ( Size i=1; i<= pose.total_residue(); ++i ) {
363  if ( !is_pep[i] ) continue;
364  bool cg_res_has_nbr( false );
365  Residue const & rsd1( pose.residue(i) );
366  for ( Size j=1; j<= pose.total_residue(); ++j ) {
367  if ( cg_res_has_nbr ) break;
368  if ( is_pep[j] ) continue;
369  Residue const & rsd2( pose.residue(j) );
370  for ( Size ii=1; ii<= rsd1.natoms(); ++ii ) {
371  if ( cg_res_has_nbr ) break;
372  for ( Size jj=1; jj<= rsd2.natoms(); ++jj ) {
373  if ( rsd1.xyz(ii).distance( rsd2.xyz(jj) ) < cutoff_cg ) {
374  cg_res_has_nbr = true;
375  break;
376  }
377  }
378  }
379  }
380  if ( cg_res_has_nbr ) ++n_pep_nbrs;
381  }
382  return n_pep_nbrs;
383 }
384 
385 
386 /*
387 void
388 dump_efactor_pdb(
389 pose::Pose & pose,
390 core::scoring::ScoreFunctionOP const scorefxn,
391 std::string const & tag
392 )
393 {
394 Size const nres( pose.total_residue() );
395 // id::AtomID_Mask const & mask;
396 // id::initialize( mask, pose );
397 
398 if( option[ pepspec::homol_csts ].user() ) scorefxn->set_weight( coordinate_constraint, option[ OptionKeys::constraints::cst_weight ] );
399 ( *scorefxn )( pose );
400 
401 char const *filename = tag.c_str();
402 std::fstream out( filename, std::ios::out );
403 
404 
405 Size number(0);
406 
407 static std::string const chains( " ABCDEFGHIJKLMNOPQRSTUVWXYZ1234567890" );
408 
409 out << "MODEL " << tag << "\n";
410 for ( Size i=1; i<= pose.total_residue(); ++i ) {
411 conformation::Residue const & rsd( pose.residue(i) );
412 Real residue_total_energy( pose.energies().residue_total_energies( i ).dot( scorefxn->weights() ) );
413 for ( Size j=1; j<= rsd.natoms(); ++j ) {
414 conformation::Atom const & atom( rsd.atom(j) );
415 
416 // if ( ! mask[ id::AtomID( j,i ) ] ) continue;
417 
418 //skip outputing virtual atom unless specified
419 // if ( !options::option[ options::OptionKeys::out::file::output_virtual ]() &&
420 // rsd.is_virtual(j) ) continue;
421 
422 ++number;
423 assert( rsd.chain() < int(chains.size()) ); // silly restriction
424 char const chain( chains[ rsd.chain() ] );
425 out << "ATOM " << I(5,number) << ' ' << rsd.atom_name(j) << ' ' <<
426 rsd.name3() << ' ' << chain << I(4,rsd.seqpos() ) << " " <<
427 F(8,3,atom.xyz()(1)) <<
428 F(8,3,atom.xyz()(2)) <<
429 F(8,3,atom.xyz()(3)) <<
430 F(6,2,1.0) << F(6,2, residue_total_energy ) << '\n';
431 }
432 }
433 out << "ENDMDL\n\n\n";
434 out << "###residue_total_energies###\n";
435 for ( Size i=1; i<= pose.total_residue(); ++i ) {
436 Real residue_total_energy( pose.energies().residue_total_energies( i ).dot( scorefxn->weights() ) );
437 out << string_of( i ) << " " << pose.residue( i ).name1() << " " << pose.energies().residue_total_energies( i ).weighted_string_of( scorefxn->weights() ) << "\t" << residue_total_energy << "\n";
438 }
439 }
440 pose::Pose
441 my_boltzmann(
442 pose::Pose pose,
443 Real score
444 )
445 {
446 
447 }
448 void
449 gen_fold_tree_for_nbr_segments(
450 pose::Pose & pose,
451 FoldTree & ftree,
452 vector1< bool > const & is_ligand,
453 vector1< bool > const & is_skipped,
454 Real const & nbr_cutoff,
455 vector1< bool > & is_nbr
456 )
457 {
458 Size nres( pose.total_residue() );
459 //define neighbors, all protein residues within cutoff from ligand, excluding is_skipped
460 set_ss_from_phipsi( pose );
461 for( Size i=1; i<= pose.total_residue(); ++i ) {
462 if ( !is_ligand[ i ] ) continue;
463 Residue const & rsd1( pose.residue(i) );
464 for ( Size j=1; j<= pose.total_residue(); ++j ) {
465 Residue const & rsd2( pose.residue(j) );
466 if( is_ligand[ j ] || is_skipped[ j ] || !rsd2.is_protein() || !( pose.secstruct( j ) == 'L' ) ) continue;
467 for ( Size ii=1; ii<= rsd1.natoms(); ++ii ) {
468 if( is_nbr[ j ] ) break;
469 for ( Size jj=1; jj<= rsd2.natoms(); ++jj ) {
470 if ( rsd1.xyz(ii).distance_squared( rsd2.xyz(jj) ) < nbr_cutoff*nbr_cutoff ) {
471 is_nbr[ j ] = true;
472 break;
473 }
474 }
475 }
476 }
477 }
478 
479 for( Size i = 2, j = 1; i <= pose.total_residue(); ++i, ++j ){
480 if( is_nbr[ i ] && !is_nbr[ i-1 ] ) j = 1;
481 if( is_nbr[ i ] && !is_nbr[ i+1 ] ){
482 ftree.new_jump( i - j, i - static_cast< int >( j / 2 ), i - j );
483 ftree.new_jump( i - j, i + 1, i );
484 }
485 }
486 
487 TR << "sel flex_prot, ";
488 for( Size i = 1; i <= pose.total_residue(); ++i ){
489 if( is_nbr[ i ] ) TR << "resi " << string_of( i ) << " or ";
490 }
491 TR << std::endl;
492 }
493 */
494 
495 /// @details This function will make a sequence mutation while trying to preserve the variants
496 void
497 make_sequence_change(
498  Size const seqpos,
499  AA const & new_aa,
500  pose::Pose & pose
501 )
502 {
503  conformation::Residue const & current_rsd( pose.residue( seqpos ) );
504  if ( current_rsd.aa() == new_aa ) return; // already done
505 
506  ResidueTypeCOP rsd_type( current_rsd.residue_type_set()->get_representative_type_aa( new_aa, current_rsd.type().variant_types() ) );
507 
508  if ( ! rsd_type ) TR << "make_sequence_change failed: new_aa= "+name_from_aa(new_aa)+" -- no matching residue types" << std::endl;
509  else {
510  conformation::ResidueOP new_rsd( ResidueFactory::create_residue( *rsd_type, current_rsd, pose.conformation() ) );
511  pose.replace_residue( seqpos, *new_rsd, false );
512  }
513 }
514 
515 
516 bool
517 has_clash(
518  pose::Pose pose,
519  vector1< bool > is_checked,
520  scoring::ScoreFunctionOP const & scorefxn,
521  Real const clash_threshold
522 )
523 {
524  using namespace scoring;
525  using namespace chemical;
526 
527  bool is_clash( false );
528  for ( Size seqpos = 1; seqpos <= pose.total_residue(); ++seqpos ) {
529  if ( !is_checked[ seqpos ] ) continue;
530 
531  ( *scorefxn )( pose );
532 
533 
534  // cached energies object
535  Energies & energies( pose.energies() );
536 
537  // the neighbor/energy links
538  EnergyGraph & energy_graph( energies.energy_graph() );
539 
540  // search upstream
541  for ( graph::Graph::EdgeListIter
542  iru = energy_graph.get_node( seqpos )->lower_edge_list_begin(),
543  irue = energy_graph.get_node( seqpos )->lower_edge_list_end();
544  iru != irue; ++iru ) {
545  EnergyEdge * edge( static_cast< EnergyEdge *> (*iru) );
546  // Size const j( edge->get_first_node_ind() );
547 
548  // the pair energies cached in the link
549  EnergyMap const & emap( edge->fill_energy_map());
550  Real const clash( emap[ fa_rep ] );
551  if ( clash > clash_threshold ) {
552  //TR<< "fa_rep: " << string_of( clash ) << " at " << pose.residue( j ).name1() << string_of( j ) << "-" << pose.residue( seqpos ).name1() << string_of( seqpos ) << std::endl;
553  is_clash = true;
554  break;
555  }
556  }
557  if ( is_clash == true ) break;
558 
559  // and downstream
560  for ( graph::Graph::EdgeListIter
561  iru = energy_graph.get_node( seqpos )->upper_edge_list_begin(),
562  irue = energy_graph.get_node( seqpos )->upper_edge_list_end();
563  iru != irue; ++iru ) {
564  EnergyEdge * edge( static_cast< EnergyEdge *> (*iru) );
565  // Size const j( edge->get_second_node_ind() );
566 
567  // the pair energies cached in the link
568  EnergyMap const & emap( edge->fill_energy_map());
569  Real const clash( emap[ fa_rep ] );
570  if ( clash > clash_threshold ) {
571  //TR<< "fa_rep: " << string_of( clash ) << " at " << pose.residue( j ).name1() << string_of( j ) << "-" << pose.residue( seqpos ).name1() << string_of( seqpos ) << std::endl;
572  is_clash = true;
573  break;
574  }
575  }
576  if ( is_clash == true ) break;
577  }
578  return is_clash;
579 }
580 
581 vector1< std::pair< Size, Size > >
582 get_clash_pairs(
583  pose::Pose pose,
584  vector1< bool > is_checked,
585  scoring::ScoreFunctionOP const & scorefxn,
586  Real const clash_threshold
587 )
588 {
589  using namespace scoring;
590  using namespace chemical;
591  vector1< std::pair< Size, Size > > clash_pairs;
592 
593  for ( Size seqpos = 1; seqpos <= pose.total_residue(); ++seqpos ) {
594  if ( !is_checked[ seqpos ] ) continue;
595  ( *scorefxn )( pose );
596  // cached energies object
597  Energies & energies( pose.energies() );
598  // the neighbor/energy links
599  EnergyGraph & energy_graph( energies.energy_graph() );
600  // search upstream
601  for ( graph::Graph::EdgeListIter
602  iru = energy_graph.get_node( seqpos )->lower_edge_list_begin(),
603  irue = energy_graph.get_node( seqpos )->lower_edge_list_end();
604  iru != irue; ++iru ) {
605  EnergyEdge * edge( static_cast< EnergyEdge *> (*iru) );
606  Size const j( edge->get_first_node_ind() );
607 
608  // the pair energies cached in the link
609  EnergyMap const & emap( edge->fill_energy_map());
610  Real const clash( emap[ fa_rep ] );
611  if ( clash > clash_threshold ) {
612  clash_pairs.push_back( std::pair< Size, Size >( seqpos, j ) );
613  }
614  }
615  // and downstream
616  for ( graph::Graph::EdgeListIter
617  iru = energy_graph.get_node( seqpos )->upper_edge_list_begin(),
618  irue = energy_graph.get_node( seqpos )->upper_edge_list_end();
619  iru != irue; ++iru ) {
620  EnergyEdge * edge( static_cast< EnergyEdge *> (*iru) );
621  Size const j( edge->get_second_node_ind() );
622 
623  // the pair energies cached in the link
624  EnergyMap const & emap( edge->fill_energy_map());
625  Real const clash( emap[ fa_rep ] );
626  if ( clash > clash_threshold ) {
627  clash_pairs.push_back( std::pair< Size, Size >( seqpos, j ) );
628  }
629  }
630  }
631  return clash_pairs;
632 }
633 
634 std::string
635 pep_rmsd_analysis(
636  pose::Pose pose
637 )
638 {
639 
640  pose::Pose ref_pose;
641  std::string ref_name( option[ in::file::native ]() );
642  import_pose::pose_from_pdb( ref_pose, ref_name );
643 
644  Size ref_pep_anchor_in( option[ pepspec::pep_anchor ] );
645  if ( option[ pepspec::native_pep_anchor ].user() ) ref_pep_anchor_in = option[ pepspec::native_pep_anchor ];
646  std::string ref_pep_chain_in( option[ pepspec::pep_chain ] );
647  if ( option[ pepspec::native_pep_chain ].user() ) ref_pep_chain_in = option[ pepspec::native_pep_chain ];
648  Size ref_pep_anchor( ref_pose.pdb_info()->pdb2pose( ref_pep_chain_in[0], ref_pep_anchor_in ) );
649  Size ref_pep_chain( ref_pose.chain( ref_pep_anchor ) );
650  Size ref_pep_begin( ref_pose.conformation().chain_begin( ref_pep_chain ) );
651  Size ref_pep_end( ref_pose.conformation().chain_end( ref_pep_chain ) );
652  Size this_pep_begin( pep_begin );
653  Size this_pep_end( pep_end );
654 
655  Size ref_prot_chain( 1 );
656  for ( Size i = 1; i <= ref_pose.conformation().num_chains(); ++i ) {
657  if ( i == ref_pep_chain ) continue;
658  if ( !( ref_pose.residue( ref_pose.conformation().chain_begin( i ) ).is_protein() ) ) continue;
659  else {
660  ref_prot_chain = i;
661  break;
662  }
663  }
664  Size ref_prot_begin( ref_pose.conformation().chain_begin( ref_prot_chain ) );
665 
666  Size pep_nterm( pep_anchor - this_pep_begin );
667  Size ref_pep_nterm( ref_pep_anchor - ref_pep_begin );
668  Size pep_cterm( this_pep_end - pep_anchor );
669  Size ref_pep_cterm( ref_pep_end - ref_pep_anchor );
670 
671  Size nterm( pep_nterm );
672  Size cterm( pep_cterm );
673  if ( pep_nterm < ref_pep_nterm ) ref_pep_begin += ref_pep_nterm - pep_nterm;
674  else if ( pep_nterm > ref_pep_nterm ) {
675  this_pep_begin += pep_nterm - ref_pep_nterm;
676  nterm = ref_pep_nterm;
677  }
678  if ( pep_cterm < ref_pep_cterm ) ref_pep_end -= ref_pep_cterm - pep_cterm;
679  else if ( pep_cterm > ref_pep_cterm ) {
680  this_pep_end -= pep_cterm - ref_pep_cterm;
681  cterm = ref_pep_cterm;
682  }
683  //superpose if needed
684  if ( option[ pepspec::native_align ] ) {
685  id::AtomID_Map< id::AtomID > atom_map;
686  pose::initialize_atomid_map( atom_map, pose, id::BOGUS_ATOM_ID );
687 
688  for ( Size i = prot_begin; i <= prot_end; ++i ) {
689  id::AtomID const id1( pose.residue( i ).atom_index( "CA" ), i );
690  id::AtomID const id2( ref_pose.residue( i + static_cast< int >( ref_prot_begin ) - static_cast< int >( prot_begin ) ).atom_index( "CA" ), i + static_cast< int >( ref_prot_begin ) - static_cast< int >( prot_begin ) );
691  atom_map[ id1 ] = id2;
692  }
693  core::scoring::ScoreFunctionOP full_scorefxn( get_score_function() );
694  core::scoring::superimpose_pose( pose, ref_pose, atom_map );
695 
696  }
697  Real total_rmsd( 0 );
698  std::string rmsd_analysis( "" );
699  for ( Size i_seq = 0; i_seq <= nterm + cterm; ++i_seq ) {
700  Size ref_pep_seqpos = ref_pep_begin + i_seq;
701  Size pep_seqpos = this_pep_begin + i_seq;
702  Real sd( ref_pose.residue( ref_pep_seqpos ).xyz( "CA" ).distance_squared( pose.residue( pep_seqpos ).xyz( "CA" ) ) );
703  total_rmsd += sd;
704  rmsd_analysis += "p" + string_of( i_seq - static_cast< int >( nterm ) ) + "_rmsd:\t" + string_of( std::sqrt( sd ) ) + "\t";
705  }
706  total_rmsd = std::sqrt( total_rmsd / ( nterm+cterm+1 ) );
707  rmsd_analysis += "total_rmsd:\t" + string_of( total_rmsd ) + "\t";
708  return rmsd_analysis;
709 
710 }
711 
712 std::string
713 pep_phipsi_analysis(
714  pose::Pose pose
715 )
716 {
717 
718  pose::Pose ref_pose;
719  std::string ref_name( option[ in::file::native ]() );
720  import_pose::pose_from_pdb( ref_pose, ref_name );
721 
722  Size ref_pep_anchor_in( option[ pepspec::pep_anchor ] );
723  if ( option[ pepspec::native_pep_anchor ].user() ) ref_pep_anchor_in = option[ pepspec::native_pep_anchor ];
724  std::string ref_pep_chain_in( option[ pepspec::pep_chain ] );
725  if ( option[ pepspec::native_pep_chain ].user() ) ref_pep_chain_in = option[ pepspec::native_pep_chain ];
726  Size ref_pep_anchor( ref_pose.pdb_info()->pdb2pose( ref_pep_chain_in[0], ref_pep_anchor_in ) );
727  Size ref_pep_chain( ref_pose.chain( ref_pep_anchor ) );
728  Size ref_pep_begin( ref_pose.conformation().chain_begin( ref_pep_chain ) );
729  Size ref_pep_end( ref_pose.conformation().chain_end( ref_pep_chain ) );
730  Size this_pep_begin( pep_begin );
731  Size this_pep_end( pep_end );
732 
733  Size pep_nterm( pep_anchor - this_pep_begin );
734  Size ref_pep_nterm( ref_pep_anchor - ref_pep_begin );
735  Size pep_cterm( this_pep_end - pep_anchor );
736  Size ref_pep_cterm( ref_pep_end - ref_pep_anchor );
737 
738  Size nterm( pep_nterm );
739  Size cterm( pep_cterm );
740  if ( pep_nterm < ref_pep_nterm ) ref_pep_begin += ref_pep_nterm - pep_nterm;
741  else if ( pep_nterm > ref_pep_nterm ) {
742  this_pep_begin += pep_nterm - ref_pep_nterm;
743  nterm = ref_pep_nterm;
744  }
745  if ( pep_cterm < ref_pep_cterm ) ref_pep_end -= ref_pep_cterm - pep_cterm;
746  else if ( pep_cterm > ref_pep_cterm ) {
747  this_pep_end -= pep_cterm - ref_pep_cterm;
748  cterm = ref_pep_cterm;
749  }
750 
751  std::string phipsi_analysis( "" );
752  for ( Size i_seq = 0; i_seq <= nterm + cterm; ++i_seq ) {
753  Size ref_pep_seqpos = ref_pep_begin + i_seq;
754  Size pep_seqpos = this_pep_begin + i_seq;
755  Real ramadev( 0 );
756 
757  //delta phi, psi, omega
758  Real phi( std::abs( ref_pose.phi( ref_pep_seqpos ) - pose.phi( pep_seqpos ) ) );
759  if ( phi > 180 ) phi = std::abs( 360 -phi );
760  ramadev += ( phi * phi );
761  phipsi_analysis += "p" + string_of( i_seq - static_cast< int >( nterm ) ) + "_phi:\t" + string_of( phi ) + "\t";
762 
763  Real psi( std::abs( ref_pose.psi( ref_pep_seqpos ) - pose.psi( pep_seqpos ) ) );
764  if ( psi > 180 ) psi = std::abs( 360 - psi );
765  ramadev += ( psi * psi );
766  phipsi_analysis += "p" + string_of( i_seq - static_cast< int >( nterm ) ) + "_psi:\t" + string_of( psi ) + "\t";
767 
768  Real omega( std::abs( ref_pose.omega( ref_pep_seqpos ) - pose.omega( pep_seqpos ) ) );
769  if ( omega > 180 ) omega = std::abs( 360 - omega );
770  phipsi_analysis += "p" + string_of( i_seq - static_cast< int >( nterm ) ) + "_omega:\t" + string_of( omega ) + "\t";
771 
772  phipsi_analysis += "p" + string_of( i_seq - static_cast< int >( nterm ) ) + "_ramadev:\t" + string_of( std::sqrt( ramadev / 2 ) ) + "\t";
773  }
774  return phipsi_analysis;
775 
776 }
777 
778 void add_termini(
779  pose::Pose & pose
780 )
781 {
782  add_lower_terminus_type_to_pose_residue( pose, pep_begin );
783  add_upper_terminus_type_to_pose_residue( pose, pep_end );
784  pose.conformation().update_polymeric_connection( pep_begin );
785  pose.conformation().update_polymeric_connection( pep_begin + 1 );
786  pose.conformation().update_polymeric_connection( pep_end );
787  pose.conformation().update_polymeric_connection( pep_end - 1 );
788 
789 }
790 
791 void add_pep_res(
792  pose::Pose & pose,
793  bool add_nterm,
794  bool add_cterm
795 )
796 {
797 
798  ResidueTypeSet const & rsd_set( *pose.residue(1).residue_type_set() );
799  ResidueOP vrt( ResidueFactory::create_residue( rsd_set.name_map( "GLY" ) ) );
800  // ResidueOP vrt( ResidueFactory::create_residue( rsd_set.name_map( "VirtBB" ) ) );
801 
802  if ( add_cterm ) {
803  pose.conformation().safely_append_polymer_residue_after_seqpos( *vrt, pep_end, false );
804  pep_end = pep_end + 1;
805  pose.set_omega( pep_end - 1, 180.0 );
806  pose.conformation().update_polymeric_connection( pep_end );
807  pose.conformation().update_polymeric_connection( pep_end - 1 );
808  add_variant_type_to_pose_residue( pose, core::chemical::VIRTUAL_BB, pep_end );
809  }
810 
811  if ( add_nterm ) {
812  pose.conformation().safely_prepend_polymer_residue_before_seqpos( *vrt, pep_begin, false );
813  pep_end = pep_end + 1;
814  pep_anchor = pep_anchor + 1;
815  pose.set_omega( pep_begin, 180.0 );
816  pose.conformation().update_polymeric_connection( pep_begin );
817  pose.conformation().update_polymeric_connection( pep_begin + 1 );
818  add_variant_type_to_pose_residue( pose, core::chemical::VIRTUAL_BB, pep_begin );
819  }
820 
821  //replace termini
822  add_termini( pose );
823 }
824 
825 void remove_pep_res(
826  pose::Pose & pose,
827  bool add_nterm,
828  bool add_cterm
829 )
830 {
831  if ( add_cterm ) {
832  pose.conformation().delete_residue_slow( pep_end );
833  pep_end = pep_end - 1;
834  pose.conformation().update_polymeric_connection( pep_end );
835  pose.conformation().update_polymeric_connection( pep_end - 1 );
836  }
837 
838  if ( add_nterm ) {
839  pose.conformation().delete_residue_slow( pep_begin );
840  pep_anchor = pep_anchor - 1;
841  pep_end = pep_end - 1;
842  pose.conformation().update_polymeric_connection( pep_begin );
843  pose.conformation().update_polymeric_connection( pep_begin + 1 );
844  }
845 
846  //replace termini
847  add_termini( pose );
848 }
849 
850 void
851 initialize_peptide(
852  pose::Pose & pose
853 )
854 {
855 
856  if ( option[ pepspec::remove_input_bb ] ) {
857  //remove non-anchor positions
858  while ( pep_anchor != pep_begin ) {
859  pose.conformation().delete_residue_slow( pep_begin );
860  pep_anchor = pep_anchor - 1;
861  pep_end = pep_end - 1;
862  pose.conformation().update_polymeric_connection( pep_begin );
863  pose.conformation().update_polymeric_connection( pep_begin + 1 );
864  }
865  while ( pep_anchor != pep_end ) {
866  pose.conformation().delete_residue_slow( pep_end );
867  pep_end = pep_end - 1;
868  pose.conformation().update_polymeric_connection( pep_end );
869  pose.conformation().update_polymeric_connection( pep_end - 1 );
870  }
871  }
872 
873  Pose start_pose( pose );
874  ResidueTypeSet const & rsd_set( *pose.residue(1).residue_type_set() );
875  //remove pep_anchor termini
876  if ( pep_begin == pep_anchor && pep_end == pep_anchor ) {
877  std::string pep_anchor_type( pose.residue( pep_anchor ).name3() );
878  Residue rsd1( pose.residue( prot_anchor ) );
879  Residue rsd2( pose.residue( pep_anchor ) );
880  StubID upstubid( AtomID( rsd1.atom_index( "N" ), prot_anchor ) ,
881  AtomID( rsd1.atom_index( "CA" ), prot_anchor ) ,
882  AtomID( rsd1.atom_index( "C" ), prot_anchor ) ) ;
883 
884  StubID downstubid( AtomID( rsd2.atom_index( "N" ), pep_anchor ) ,
885  AtomID( rsd2.atom_index( "CA" ), pep_anchor ) ,
886  AtomID( rsd2.atom_index( "C" ), pep_anchor ) ) ;
887  RT const rt( pose.conformation().get_stub_transform( upstubid, downstubid ) );
888 
889  pose = *pose.split_by_chain( prot_chain );
890 
891  ResidueOP pep_anchor_res_ptr( ResidueFactory::create_residue( pose.residue( 1 ).residue_type_set()->name_map( pep_anchor_type ) ) );
892  pose.append_residue_by_jump( *pep_anchor_res_ptr, prot_anchor, "", "", true );
893  pose.conformation().set_stub_transform( upstubid, downstubid, rt );
894  for ( Size i_chi = 1; i_chi <= pose.residue( pep_anchor ).nchi(); ++i_chi ) {
895  pose.set_chi( i_chi, pep_anchor, start_pose.residue( pep_anchor ).chi( i_chi ) );
896  }
897  }
898  Residue pep_anchor_res = pose.residue( pep_anchor );
899 
900  //must maintain fold tree!
901  pose.fold_tree( start_pose.fold_tree() );
902 
903  //if using input backbone and designing, mutate all res to cg_res_type
904  std::string cg_res_type( option[ pepspec::cg_res_type ] );
905  ResidueOP res( ResidueFactory::create_residue( rsd_set.name_map( cg_res_type ) ) );
906  if ( option[ pepspec::use_input_bb ] && !option[ pepspec::no_design ] && !option[ pepspec::input_seq ].user() ) {
907  for ( Size mut_site = pep_begin; mut_site <= pep_end; mut_site++ ) {
908  if ( mut_site==pep_anchor ) continue;
909  make_sequence_change( mut_site, chemical::AA( aa_from_oneletter_code( res->name1() ) ) , pose );
910  }
911  }
912 
913 }
914 
915 void
916 set_pep_csts(
917  pose::Pose & pose
918 )
919 {
920  using namespace core::scoring::constraints;
921  for ( Size i_cst = 1; i_cst <= pep_coord_csts.size(); ++i_cst ) {
922  pep_coord_cst pep_cst( pep_coord_csts[ i_cst ] );
923  Size seqpos( pep_cst.pep_pos + pep_anchor );
924  if ( seqpos < pep_begin ) continue;
925  if ( seqpos > pep_end ) break;
926  if ( numeric::random::rg().uniform() > option[ pepspec::p_homol_csts ] ) continue;
927  Vector pep_cst_vector;
928  pep_cst_vector.x( pep_cst.x );
929  pep_cst_vector.y( pep_cst.y );
930  pep_cst_vector.z( pep_cst.z );
931  ConstraintCOP this_cst( ConstraintOP( new CoordinateConstraint( AtomID( pose.residue( seqpos ).atom_index( pep_cst.atom_name ), seqpos ), AtomID( pose.residue( prot_anchor ).atom_index( "CA" ), prot_anchor ), pep_cst_vector, core::scoring::func::FuncOP( new core::scoring::func::FlatHarmonicFunc( pep_cst.x0, pep_cst.sd, pep_cst.tol ) ) ) ) );
932  // ConstraintCOP cst( new CoordinateConstraint( AtomID( pose.residue( i ).atom_index( "CA" ), i ), AtomID( pose.residue( pep_anchor ).atom_index( "CA" ), pep_anchor ), pose.residue( i ).xyz( "CA" ), new FlatHarmonicFunc( 0.0, 0.1, 2.0 ) ) );
933  pose.add_constraint( this_cst );
934  }
935 }
936 
937 /*
938 void
939 refine_fa_pep_bb(
940 Pose & pose,
941 vector1< bool > is_pep,
942 scoring::ScoreFunctionOP cen_scorefxn
943 )
944 {
945 
946 kinematics::MoveMapOP mm ( new kinematics::MoveMap );
947 mm->set_bb( is_pep );
948 mm->set_chi( pep_anchor, true );
949 
950 rigid::RigidBodyPerturbMoverOP rb_mover = new rigid::RigidBodyPerturbMover( pep_jump, 0.1, 0.0 );
951 rb_mover->apply( pose );
952 
953 protocols::simple_moves::SmallMoverOP rep_small_mover( new protocols::simple_moves::SmallMover( mm, 2.0, 10 ) );
954 rep_small_mover->angle_max( 'H', 1.5 );
955 rep_small_mover->angle_max( 'E', 1.5 );
956 rep_small_mover->angle_max( 'L', 1.5 );
957 
958 protocols::simple_moves::MinMoverOP min_mover = new protocols::simple_moves::MinMover( mm, cen_scorefxn, "dfpmin", 0.001, true );
959 
960 RandomMoverOP rand_mover( new protocols::moves::RandomMover() );
961 rand_mover->add_mover( rep_small_mover, 9 );
962 // rand_mover->add_mover( rb_mover, 1 );
963 // rand_mover->add_mover( min_mover, 1 );
964 
965 MonteCarloOP mc_rep ( new MonteCarlo( pose, *cen_scorefxn, 0.8 ) );
966 TrialMoverOP rep_trial = new TrialMover( rand_mover, mc_rep );
967 for( Size i = 1; i <= 1000; ++i ){
968 rep_trial->apply( pose );
969 mc_rep->boltzmann( pose );
970 }
971 // mc_rep->recover_low( pose );
972 }
973 */
974 
975 /// @brief helper code for fragments generation, copied from S.M.Lewis
976 core::fragment::FragSetCOP
977 make_frags(
978  core::Size const start,
979  core::Size const stop,
980  std::string const & seq
981 )
982 {
983  core::Size const frags_length( 1 ); //magic number: 3mer fragments!!
984  core::fragment::FragSetOP fragset( new core::fragment::ConstantLengthFragSet( frags_length ) );
985  std::string ss_string( frags_length, 'L' );
986  core::fragment::FragDataOPs list;
987 
988  for ( core::Size j = start; j <= stop - frags_length + 1; ++j ) {
989  std::string const seqsubstr( seq, j-1, frags_length ); //j-1 accounts for string [] from 0
990  list = core::fragment::picking_old::vall::pick_fragments_by_ss_plus_aa( ss_string, seqsubstr, 200, false ); //magic number: 200 fragments per position (not duplicated - this will be like robetta server fragments)
991  core::fragment::FrameOP frame;
992  frame = core::fragment::FrameOP( new core::fragment::Frame( j ) );
993  frame->add_fragment( list );
994  fragset->add( frame );
995  }
996  return fragset;
997 }
998 
999 //@brief helper code for fragments generation, copied from S.M.Lewis
1000 core::fragment::FragSetCOP
1001 make_1mer_frags(
1002  core::Size const seqpos_start,
1003  core::Size const seqpos_stop,
1004  std::string const & seq,
1005  Size const nfrags
1006 )
1007 {
1008  core::Size const frags_length( 1 );
1009  core::fragment::FragSetOP fragset( new core::fragment::ConstantLengthFragSet( frags_length ) );
1010  //80% L, 10% H, 10% S
1011  core::fragment::FragDataOPs list;
1012 
1013  for ( core::Size j = seqpos_start; j <= seqpos_stop - frags_length + 1; ++j ) {
1014  //get substr
1015  //get ss str
1016  std::string ss_string( "L" );
1017  if ( option[ pepspec::ss_type ].user() ) ss_string = std::string( option[ pepspec::ss_type ] );
1018  else {
1019  if ( numeric::random::rg().uniform() < 0.1 ) ss_string = "H";
1020  else if ( numeric::random::rg().uniform() > 0.9 ) ss_string = "S";
1021  }
1022  //get fraglist from ss or ss and seq
1023  if ( option[ pepspec::input_seq ].user() && seq[ j - 1 ] != 'X' ) {
1024  std::string const seqsubstr( seq, j-1, frags_length ); //j-1 accounts for string [] from 0
1025  list = core::fragment::picking_old::vall::pick_fragments_by_ss_plus_aa( ss_string, seqsubstr, nfrags, true );
1026  } else list = core::fragment::picking_old::vall::pick_fragments_by_ss( ss_string, nfrags, true );
1027  core::fragment::FrameOP frame;
1028  frame = core::fragment::FrameOP( new core::fragment::Frame( j ) );
1029  frame->add_fragment( list );
1030  fragset->add( frame );
1031  }
1032  return fragset;
1033 }
1034 
1035 void
1036 gen_pep_bb_sequential(
1037  pose::Pose & pose,
1038  scoring::ScoreFunctionOP cen_scorefxn
1039 )
1040 {
1041  using namespace scoring;
1042  using namespace chemical;
1043  using namespace scoring::methods;
1044  using namespace scoring::etable;
1045  using namespace scoring::etable::count_pair;
1046  using core::pack::task::operation::TaskOperationCOP;
1047 
1048  Size n_build_loop( option[ pepspec::n_build_loop ] );
1049  Size n_prepend( option[ pepspec::n_prepend ] );
1050  Size n_append( option[ pepspec::n_append ] );
1051 
1052  ( *cen_scorefxn )( pose );
1053  Real min_score( pose.energies().total_energies()[ vdw ] );
1054  Real clash_cutoff( option[ pepspec::clash_cutoff ] );
1055  clash_cutoff = std::max( clash_cutoff, pose.energies().total_energies()[ fa_rep ] );
1056 
1057  Pose replace_res_pose( pose );
1058  ResidueTypeSet const & rsd_set( *pose.residue(1).residue_type_set() );
1059  std::string cg_res_type( option[ pepspec::cg_res_type ] );
1060  ResidueOP res( ResidueFactory::create_residue( rsd_set.name_map( cg_res_type ) ) );
1061  Size break_loop( 0 );
1062  Size n_prepended( 0 );
1063  Size n_appended( 0 );
1064  //Size iqq( 1 ); // unused ~Labonte
1065  while ( ( n_prepended < n_prepend || n_appended < n_append ) ) {
1066  Size this_prepend( n_prepend );
1067  Size this_append( n_append );
1068  if ( option[ pepspec::gen_pep_bb_sequential ] || !option[ pepspec::homol_csts ].user() ) {
1069  this_prepend = static_cast< int >( numeric::random::rg().uniform() * ( n_prepend - n_prepended + 1 ) );
1070  this_append = static_cast< int >( numeric::random::rg().uniform() * ( n_append - n_appended + 1 ) );
1071  }
1072  if ( this_prepend + this_append == 0 ) continue;
1073  //std::cout << "Prepending/Appending " + string_of( this_prepend ) + " / " + string_of( this_append ) << std::endl;
1074 
1075  std::string this_input_seq( input_seq );
1076  for ( Size ii = 1; ii <= this_prepend; ++ii ) {
1077  // pose.prepend_polymer_residue_before_seqpos( *res, pep_begin, false );
1078  pose.prepend_polymer_residue_before_seqpos( *res, pep_begin, true );
1079  pep_end = pep_end + 1;
1080  pep_anchor = pep_anchor + 1;
1081  pose.set_omega( pep_begin, 180.0 );
1082  pose.conformation().update_polymeric_connection( pep_begin );
1083  pose.conformation().update_polymeric_connection( pep_begin + 1 );
1084  }
1085  for ( Size ii = 1; ii <= this_append; ++ii ) {
1086  // pose.append_polymer_residue_after_seqpos( *res, pep_end, false );
1087  pose.append_polymer_residue_after_seqpos( *res, pep_end, true );
1088  pep_end = pep_end + 1;
1089  pose.set_omega( pep_end - 1, 180.0 );
1090  pose.conformation().update_polymeric_connection( pep_end );
1091  pose.conformation().update_polymeric_connection( pep_end - 1 );
1092  }
1093  //use input seq?
1094  this_input_seq.erase( n_prepended + this_prepend + 1, n_append - n_appended - this_append );
1095  this_input_seq.erase( 0, n_prepend - n_prepended - this_prepend );
1096  if ( option[ pepspec::input_seq ].user() ) {
1097  for ( Size ii = pep_begin; ii <= pep_end; ++ii ) {
1098  Size this_input_seq_pos( ii - pep_begin );
1099  char this_input_seq_aa( this_input_seq[ this_input_seq_pos ] );
1100  if ( this_input_seq_aa != 'X' ) {
1101  make_sequence_change( ii, chemical::AA( aa_from_oneletter_code( this_input_seq_aa ) ), pose );
1102  }
1103  }
1104  }
1105 
1106  pose.update_residue_neighbors();
1107 
1108  //turn on homol_csts
1109  if ( option[ pepspec::homol_csts ].user() ) {
1110  cen_scorefxn->set_weight( coordinate_constraint, option[ OptionKeys::constraints::cst_weight ] );
1111  set_pep_csts( pose );
1112  }
1113 
1114  //setup pep vectors
1115  Size nres_pep( pep_end - pep_begin + 1 );
1116  vector1< bool > is_pep( pose.total_residue(), false );
1117  for ( Size i = 1; i <= pose.total_residue(); ++i ) is_pep[ i ] = ( i >= pep_begin && i <= pep_end );
1118  // vector1< bool > is_insert( nres_pep, false );
1119  kinematics::MoveMapOP mm_frag( new kinematics::MoveMap );
1120  for ( Size ii = 1; ii <= this_prepend + 1; ++ii ) {
1121  Size moveable_seqpos( pep_begin + ii - 1 );
1122  mm_frag->set_bb( moveable_seqpos, true );
1123  // is_insert[ ii ] = true;
1124  }
1125  for ( Size ii = nres_pep; ii >= nres_pep - this_append; --ii ) {
1126  Size moveable_seqpos( pep_begin + ii - 1 );
1127  mm_frag->set_bb( moveable_seqpos, true );
1128  // is_insert[ ii ] = true;
1129  }
1130 
1131  ( *cen_scorefxn )( pose );
1132  Real best_score( pose.energies().total_energies()[ total_score ] );
1133  MonteCarloOP mc_frag( new MonteCarlo( pose, *cen_scorefxn, 2.0 ) );
1134 
1135  //replace prot w/ A/G/P?
1136  if ( option[ pepspec::no_cen_rottrials ] ) {
1137  for ( Size seqpos = prot_begin; seqpos <= prot_end; ++seqpos ) {
1138  if ( pose.residue( seqpos ).name1() != 'G' && pose.residue( seqpos ).name1() != 'P' ) {
1139  make_sequence_change( seqpos, chemical::aa_ala, pose );
1140  }
1141  }
1142  }
1143 
1144  //phi/psi insertion loop
1145  if ( !option[ pepspec::use_input_bb ] && n_build_loop > 0 ) {
1146  // Size this_build_loop( static_cast< int >( n_build_loop * std::max( this_prepend, this_append ) / std::max( n_prepend, n_append ) ) );
1147  //do at least 10 frags/phipsi
1148  Size this_build_loop( std::max(
1149  Size( 10 ),
1150  Size( n_build_loop * std::max( this_prepend, this_append ) / std::max( n_prepend, n_append ) ) ) );
1151 
1152  //get frags
1153  Size const nfrags( this_build_loop );
1154  core::fragment::FragSetCOP fragset( make_1mer_frags( pep_begin, pep_end, this_input_seq, nfrags ) );
1155  protocols::simple_moves::ClassicFragmentMover frag_mover( fragset, mm_frag );
1156 
1157  for ( Size build_loop_inner = 1; build_loop_inner <= this_build_loop; ++build_loop_inner ) {
1158  /*
1159  //choose an insertion position = seqpos
1160  Size insert_peppos = 0;
1161  while( true ){
1162  insert_peppos = static_cast< int > ( nres_pep * numeric::random::rg().uniform() ) + 1;
1163  if( is_insert[ insert_peppos ] ) break;
1164  }
1165  Size insert_seqpos( pep_begin + insert_peppos - 1 );
1166  */
1167  frag_mover.apply( pose );
1168  //random perturb angles +/- 1 degree
1169  // pose.set_phi( insert_seqpos, pose.phi( insert_seqpos ) + 1 * ( 2 * numeric::random::rg().uniform() - 1 ) );
1170  // pose.set_psi( insert_seqpos, pose.psi( insert_seqpos ) + 1 * ( 2 * numeric::random::rg().uniform() - 1 ) );
1171 
1172  Real test_score( pose.energies().total_energies().dot( cen_scorefxn->weights() ) );
1173  if ( mc_frag->boltzmann( pose ) ) {
1174  //need to do my own eval for <= because MC mover is < only
1175  if ( test_score <= best_score ) {
1176  best_score = test_score;
1177  mc_frag->reset( pose );
1178  }
1179  }
1180  }
1181  mc_frag->recover_low( pose );
1182  }
1183  ( *cen_scorefxn )( pose );
1184 
1185  bool this_has_clash( false );
1186 
1187  //check no_cen clash
1188  //break if has clash already
1189  if ( has_clash( pose, is_pep, cen_scorefxn, clash_cutoff ) ) {
1190  this_has_clash = true;
1191  //replace prot residues
1192  for ( Size seqpos = prot_begin; seqpos <= prot_end; ++seqpos ) {
1193  if ( pose.residue( seqpos ).name1() != 'G' && pose.residue( seqpos ).name1() != 'P' ) {
1194  pose.replace_residue( seqpos, replace_res_pose.residue( seqpos ), true );
1195  }
1196  }
1197  } else if ( option[ pepspec::no_cen_rottrials ] ) {
1198  //else replace rotamers and repack if necessary
1199  //ignore preexisting clashes
1200  vector1< bool > ignore_clash( pose.total_residue(), false );
1201  for ( Size i = 1; i <= pose.total_residue() - 1; ++i ) {
1202  vector1< bool > this_clash( pose.total_residue(), false );
1203  this_clash[ i ] = true;
1204  if ( has_clash( pose, this_clash, cen_scorefxn, clash_cutoff ) ) {
1205  ignore_clash[ i ] = true;
1206  }
1207  }
1208  //replace prot residues
1209  for ( Size seqpos = prot_begin; seqpos <= prot_end; ++seqpos ) {
1210  pose.replace_residue( seqpos, replace_res_pose.residue( seqpos ), true );
1211  }
1212  //find clashes and repack
1213  if ( has_clash( pose, is_pep, cen_scorefxn, clash_cutoff ) ) {
1214  vector1< bool > repack_this( pose.total_residue(), false );
1215  vector1< std::pair< Size, Size > > clash_pairs( get_clash_pairs( pose, is_pep, cen_scorefxn, clash_cutoff ) );
1216  //get all prot clash res and their nbrs
1217  for ( Size i = 1; i <= clash_pairs.size(); ++i ) {
1218  std::pair< Size, Size > clash_pair( clash_pairs[ i ] );
1219  Size this_clash( clash_pair.second );
1220  EnergyGraph const & energy_graph( pose.energies().energy_graph() );
1221  for ( graph::Graph::EdgeListConstIter
1222  ir = energy_graph.get_node( this_clash )->const_edge_list_begin(),
1223  ire = energy_graph.get_node( this_clash )->const_edge_list_end();
1224  ir != ire; ++ir ) {
1225  Size this_nbr( (*ir)->get_other_ind( this_clash ) );
1226  repack_this[ this_nbr ] = true;
1227  }
1228  }
1229  //and repack
1230  pack::task::TaskFactoryOP rp_task_factory( new pack::task::TaskFactory );
1231  pack::task::operation::RestrictResidueToRepackingOP restrict_to_repack_taskop( new pack::task::operation::RestrictResidueToRepacking() );
1232  pack::task::operation::PreventRepackingOP prevent_repack_taskop( new pack::task::operation::PreventRepacking() );
1233  for ( Size i=1; i<= pose.total_residue(); ++i ) {
1234  if ( repack_this[ i ] ) restrict_to_repack_taskop->include_residue( i );
1235  else prevent_repack_taskop->include_residue( i );
1236  }
1237  rp_task_factory->push_back( TaskOperationCOP( new pack::task::operation::IncludeCurrent() ) );
1238  rp_task_factory->push_back( restrict_to_repack_taskop );
1239  rp_task_factory->push_back( prevent_repack_taskop );
1240  protocols::simple_moves::RotamerTrialsMoverOP dz_rottrial( new protocols::simple_moves::RotamerTrialsMover( cen_scorefxn, rp_task_factory ) );
1241  dz_rottrial->apply( pose );
1242  //check clashes one last time
1243  vector1< bool > check_clash( pose.total_residue(), false );
1244  for ( Size i = 1; i <= pose.total_residue() - 1; ++i ) if ( is_pep[ i ] || ( repack_this[ i ] && !ignore_clash[ i ] ) ) check_clash[ i ] = true;
1245 
1246  this_has_clash = has_clash( pose, check_clash, cen_scorefxn, clash_cutoff );
1247  }
1248  } else {
1249  //check cen clash
1250  Real clash_score( pose.energies().total_energies()[ vdw ] );
1251  this_has_clash = ( clash_score > min_score + 0.03 );
1252  }
1253 
1254  //if gen clash, step back a step
1255  if ( break_loop < 5 && this_has_clash ) {
1256  ++break_loop;
1257  //std::cout << "Removing " + string_of( this_prepend ) + " / " + string_of( this_append ) << std::endl;
1258  for ( Size ii = 1; ii <= this_prepend; ++ii ) {
1259  pose.conformation().delete_residue_slow( pep_begin );
1260  pep_anchor = pep_anchor - 1;
1261  pep_end = pep_end - 1;
1262  pose.conformation().update_polymeric_connection( pep_begin );
1263  pose.conformation().update_polymeric_connection( pep_begin + 1 );
1264  }
1265  for ( Size ii = 1; ii <= this_append; ++ii ) {
1266  pose.conformation().delete_residue_slow( pep_end );
1267  pep_end = pep_end - 1;
1268  pose.conformation().update_polymeric_connection( pep_end );
1269  pose.conformation().update_polymeric_connection( pep_end - 1 );
1270  }
1271  } else {
1272  n_prepended += this_prepend;
1273  n_appended += this_append;
1274  }
1275  }
1276 }
1277 
1278 
1279 void
1280 perturb_pep_bb(
1281  pose::Pose & pose,
1282  kinematics::MoveMapOP mm_move,
1283  scoring::ScoreFunctionOP cen_scorefxn,
1284  Size n_iter
1285 )
1286 {
1287  ( *cen_scorefxn )( pose );
1288  protocols::simple_moves::SmallMoverOP cg_small( new protocols::simple_moves::SmallMover( mm_move, 5.0, 1 ) );
1289  cg_small->angle_max( 'H', 10.0 );
1290  cg_small->angle_max( 'E', 10.0 );
1291  cg_small->angle_max( 'L', 10.0 );
1292 
1293  MonteCarloOP mc_cg( new MonteCarlo( pose, *cen_scorefxn, 2.0 ) );
1294  Real best_score( pose.energies().total_energies().dot( cen_scorefxn->weights() ) );
1295  for ( Size cgmove_loop = 1; cgmove_loop <= n_iter; cgmove_loop++ ) {
1296  cg_small->apply( pose );
1297  if ( mc_cg->boltzmann( pose ) ) {
1298  Real test_score( pose.energies().total_energies().dot( cen_scorefxn->weights() ) );
1299  if ( test_score <= best_score ) {
1300  best_score = test_score;
1301  mc_cg->reset( pose );
1302  }
1303  }
1304  }
1305  mc_cg->recover_low( pose );
1306 }
1307 
1308 void
1309 mutate_random_residue(
1310  Pose & pose,
1311  vector1< bool > is_mutable,
1312  ScoreFunctionOP soft_scorefxn,
1313  ScoreFunctionOP full_scorefxn
1314 )
1315 {
1316  //pick rand seqpos
1317  assert( is_mutable.size() == pose.total_residue() );
1318  bool mutate( false );
1319  Size seqpos( 0 );
1320  while ( !mutate ) {
1321  seqpos = static_cast< int >( numeric::random::rg().uniform() * is_mutable.size() + 1 );
1322  mutate = is_mutable[ seqpos ];
1323  }
1324 
1325  //pick rand aa and mutate
1326  make_sequence_change( seqpos, chemical::AA( static_cast< int > ( 20 * numeric::random::rg().uniform() + 1 ) ), pose );
1327 
1328  //rottrial seqpos only
1329  using core::pack::task::operation::TaskOperationCOP;
1330  pack::task::TaskFactoryOP mut_task_factory( new pack::task::TaskFactory );
1331  {
1332  pack::task::operation::RestrictResidueToRepackingOP restrict_to_repack_taskop( new pack::task::operation::RestrictResidueToRepacking() );
1333  pack::task::operation::PreventRepackingOP prevent_repack_taskop( new pack::task::operation::PreventRepacking() );
1334  for ( Size i=1; i<= pose.total_residue(); ++i ) {
1335  if ( i == seqpos ) restrict_to_repack_taskop->include_residue( i );
1336  else prevent_repack_taskop->include_residue( i );
1337  }
1338  mut_task_factory->push_back( TaskOperationCOP( new pack::task::operation::InitializeFromCommandline() ) );
1339  mut_task_factory->push_back( restrict_to_repack_taskop );
1340  mut_task_factory->push_back( prevent_repack_taskop );
1341  }
1342  protocols::simple_moves::RotamerTrialsMoverOP mut_rottrial( new protocols::simple_moves::RotamerTrialsMover( soft_scorefxn, mut_task_factory ) );
1343  mut_rottrial->apply( pose );
1344 
1345  //get seqpos nbrs from energy map
1346  vector1< bool > is_nbr( pose.total_residue(), false );
1347  EnergyGraph const & energy_graph( pose.energies().energy_graph() );
1348  for ( graph::Graph::EdgeListConstIter
1349  ir = energy_graph.get_node( seqpos )->const_edge_list_begin(),
1350  ire = energy_graph.get_node( seqpos )->const_edge_list_end();
1351  ir != ire; ++ir ) {
1352  Size i_nbr( (*ir)->get_other_ind( seqpos ) );
1353  is_nbr[ i_nbr ] = true && ( i_nbr != pep_anchor );
1354  }
1355 
1356  //rottrial seqpos nbrs too
1357  pack::task::TaskFactoryOP rp_task_factory( new pack::task::TaskFactory );
1358  {
1359  pack::task::operation::RestrictResidueToRepackingOP restrict_to_repack_taskop( new pack::task::operation::RestrictResidueToRepacking() );
1360  pack::task::operation::PreventRepackingOP prevent_repack_taskop( new pack::task::operation::PreventRepacking() );
1361  for ( Size i=1; i<= pose.total_residue(); ++i ) {
1362  if ( ( i == seqpos || is_nbr[ i ] ) && i != pep_anchor ) restrict_to_repack_taskop->include_residue( i );
1363  else prevent_repack_taskop->include_residue( i );
1364  }
1365  rp_task_factory->push_back( TaskOperationCOP( new pack::task::operation::InitializeFromCommandline() ) );
1366  rp_task_factory->push_back( TaskOperationCOP( new pack::task::operation::IncludeCurrent() ) );
1367  rp_task_factory->push_back( restrict_to_repack_taskop );
1368  rp_task_factory->push_back( prevent_repack_taskop );
1369  }
1370  protocols::simple_moves::RotamerTrialsMoverOP rp_rottrial( new protocols::simple_moves::RotamerTrialsMover( soft_scorefxn, rp_task_factory ) );
1371  rp_rottrial->apply( pose );
1372 
1373  /*
1374  if( pose.residue( seqpos ).name3() != "PRO" && pose.residue( seqpos ).name3() != "GLY" && pose.residue( seqpos ).name3() != "ALA" ){
1375  MonteCarloOP mc( new MonteCarlo( pose, *soft_scorefxn, 1.0 ) );
1376  protocols::simple_moves::sidechain_moves::SidechainMoverOP sidechain_mover( new protocols::simple_moves::sidechain_moves::SidechainMover() );
1377  sidechain_mover->set_task_factory( rp_task_factory );
1378  sidechain_mover->set_prob_uniform( 0.05 );
1379  for( Size ii = 1; ii <= 20; ++ii ){
1380  sidechain_mover->apply( pose );
1381  mc->boltzmann( pose );
1382  }
1383  mc->recover_low( pose );
1384  }
1385  */
1386 
1387  kinematics::MoveMapOP mm_min( new kinematics::MoveMap );
1388  mm_min->set_chi( seqpos );
1389  mm_min->set_chi( is_nbr );
1390  protocols::simple_moves::MinMoverOP min_mover( new protocols::simple_moves::MinMover( mm_min, full_scorefxn, "dfpmin", 0.001, true ) );
1391  min_mover->apply( pose );
1392 }
1393 
1394 void
1395 packmin_unbound_pep(
1396  pose::Pose & pose,
1397  scoring::ScoreFunctionOP full_scorefxn
1398 )
1399 {
1400  pose.update_residue_neighbors();
1401  ( *full_scorefxn )( pose );
1402  pack::task::PackerTaskOP task( pack::task::TaskFactory::create_packer_task( pose ));
1403  task->initialize_from_command_line().or_include_current( true );
1404  task->restrict_to_repacking();
1405  protocols::simple_moves::PackRotamersMoverOP pack( new protocols::simple_moves::PackRotamersMover( full_scorefxn, task, 1 ) );
1406  pack->apply( pose );
1407  kinematics::MoveMapOP mm( new kinematics::MoveMap );
1408  mm->set_chi( true );
1409  protocols::simple_moves::MinMoverOP min_mover( new protocols::simple_moves::MinMover( mm, full_scorefxn, "dfpmin", 0.001, true ) );
1410  min_mover->apply( pose );
1411 }
1412 
1413 Real
1414 get_binding_score(
1415  Pose pose,
1416  Size pep_chain,
1417  ScoreFunctionOP full_scorefxn
1418 )
1419 {
1420  ( *full_scorefxn )( pose );
1421  Real total_score( pose.energies().total_energies().dot( full_scorefxn->weights() ) );
1422  Pose pep_pose( *pose.split_by_chain( pep_chain ) );
1423  packmin_unbound_pep( pep_pose, full_scorefxn );
1424  ( *full_scorefxn )( pep_pose );
1425  Real pep_score = pep_pose.energies().total_energies().dot( full_scorefxn->weights() );
1426  Real bind_score = total_score - pep_score;
1427  return bind_score;
1428 }
1429 
1430 
1431 void
1432 print_pep_analysis(
1433  std::string pdb_name,
1434  std::fstream & out_file,
1435  Pose pose,
1436  Real prot_score,
1437  ScoreFunctionOP full_scorefxn,
1438  bool dump_pdb
1439 )
1440 {
1441  //must remove constraints!
1442  ( *full_scorefxn )( pose );
1443  std::string output_seq;
1444  for ( Size i = pep_begin; i <= pep_end; i++ ) {
1445  output_seq.append( 1, pose.residue( i ).name1() );
1446  }
1447  Real total_score( pose.energies().total_energies().dot( full_scorefxn->weights() ) );
1448  out_file << pdb_name + "\t" << output_seq << "\t" << pose.energies().total_energies().weighted_string_of( full_scorefxn->weights() );
1449  out_file<<"\ttotal_score:\t"<<total_score<<"\t";
1450  out_file<<"total-prot_score:\t"<<total_score - prot_score<<"\t";
1451  if ( option[ pepspec::homol_csts ].user() ) {
1452  full_scorefxn->set_weight( coordinate_constraint, option[ OptionKeys::constraints::cst_weight ] );
1453  ( *full_scorefxn )( pose );
1454  out_file << "coordinate_constraint:\t" << pose.energies().total_energies()[ coordinate_constraint ] << "\t";
1455  out_file << "total_cst_score:\t" << pose.energies().total_energies().dot( full_scorefxn->weights() ) << "\t";
1456  // pose.constraint_set( 0 );
1457  full_scorefxn->set_weight( coordinate_constraint, 0 );
1458  ( *full_scorefxn )( pose );
1459  }
1460  //score pep unbound state//
1461  if ( option[ pepspec::binding_score ] ) {
1462  /*
1463  Pose pep_pose( pose.split_by_chain( pep_chain ) );
1464  packmin_unbound_pep( pep_pose, full_scorefxn );
1465  ( *full_scorefxn )( pep_pose );
1466  Real pep_score = pep_pose.energies().total_energies().dot( full_scorefxn->weights() );
1467  Real bind_score = total_score - pep_score;
1468  */
1469  Real bind_score( get_binding_score( pose, pep_chain, full_scorefxn ) );
1470  out_file<<"binding_score:\t"<<bind_score<<"\t";
1471  out_file<<"binding-prot_score:\t"<<bind_score - prot_score<<"\t";
1472  // TR << pdb_name + ".pep\t" << output_seq << "\t" << pep_pose.energies().total_energies().weighted_string_of( full_scorefxn->weights() ) + "\ttotal_score:\t" << pep_score << "\n";
1473  }
1474  if ( option[ pepspec::upweight_interface ] ) {
1475  basic::MetricValue< Real > weighted_interface_score;
1476  pose.metric("energy_decomposition", "weighted_total", weighted_interface_score);
1477  out_file << "reweighted_score:\t" + string_of( weighted_interface_score.value() ) + "\t";
1478  }
1479  basic::MetricValue< Real > interface_score;
1480  pose.metric( "interface_delta_energies", "weighted_total", interface_score );
1481  out_file << "interface_score:\t" + string_of( interface_score.value() ) + "\t";
1482  if ( option[ pepspec::calc_sasa ] ) {
1483  basic::MetricValue< Real > delta_sasa;
1484  pose.metric( "sasa_interface", "delta_sasa", delta_sasa );
1485  out_file << "interface_sasa:\t" + string_of( delta_sasa.value() ) + "\t";
1486  }
1487  // if( option[ pepspec::rmsd_analysis ] ) out_file << pep_rmsd_analysis( pose );
1488  // if( option[ pepspec::phipsi_analysis ] ) out_file << pep_phipsi_analysis( pose );
1489  out_file<<std::endl;
1490  if ( dump_pdb ) pose.dump_scored_pdb( pdb_name, *full_scorefxn );
1491 
1492 }
1493 
1494 void
1495 RunPepSpec()
1496 {
1497  //load loop lengths//
1498  Size n_peptides( option[ pepspec::n_peptides ] );
1499  Size n_cgrelax_loop( option[ pepspec::n_cgrelax_loop ] );
1500 
1501  //data out
1502  std::string out_nametag( "data" );
1503  if ( option[ out::file::o ].user() ) out_nametag = option[ out::file::o ];
1504  std::string pdb_dir( out_nametag + ".pdbs" );
1505  std::string pdb_name;
1506  utility::file::create_directory( pdb_dir );
1507  std::string out_file_name_str( out_nametag + ".spec" );
1508  char const *out_file_name = out_file_name_str.c_str();
1509  std::fstream out_file( out_file_name, std::ios::out );
1510 
1511  //load input structure pdbs/csts//
1512  pose::Pose pose;
1513  vector1< std::string > pdb_filenames;
1514  if ( option[ pepspec::pdb_list ].user() ) {
1515  std::string pdb_list_filename( option[ pepspec::pdb_list ] );
1516  std::ifstream pdb_list_data( pdb_list_filename.c_str() );
1517  if ( !pdb_list_data.good() ) {
1518  utility_exit_with_message( "Unable to open file: " + pdb_list_filename + '\n' );
1519  }
1520  std::string pdb_list_line;
1521  while ( !getline( pdb_list_data, pdb_list_line, '\n' ).eof() ) {
1522  std::string this_filename( pdb_list_line );
1523  pdb_filenames.push_back( this_filename );
1524  }
1525  } else {
1526  pdb_filenames.push_back( basic::options::start_file() );
1527 
1528  }
1529 
1530  Size n_prepend( option[ pepspec::n_prepend ] );
1531  Size n_append( option[ pepspec::n_append ] );
1532 
1533  bool const save_low_pdbs( option[ pepspec::save_low_pdbs ] );
1534  bool const save_all_pdbs( option[ pepspec::save_all_pdbs ] );
1535 
1536  //load homol_csts?
1537  if ( option[ pepspec::homol_csts ].user() ) {
1538  //parse cst line as: atom_name pep_pos x y z x0 sd tol
1539  std::string cst_filename( option[ pepspec::homol_csts ] );
1540  std::ifstream cst_file( cst_filename.c_str() );
1541  std::string str_buffer;
1542  while ( !getline( cst_file, str_buffer, '\t' ).eof() ) {
1543  pep_coord_cst pep_cst;
1544 
1545  pep_cst.atom_name = str_buffer;
1546  getline( cst_file, str_buffer, '\t' );
1547  {
1548  std::istringstream istr_buffer( str_buffer );
1549  istr_buffer >> pep_cst.pep_pos;
1550  }
1551  getline( cst_file, str_buffer, '\t' );
1552  {
1553  std::istringstream istr_buffer( str_buffer );
1554  istr_buffer >> pep_cst.x;
1555  }
1556  getline( cst_file, str_buffer, '\t' );
1557  {
1558  std::istringstream istr_buffer( str_buffer );
1559  istr_buffer >> pep_cst.y;
1560  }
1561  getline( cst_file, str_buffer, '\t' );
1562  {
1563  std::istringstream istr_buffer( str_buffer );
1564  istr_buffer >> pep_cst.z;
1565  }
1566  getline( cst_file, str_buffer, '\t' );
1567  {
1568  std::istringstream istr_buffer( str_buffer );
1569  istr_buffer >> pep_cst.x0;
1570  }
1571  getline( cst_file, str_buffer, '\t' );
1572  {
1573  std::istringstream istr_buffer( str_buffer );
1574  istr_buffer >> pep_cst.sd;
1575  }
1576  getline( cst_file, str_buffer, '\n' );
1577  {
1578  std::istringstream istr_buffer( str_buffer );
1579  istr_buffer >> pep_cst.tol;
1580  }
1581  pep_coord_csts.push_back( pep_cst );
1582  }
1583  }
1584 
1585 
1586  //define scoring functions//
1587  core::scoring::ScoreFunctionOP full_scorefxn( get_score_function() );
1588  core::scoring::ScoreFunctionOP soft_scorefxn( ScoreFunctionFactory::create_score_function( option[ pepspec::soft_wts ] ) );
1589  core::scoring::ScoreFunctionOP cen_scorefxn( new core::scoring::ScoreFunction() );
1590  if ( option[ pepspec::cen_wts ].user() ) {
1591  cen_scorefxn = ScoreFunctionFactory::create_score_function( option[ pepspec::cen_wts ] );
1592  } else {
1593  cen_scorefxn->set_weight( fa_rep, 0.5 );
1594  }
1595  // Real const max_interchain_wt( cen_scorefxn->get_weight( interchain_contact ) );
1596 
1597  //////////////fragment library/////////////////
1598  /*
1599  VallLibrarian librarian;
1600  librarian.add_fragment_gen( new LengthGen( 1 ) );
1601  librarian.catalog( FragmentLibraryManager::get_instance()->get_Vall() );
1602  */
1603 
1604 
1605  // set up the pose_metric_calculators
1606  core::pose::metrics::CalculatorFactory & calculator_factory(core::pose::metrics::CalculatorFactory::Instance());
1607  protocols::toolbox::pose_metric_calculators::ResidueDecompositionByChainCalculatorOP res_decomp_calculator( new protocols::toolbox::pose_metric_calculators::ResidueDecompositionByChainCalculator() );
1608  calculator_factory.register_calculator("residue_decomposition", res_decomp_calculator);
1609  protocols::toolbox::pose_metric_calculators::DecomposeAndReweightEnergiesCalculatorOP decomp_reweight_calculator( new protocols::toolbox::pose_metric_calculators::DecomposeAndReweightEnergiesCalculator( "residue_decomposition" ) );
1610  vector1< Real > reweight_vector( 1, 3.0 ); reweight_vector.push_back( 2.0 );
1611  decomp_reweight_calculator->master_weight_vector( reweight_vector );
1612  decomp_reweight_calculator->num_sets( ( Size ) 2 );
1613  calculator_factory.register_calculator("energy_decomposition", decomp_reweight_calculator);
1614  core::pose::metrics::PoseMetricCalculatorOP int_calculator( new core::pose::metrics::simple_calculators::InterfaceNeighborDefinitionCalculator( (Size)1, (Size)2 ) );
1615  calculator_factory.register_calculator( "interface", int_calculator );
1616  core::pose::metrics::PoseMetricCalculatorOP int_delta_energy_calculator( new core::pose::metrics::simple_calculators::InterfaceDeltaEnergeticsCalculator( "interface" ) );
1617  calculator_factory.register_calculator( "interface_delta_energies", int_delta_energy_calculator );
1618  core::pose::metrics::PoseMetricCalculatorOP int_sasa_calculator( new core::pose::metrics::simple_calculators::InterfaceSasaDefinitionCalculator( (Size)1, (Size)2 ) );
1619  calculator_factory.register_calculator( "sasa_interface", int_sasa_calculator );
1620 
1621  int pose_index( 0 );
1622  if ( option[ pepspec::run_sequential ] ) n_peptides = pdb_filenames.size();
1623  for ( Size peptide_loop = 1; peptide_loop <= n_peptides; ++peptide_loop ) {
1624 
1625  /*
1626  VallData vall( option[ pepspec::frag_file ] );
1627  TorsionFragmentLibrary lib;
1628  lib.resize( 20 );
1629  //lib indexed same as AA enum
1630  std::string ss_seq( "L" );
1631  if( option[ pepspec::ss_type ].user() ) ss_seq = option[ pepspec::ss_type ];
1632  for( Size i = 1; i <= 20; ++i ){
1633  char aa( oneletter_code_from_aa( chemical::AA( i ) ) );
1634  std::string frag_seq;
1635  frag_seq.push_back( aa );
1636  Real seq_wt( 0.1 );
1637  Real ss_wt( 0.1 );
1638  vall.get_frags( 1000, frag_seq, ss_seq, seq_wt, ss_wt, false, false, true, lib[ i ] );
1639  }
1640  */
1641 
1642  //load random start pdb//
1643  if ( option[ pepspec::run_sequential ] ) ++pose_index;
1644  else pose_index = static_cast< int >( numeric::random::rg().uniform() * pdb_filenames.size() + 1 );
1645  std::string pdb_filename( pdb_filenames[ pose_index ] );
1646  TR<<"Initializing "<< out_nametag + "_" + string_of( peptide_loop ) + " with " + pdb_filename << std::endl;
1647  import_pose::pose_from_pdb( pose, pdb_filename );
1648 
1649  // ResidueTypeSet const & rsd_set( pose.residue(1).residue_type_set() );
1650  //convert user input to internal values//
1651  if ( option[ pepspec::pep_anchor ].user() ) {
1652  Size const pep_anchor_in( option[ pepspec::pep_anchor ] );
1653  std::string const pep_chain_in( option[ pepspec::pep_chain ] );
1654  pep_anchor = pose.pdb_info()->pdb2pose( pep_chain_in[ 0 ], pep_anchor_in );
1655  pep_chain = pose.chain( pep_anchor );
1656  pep_begin = pose.conformation().chain_begin( pep_chain );
1657  pep_end = pose.conformation().chain_end( pep_chain );
1658  } else {
1659  pep_chain = 2;
1660  pep_begin = pose.conformation().chain_begin( pep_chain );
1661  pep_end = pose.conformation().chain_end( pep_chain );
1662  pep_anchor = pep_begin + static_cast< int >( ( pose.conformation().chain_end( pep_chain ) - pose.conformation().chain_begin( pep_chain ) ) * numeric::random::rg().uniform() );
1663  }
1664 
1665  if ( pep_anchor == 0 ) utility_exit_with_message( "pep_chain / pep_anchor combo not found\n" );
1666  for ( Size i = 1; i <= pose.conformation().num_chains(); ++i ) {
1667  if ( i == pep_chain ) continue;
1668  if ( !( pose.residue( pose.conformation().chain_begin( i ) ).is_protein() ) ) continue;
1669  else {
1670  prot_chain = i;
1671  break;
1672  }
1673  }
1674  prot_begin = pose.conformation().chain_begin( prot_chain );
1675  prot_end =pose.conformation().chain_end( prot_chain );
1676  prot_anchor =prot_begin;
1677 
1678  Real cutoff( option[ pepspec::interface_cutoff ] );
1679  //std::string output_seq;
1680 
1681  //gen fold tree//
1682  FoldTree f( pose.total_residue() );
1683 
1684  std::string pep_anchor_root( "CB" );
1685  if ( pose.residue( pep_anchor ).name1() == 'G' ) pep_anchor_root = "CA";
1686  pep_jump = 2;
1687  if ( prot_chain < pep_chain ) {
1688  pep_jump = f.new_jump( prot_anchor, pep_anchor, pep_begin - 1 );
1689  f.set_jump_atoms( pep_jump, "CA", pep_anchor_root );
1690  if ( pep_end != pose.total_residue() ) f.new_jump( prot_anchor, pep_end + 1, pep_end );
1691  if ( prot_end + 1 != pep_begin ) f.new_jump( prot_anchor, prot_end + 1, prot_end );
1692  } else {
1693  pep_jump = f.new_jump( pep_anchor, prot_anchor, prot_begin - 1 );
1694  f.set_jump_atoms( pep_jump, pep_anchor_root, "CA" );
1695  if ( prot_end != pose.total_residue() ) f.new_jump( prot_anchor, prot_end + 1, prot_end );
1696  if ( pep_end + 1 != prot_begin ) f.new_jump( pep_anchor, pep_end + 1, pep_end );
1697  }
1698  f.reorder( prot_anchor );
1699  TR << f << std::endl;
1700 
1701  //set foldtree in the pose//
1702  pose.fold_tree( f );
1703 
1704  pose::Pose start_pose( pose );
1705  Pose prot_pose = *pose.split_by_chain( prot_chain );
1706  ( *full_scorefxn )( prot_pose );
1707  packmin_unbound_pep( prot_pose, full_scorefxn );
1708  Real prot_score( prot_pose.energies().total_energies().dot( full_scorefxn->weights() ) );
1709 
1710  protocols::viewer::add_conformation_viewer( pose.conformation(), "pepspec_pose" );
1711 
1712  //init peptide
1713  initialize_peptide( pose );
1714 
1715  //load input seq?
1716  input_seq = std::string( n_prepend + n_append + 1, 'X' );
1717  if ( option[ pepspec::input_seq ].user() ) input_seq = std::string( option[ pepspec::input_seq ] );
1718 
1719  //bool add_nterm( true );
1720  //if( pep_anchor == pep_begin ) add_nterm = false; // unused ~Labonte
1721  //bool add_cterm( true );
1722  //if( pep_anchor == pep_end ) add_cterm = false; // unused ~Labonte
1723 
1724  Residue pep_anchor_res = pose.residue( pep_anchor );
1725 
1726  //convert to CG residues//
1727  if ( !option[ pepspec::no_cen ] ) core::util::switch_to_residue_type_set( pose, core::chemical::CENTROID );
1728  ( *cen_scorefxn )( pose );
1729 
1730  Real min_vdw( pose.energies().total_energies()[ vdw ] );
1731 
1732  gen_pep_bb_sequential( pose, cen_scorefxn );
1733 
1734  if ( option[ pepspec::homol_csts ].user() ) {
1735  cen_scorefxn->set_weight( coordinate_constraint, option[ OptionKeys::constraints::cst_weight ] );
1736  set_pep_csts( pose );
1737  }
1738 
1739  ( *cen_scorefxn )( pose );
1740  if ( peptide_loop > 1 && pose.energies().total_energies()[ vdw ] > min_vdw + 0.1 ) {
1741  --peptide_loop;
1742  continue;
1743  }
1744 
1745  //define prot, pep, and anchor residues//
1746  vector1< bool > is_pep( pose.total_residue(), false ), is_mutable( pose.total_residue(), false ), is_prot( pose.total_residue(), false ), is_anchor( pose.total_residue(), false );
1747  vector1< Size > pep_res_vec;
1748  for ( Size i=1; i<= pose.total_residue(); ++i ) {
1749  is_pep[i] = ( i >= pep_begin && i <= pep_end );
1750  if ( is_pep[ i ] ) pep_res_vec.push_back( i );
1751  is_prot[i] = ( i >= prot_begin && i <= prot_end );
1752  is_anchor[i] = ( ( i == pep_anchor || i == prot_anchor ) );
1753  is_mutable[i] = ( i >= pep_begin && i <= pep_end && !(i==pep_anchor) );
1754  if ( option[ pepspec::input_seq ].user() && i >= pep_begin ) is_mutable[i] = ( input_seq[ i - pep_begin ] == 'X' && i != pep_anchor );
1755  }
1756 
1757  //make peptide from CG backbone//
1758  Pose restart_cgrelax_pose( pose );
1759  Size restart_pep_begin( pep_begin );
1760  Size restart_pep_anchor( pep_anchor );
1761  Size restart_pep_end( pep_end );
1762  for ( Size cgrelax_loop = 1; cgrelax_loop <= n_cgrelax_loop; cgrelax_loop++ ) {
1763  pose = restart_cgrelax_pose;
1764  pep_begin = restart_pep_begin;
1765  pep_anchor = restart_pep_anchor;
1766  pep_end = restart_pep_end;
1767 
1768  //min movemap//
1769  kinematics::MoveMapOP mm_min( new kinematics::MoveMap );
1770  mm_min->set_chi( is_pep );
1771 
1772  //small-small movemap//
1773  kinematics::MoveMapOP mm_move( new kinematics::MoveMap );
1774  mm_move->set_bb( is_pep );
1775 
1776  //small CG backbone moves//
1777  Size n_cgperturb_iter( 100 );
1778  if ( option[ pepspec::use_input_bb ] ) n_cgperturb_iter = option[ pepspec::n_build_loop ];
1779  if ( cgrelax_loop > 1 || option[ pepspec::use_input_bb ] ) {
1780  perturb_pep_bb( pose, mm_move, cen_scorefxn, n_cgperturb_iter );
1781  }
1782 
1783  //debug dump CG data and pdb//
1784  if ( option[ pepspec::dump_cg_bb ] ) {
1785  //remove buffer res
1786  std::string pdb_name( pdb_dir + "/" + out_nametag + "_" + string_of( peptide_loop ) + ".pdb" );
1787  if ( n_cgrelax_loop > 1 ) pdb_name = pdb_dir + "/" + out_nametag + "_" + string_of( peptide_loop ) + "_" + string_of( cgrelax_loop ) + ".pdb";
1788  if ( option[ pepspec::save_low_pdbs ] ) pose.dump_scored_pdb( pdb_name, *full_scorefxn );
1789  ( *cen_scorefxn )( pose );
1790  out_file << pdb_name + "\t"<<pose.energies().total_energies().weighted_string_of( cen_scorefxn->weights() );
1791  out_file<<"\ttotal_score:\t"<<pose.energies().total_energies()[ total_score ]<<"\t";
1792  // if( option[ pepspec::rmsd_analysis ] ) out_file << pep_rmsd_analysis( pose );
1793  // if( option[ pepspec::phipsi_analysis ] ) out_file << pep_phipsi_analysis( pose );
1794  out_file<<std::endl;
1795  continue;
1796  }
1797 
1798  //switch back to fullatom
1799  if ( option[ pepspec::add_buffer_res ] ) {
1800  core::pose::remove_variant_type_from_pose_residue( pose, core::chemical::VIRTUAL_BB, pep_begin );
1801  core::pose::remove_variant_type_from_pose_residue( pose, core::chemical::VIRTUAL_BB, pep_end );
1802  }
1803  core::util::switch_to_residue_type_set( pose, core::chemical::FA_STANDARD );
1804  if ( option[ pepspec::add_buffer_res ] ) {
1805  core::pose::add_variant_type_to_pose_residue( pose, core::chemical::VIRTUAL_BB, pep_begin );
1806  core::pose::add_variant_type_to_pose_residue( pose, core::chemical::VIRTUAL_BB, pep_end );
1807  }
1808 
1809  //replace prot residues w/ original rotamers
1810  for ( Size resnum = prot_begin; resnum <= prot_end; ++resnum ) {
1811  pose.replace_residue( resnum, start_pose.residue( resnum ), true );
1812  }
1813  pose.replace_residue( pep_anchor, pep_anchor_res, true );
1814 
1815  //replace termini
1816  add_termini( pose );
1817  //turn off csts//
1818  full_scorefxn->set_weight( coordinate_constraint, 0 );
1819 
1820  //randomize pep sequence//
1821  if ( !( option[ pepspec::no_design ] || option[ pepspec::input_seq ].user() ) ) {
1822  for ( Size mut_site = pep_begin; mut_site <= pep_end; ++mut_site ) {
1823  if ( mut_site==pep_anchor ) continue;
1824  int resindex;
1825  resindex = static_cast< int > ( 20 * numeric::random::rg().uniform() + 1 );
1826  make_sequence_change( mut_site, chemical::AA(resindex), pose );
1827  }
1828  }
1829 
1830 
1831  //define neighbors
1832  vector1< bool > is_pep_nbr( pose.total_residue(), false );
1833  vector1< Size > nbr_res_vec;
1834  for ( Size i=1; i<= pose.total_residue(); ++i ) {
1835  Residue const & rsd1( pose.residue(i) );
1836  if ( is_pep[i] ) continue;
1837  for ( Size j=1; j<= pose.total_residue(); ++j ) {
1838  Residue const & rsd2( pose.residue(j) );
1839  if ( !is_pep[j] ) continue;
1840  if ( is_pep_nbr[i] ) break;
1841  for ( Size ii=1; ii<= rsd1.natoms(); ++ii ) {
1842  if ( is_pep_nbr[i] ) break;
1843  for ( Size jj=1; jj<= rsd2.natoms(); ++jj ) {
1844  if ( rsd1.xyz(ii).distance( rsd2.xyz(jj) ) < cutoff ) {
1845  is_pep_nbr[i] = true;
1846  nbr_res_vec.push_back( i );
1847  break;
1848  }
1849  }
1850  }
1851  }
1852  }
1853  mm_min->set_chi( is_pep_nbr );
1854 
1855  //define movers//
1856  protocols::simple_moves::MinMoverOP min_mover( new protocols::simple_moves::MinMover( mm_min, full_scorefxn, "dfpmin", 0.001, true ) );
1857 
1858  MonteCarloOP mc_relax( new MonteCarlo( pose, *full_scorefxn, 1.0 ) );
1859 
1860  Real bind_score( get_binding_score( pose, pep_chain, full_scorefxn ) );
1861  myMC mc_bind( pose, bind_score, 1.0 );
1862 
1863  //define design task and repack task
1864  pack::task::TaskFactoryOP dz_task_factory( new pack::task::TaskFactory );
1865  {
1866  pack::task::operation::RestrictResidueToRepackingOP restrict_to_repack_taskop( new pack::task::operation::RestrictResidueToRepacking() );
1867  pack::task::operation::PreventRepackingOP prevent_repack_taskop( new pack::task::operation::PreventRepacking() );
1868  for ( Size i=1; i<= pose.total_residue(); ++i ) {
1869  if ( is_pep[ i ] && i != pep_anchor ) {
1870  if ( option[ pepspec::no_design ] || ( option[ pepspec::input_seq ].user() && input_seq[ i - pep_begin ] != 'X' ) ) {
1871  restrict_to_repack_taskop->include_residue( i );
1872  }
1873  } else if ( is_pep_nbr[i] && is_prot[i] ) {
1874  restrict_to_repack_taskop->include_residue( i );
1875  } else {
1876  prevent_repack_taskop->include_residue( i );
1877  }
1878  }
1879  using core::pack::task::operation::TaskOperationCOP;
1880  if ( !option[ pepspec::diversify_pep_seqs ].user() ) dz_task_factory->push_back( TaskOperationCOP( new pack::task::operation::InitializeFromCommandline() ) );
1881  dz_task_factory->push_back( TaskOperationCOP( new pack::task::operation::IncludeCurrent() ) );
1882  dz_task_factory->push_back( restrict_to_repack_taskop );
1883  dz_task_factory->push_back( prevent_repack_taskop );
1884  }
1885  std::string wts_name( option[ score::weights ] ), soft_wts_name( option[ pepspec::soft_wts ] );
1886  //full_wts
1887  {
1888  pack::task::PackerTaskOP dz_task( dz_task_factory->create_task_and_apply_taskoperations( pose ));
1889 
1890  //upweight interface?
1891  if ( option[ pepspec::upweight_interface ] ) {
1892  Real upweight( 2.0 );
1893  core::pack::task::IGEdgeReweightContainerOP IGreweight = dz_task->set_IGEdgeReweights();
1894  core::pack::task::IGEdgeReweighterOP upweighter( new protocols::toolbox::ResidueGroupIGEdgeUpweighter( upweight, pep_res_vec, nbr_res_vec ) );
1895  IGreweight->add_reweighter( upweighter );
1896  }
1897 
1898  protocols::simple_moves::PackRotamersMoverOP dz_pack( new protocols::simple_moves::PackRotamersMover( full_scorefxn, dz_task, 1 ) );
1899  protocols::simple_moves::RotamerTrialsMoverOP dz_rottrial( new protocols::simple_moves::RotamerTrialsMover( full_scorefxn, dz_task_factory ) );
1900  SequenceMoverOP design_seq( new SequenceMover );
1901  design_seq->add_mover( dz_pack );
1902  design_seq->add_mover( dz_rottrial );
1903  design_seq->add_mover( min_mover );
1904  ( *full_scorefxn )( pose );
1905  design_seq->apply( pose );
1906  if ( option[ pepspec::binding_score ] ) {
1907  bind_score = get_binding_score( pose, pep_chain, full_scorefxn );
1908  mc_bind.roll( pose, bind_score );
1909  } else mc_relax->boltzmann( pose );
1910  pdb_name = pdb_dir + "/" + out_nametag + "_" + string_of( peptide_loop ) + "_full.pdb";
1911  print_pep_analysis( pdb_name, out_file, pose, prot_score, full_scorefxn, save_all_pdbs );
1912  }
1913  //design again with soft_wts
1914  if ( option[ pepspec::soft_wts ].user() ) {
1915  pack::task::PackerTaskOP dz_task( dz_task_factory->create_task_and_apply_taskoperations( pose ));
1916 
1917  //upweight interface?
1918  if ( option[ pepspec::upweight_interface ] ) {
1919  Real upweight( 2.0 );
1920  core::pack::task::IGEdgeReweightContainerOP IGreweight = dz_task->set_IGEdgeReweights();
1921  core::pack::task::IGEdgeReweighterOP upweighter( new protocols::toolbox::ResidueGroupIGEdgeUpweighter( upweight, pep_res_vec, nbr_res_vec ) );
1922  IGreweight->add_reweighter( upweighter );
1923  }
1924 
1925  protocols::simple_moves::PackRotamersMoverOP dz_pack( new protocols::simple_moves::PackRotamersMover( soft_scorefxn, dz_task, 1 ) );
1926  protocols::simple_moves::RotamerTrialsMoverOP dz_rottrial( new protocols::simple_moves::RotamerTrialsMover( soft_scorefxn, dz_task_factory ) );
1927  SequenceMoverOP design_seq( new SequenceMover );
1928  design_seq->add_mover( dz_pack );
1929  design_seq->add_mover( dz_rottrial );
1930  design_seq->add_mover( min_mover );
1931  ( *soft_scorefxn )( pose );
1932  design_seq->apply( pose );
1933  if ( option[ pepspec::binding_score ] ) {
1934  bind_score = get_binding_score( pose, pep_chain, full_scorefxn );
1935  mc_bind.roll( pose, bind_score );
1936  } else mc_relax->boltzmann( pose );
1937  pdb_name = pdb_dir + "/" + out_nametag + "_" + string_of( peptide_loop ) + "_soft.pdb";
1938  print_pep_analysis( pdb_name, out_file, pose, prot_score, full_scorefxn, save_all_pdbs );
1939  }
1940  //now try MCM desin against binding score
1941  if ( option[ pepspec::diversify_pep_seqs ] ) {
1942  Size iter_per_res( option[ pepspec::diversify_lvl ] );
1943  for ( Size ii = 1; ii <= iter_per_res * ( pep_end - pep_begin + 1 ); ++ii ) {
1944  mutate_random_residue( pose, is_mutable, soft_scorefxn, full_scorefxn );
1945  pdb_name = pdb_dir + "/" + out_nametag + "_" + string_of( peptide_loop ) + "_mc" + string_of( ii ) + ".pdb";
1946  print_pep_analysis( pdb_name, out_file, pose, prot_score, full_scorefxn, save_all_pdbs );
1947  if ( option[ pepspec::binding_score ] ) {
1948  bind_score = get_binding_score( pose, pep_chain, full_scorefxn );
1949  mc_bind.roll( pose, bind_score );
1950  } else mc_relax->boltzmann( pose );
1951  }
1952  }
1953  //recover best
1954  if ( option[ pepspec::binding_score ] ) mc_bind.recover_low( pose );
1955  else mc_relax->recover_low( pose );
1956  pdb_name = pdb_dir + "/" + out_nametag + "_" + string_of( peptide_loop ) + ".pdb";
1957  print_pep_analysis( pdb_name, out_file, pose, prot_score, full_scorefxn, save_low_pdbs );
1958 
1959  }
1960  }
1961 }
1962 
1963 
1964 void*
1965 my_main( void*)
1966 {
1967 
1968  RunPepSpec();
1969  exit(0);
1970 
1971 }
1972 
1973 int
1974 main( int argc, char * argv [] )
1975 {
1976  try {
1977  using namespace basic::options;
1978  using namespace basic::options::OptionKeys;
1979 
1980  protocols::init::init(argc, argv);
1981 
1982  protocols::viewer::viewer_main( my_main );
1983  TR.flush();
1984  } catch ( utility::excn::EXCN_Base const & e ) {
1985  std::cout << "caught exception " << e.msg() << std::endl;
1986  return -1;
1987  }
1988 }
get_binding_score
Real get_binding_score(Pose pose, Size pep_chain, ScoreFunctionOP full_scorefxn)
Definition: pepspec.cc:1414
utility_exit_with_message
#define utility_exit_with_message(m)
Exit with file + line + message.
Definition: exit.hh:47
make_frags
core::fragment::FragSetCOP make_frags(core::Size const start, core::Size const stop, std::string const &seq)
helper code for fragments generation, copied from S.M.Lewis
Definition: pepspec.cc:977
MetricValue.hh
prot_begin
Size prot_begin(0)
prot_anchor
Size prot_anchor(0)
basic::database::schema_generator::ConstraintOP
utility::pointer::shared_ptr< Constraint > ConstraintOP
Definition: Constraint.fwd.hh:25
THREAD_LOCAL
#define THREAD_LOCAL
Definition: backwards_thread_local.hh:17
remove_pep_res
void remove_pep_res(pose::Pose &pose, bool add_nterm, bool add_cterm)
Definition: pepspec.cc:825
utility::excn::EXCN_Base
Definition: EXCN_Base.hh:49
utility::excn::EXCN_Base::msg
virtual std::string const msg() const
Definition: EXCN_Base.hh:70
myMC::myMC
myMC(Pose pose, Real score, Real temp_init)
Definition: pepspec.cc:225
pep_begin
Size pep_begin(0)
basic::MetricValue::value
T const & value() const
Definition: MetricValue.hh:60
numeric::random::rg
RandomGenerator & rg()
Return the one-per-thread "singleton" random generator.
Definition: random.cc:46
pep_chain
Size pep_chain(0)
pep_coord_cst::sd
Real sd
Definition: pepspec.cc:197
myMC::low_score_
Real low_score_
Definition: pepspec.cc:212
app.surface_docking.surface_docking.constraints
constraints
Definition: surface_docking.py:46
input_seq
std::string input_seq
Definition: pepspec.cc:186
pep_coord_cst::x0
Real x0
Definition: pepspec.cc:196
mutate_random_residue
void mutate_random_residue(Pose &pose, vector1< bool > is_mutable, ScoreFunctionOP soft_scorefxn, ScoreFunctionOP full_scorefxn)
Definition: pepspec.cc:1309
myMC::low_pose_
Pose low_pose_
Definition: pepspec.cc:211
myMC::temp_
Real temp_
Definition: pepspec.cc:213
ObjexxFCL::getline
std::istream & getline(std::istream &stream, Fstring &s)
Get Line from Stream.
Definition: Fstring.cc:1610
pep_end
Size pep_end(0)
app.delta_score_per_mutation.mm
tuple mm
Definition: delta_score_per_mutation.py:53
utility::options::OptionKeys::options::user
BooleanOptionKey const user("options:user")
Definition: OptionKeys.hh:40
Pose
core::pose::Pose Pose
Definition: supercharge.cc:101
get_clash_pairs
vector1< std::pair< Size, Size > > get_clash_pairs(pose::Pose pose, vector1< bool > is_checked, scoring::ScoreFunctionOP const &scorefxn, Real const clash_threshold)
Definition: pepspec.cc:582
test.T200_Scoring.weights
tuple weights
Definition: T200_Scoring.py:51
file_sys_util.hh
Platform independent operations on files (except I/O)
aa2index
Size aa2index(chemical::AA aa)
Definition: pepspec.cc:282
has_clash
bool has_clash(pose::Pose pose, vector1< bool > is_checked, scoring::ScoreFunctionOP const &scorefxn, Real const clash_threshold)
Definition: pepspec.cc:517
myMC::recover_low
void recover_low(Pose &pose)
Definition: pepspec.cc:267
pep_coord_cst::x
Real x
Definition: pepspec.cc:193
pep_coord_cst
Definition: pepspec.cc:190
app.surface_docking.surface_docking.pose
tuple pose
Definition: surface_docking.py:203
random.hh
Random number generator system.
Exceptions.hh
common derived classes for thrown exceptions
basic::options::start_file
std::string start_file()
kind of like the old -s – just one file
Definition: util.cc:41
add_pep_res
void add_pep_res(pose::Pose &pose, bool add_nterm, bool add_cterm)
Definition: pepspec.cc:791
prot_end
Size prot_end(0)
utility::file::create_directory
bool create_directory(std::string const &dir_path)
Create a directory if it doesn't already exist.
Definition: file_sys_util.cc:295
PyMOL_demo.scorefxn
tuple scorefxn
Definition: PyMOL_demo.py:63
ObjexxFCL::index
Fstring::size_type index(Fstring const &s, Fstring const &ss)
First Index Position of a Substring in an Fstring.
Definition: Fstring.hh:2180
myMC::roll
void roll(Pose &pose, Real &score)
Definition: pepspec.cc:240
make_1mer_frags
core::fragment::FragSetCOP make_1mer_frags(core::Size const seqpos_start, core::Size const seqpos_stop, std::string const &seq, Size const nfrags)
Definition: pepspec.cc:1001
gen_pep_bb_sequential
void gen_pep_bb_sequential(pose::Pose &pose, scoring::ScoreFunctionOP cen_scorefxn)
Definition: pepspec.cc:1036
add_termini
void add_termini(pose::Pose &pose)
Definition: pepspec.cc:778
packmin_unbound_pep
void packmin_unbound_pep(pose::Pose &pose, scoring::ScoreFunctionOP full_scorefxn)
Definition: pepspec.cc:1395
test.T200_Scoring.rsd1
tuple rsd1
Definition: T200_Scoring.py:58
pep_coord_cst::pep_pos
int pep_pos
Definition: pepspec.cc:192
demo.D090_Ala_scan.interface_cutoff
tuple interface_cutoff
Definition: D090_Ala_scan.py:619
myMC::last_score_
Real last_score_
Definition: pepspec.cc:210
test.T200_Scoring.rsd2
tuple rsd2
Definition: T200_Scoring.py:59
docking.p
tuple p
Definition: docking.py:9
open.hh
Functions for opening database files.
main
int main(int argc, char *argv[])
Definition: pepspec.cc:1974
myMC::accept
bool accept()
Definition: pepspec.cc:275
numeric::Vector
xyzVector< Real > Vector
Definition: cyclic_coordinate_descent.cc:22
set_pep_csts
void set_pep_csts(pose::Pose &pose)
Definition: pepspec.cc:916
Tracer.hh
Tracer IO system.
prot_chain
Size prot_chain(0)
ObjexxFCL::abs
T abs(T const &x)
std::abs( x ) == | x |
Definition: Fmath.hh:295
test.T200_Scoring.ii
tuple ii
Definition: T200_Scoring.py:70
utility::vector1
std::vector with 1-based indexing
Definition: vector1.fwd.hh:44
test.T200_Scoring.emap
tuple emap
Definition: T200_Scoring.py:61
test.G100_Functions_Input.score
tuple score
Definition: G100_Functions_Input.py:27
app.delta_score_per_mutation.task
tuple task
Definition: delta_score_per_mutation.py:60
myMC
Definition: pepspec.cc:202
options
rule< Scanner, options_closure::context_t > options
Definition: Tag.cc:377
make_sequence_change
void make_sequence_change(Size const seqpos, AA const &new_aa, pose::Pose &pose)
Definition: pepspec.cc:497
RunPepSpec
void RunPepSpec()
Definition: pepspec.cc:1495
my_main
void * my_main(void *)
Definition: pepspec.cc:1965
clean_pdb_keep_ligand.resnum
list resnum
if line_edit[13:14]=='P': #Nucleic acid? Skip.
Definition: clean_pdb_keep_ligand.py:174
ContactMap.native
list native
Definition: ContactMap.py:108
TR
static THREAD_LOCAL basic::Tracer TR("apps.public.pepspec")
numeric::Real
double Real
Definition: types.hh:39
pep_jump
Size pep_jump(2)
myMC::accept_
bool accept_
Definition: pepspec.cc:214
perturb_pep_bb
void perturb_pep_bb(pose::Pose &pose, kinematics::MoveMapOP mm_move, scoring::ScoreFunctionOP cen_scorefxn, Size n_iter)
Definition: pepspec.cc:1280
pep_coord_cst::atom_name
std::string atom_name
Definition: pepspec.cc:191
initialize_peptide
void initialize_peptide(pose::Pose &pose)
Definition: pepspec.cc:851
app.delta_score_per_mutation.min_mover
tuple min_mover
Definition: delta_score_per_mutation.py:66
loops.pack
tuple pack
Definition: loops.py:50
utility::io::oc::cout
ocstream cout(std::cout)
Wrapper around std::cout.
Definition: ocstream.hh:287
pep_coord_cst::z
Real z
Definition: pepspec.cc:195
basic::database::schema_generator::ConstraintCOP
utility::pointer::shared_ptr< Constraint const > ConstraintCOP
Definition: Constraint.fwd.hh:27
utility::options::OptionKeys::options::exit
BooleanOptionKey const exit("options:exit")
Definition: OptionKeys.hh:51
vector1.hh
vector1: std::vector with 1-based indexing
get_n_pep_nbrs
Size get_n_pep_nbrs(pose::Pose const &pose, vector1< bool > const is_pep, Real const cutoff_cg)
Definition: pepspec.cc:355
pep_phipsi_analysis
std::string pep_phipsi_analysis(pose::Pose pose)
Definition: pepspec.cc:713
test.G100_Functions_Input.cst_file
tuple cst_file
Definition: G100_Functions_Input.py:34
app.surface_docking.surface_docking.start_pose
tuple start_pose
Definition: surface_docking.py:145
basic::Tracer
Class for handling user debug/warnings/errors. Use instance of this class instead of 'std::cout' for ...
Definition: Tracer.hh:134
app.delta_score_per_mutation.res
tuple res
Definition: delta_score_per_mutation.py:93
demo.D010_Pose_structure.pdb_filename
pdb_filename
Definition: D010_Pose_structure.py:175
myMC::last_pose_
Pose last_pose_
Definition: pepspec.cc:209
pep_anchor
Size pep_anchor(0)
pep_coord_cst::tol
Real tol
Definition: pepspec.cc:198
numeric::random::uniform
double uniform()
Generate a random number between 0 and 1. Threadsafe since each thread uses its own random generator...
Definition: random.cc:56
PyMOL_demo.seq
tuple seq
Definition: PyMOL_demo.py:72
ObjexxFCL::string_of
std::string string_of(Fstring const &s)
string of an Fstring
Definition: Fstring.hh:2889
print_pep_analysis
void print_pep_analysis(std::string pdb_name, std::fstream &out_file, Pose pose, Real prot_score, ScoreFunctionOP full_scorefxn, bool dump_pdb)
Definition: pepspec.cc:1432
basic::init
void init()
set global 'init_was_called' to true
Definition: init.cc:26
pep_coord_csts
vector1< pep_coord_cst > pep_coord_csts
Definition: pepspec.cc:200
xyzVector.hh
Fast (x,y,z)-coordinate numeric vector.
max
static T max(T x, T y)
Definition: Svm.cc:19
pep_rmsd_analysis
std::string pep_rmsd_analysis(pose::Pose pose)
Definition: pepspec.cc:635
numeric::kinematic_closure::frame
void frame(const utility::vector1< utility::vector1< Real > > &R, utility::vector1< utility::vector1< Real > > &U)
Definition: bridgeObjects.cc:244
numeric::random::Size
platform::Size Size
Definition: random.fwd.hh:30
pep_coord_cst::y
Real y
Definition: pepspec.cc:194
demo.D060_Folding.f
tuple f
Definition: D060_Folding.py:445
option
rule< Scanner, option_closure::context_t > option
Definition: Tag.cc:378
Generated on Sun Jan 17 2016 15:21:12 for Rosetta by   doxygen 1.8.7