Manual active models

SCRIPTING.md

@@ -448,6 +448,22 @@ in exactly this way by transparently replacing `AND` with `IF` behind the scenes
 When using `IF` with `GOTO`, it is possible to create loops as seen in the last example above.
 
 
+# INTC
+Examples:
+```
+INTC
+INTC ~
+INTC &clashes
+INTC A &clashes
+```
+
+Obtains the total internal clashes of the specified strand, or the working strand if no strand ID letter is given.
+
+If a float variable is specified, that variable will receive the internal clash total. Otherwise, it will be output to the command line.
+
+In the absence of an output variable, a tilde will prevent outputting a newline after the clash total.
+
+
 # LET
 Examples:
 ```
@@ -744,3 +760,22 @@ UNLIG
 ```
 
 Deletes any ligands that were loaded from a PDB by a `LOAD` command. Takes no arguments.
+
+
+# WORST
+Example:
+```
+WORST ~
+WORST %resno1 %resno2
+WORST $res1 $res2
+```
+
+Used after `INTC` to retrieve the two residues with the largest clashes.
+
+If no output variables are specified, the residue names (e.g. Asn152) will be written to the command line, separated by a space.
+The optional tilde means no newline at the end of the output.
+
+If integer variables are specified, the variables will be set to the residue numbers of the clashing residues. If string variables are specified,
+they will be set to the names of the clashing residues.
+
+Using `WORST` without an `INTC` causes undefined behavior.

predict/manual_tne.pepd

@@ -0,0 +1,124 @@
+LOAD "pdbs/OR5/OR5P3.upright.pdb"
+
+# Configurable variables.
+LET &vrot6 = -5.3
+LET @xlat6 = [0,3.0,0]
+LET &tilt6a = -10
+LET &tilt6b = -5
+LET &tilt6c = -6
+LET &tilt5 = -0.3
+LET &mov5 = 1.0
+LET &mov5x = 3.3
+LET &mov5z = 2.0
+LET &tilt57 = 5
+LET &bend5 = 7
+LET &room76a = 7
+LET &room76b = 6
+LET &room7E = 5
+LET &room72 = 3
+LET &tilt1 = -7
+
+# Internal variables.
+IF $45.51 = "D" THEN LET $ExO1 = "OD1"
+ELSE IF $45.51 = "E" THEN LET $ExO1 = "OE1"
+LET $ExO2 = $ExO1 FROM 1 FOR 2
+LET $ExO2 += "2"
+
+
+# Initial measurements.
+BRIDGE %45.51 %6.55
+MEASURE %45.51 $ExO1 %6.55 "OH" &d6
+MEASURE %45.51 $ExO2 %6.55 "OH" &d6a
+IF &d6a < &d6 THEN LET &d6 = &d6a
+ECHO "Old 45-6 Distance: " &d6
+MEASURE %5.58 "OH" %7.53 "OH" &d57
+ECHO "Old 5-7 Distance: " &d57
+
+
+# Get TMR6 and EXR2 bridged.
+PTALIGN @45.51 @6.55 @6.48 @axis &angle
+ROTATE @axis &vrot6 %6.48 %6.28 %6.58
+MOVEREL %6.28 %6.58 @xlat6
+
+# Trip switch and CGSHL histidine.
+ATOMTO %6.48 "OH" @3.43
+ATOMTO %6.47 "CZ" @3.37
+ATOMTO %6.40 "NE2" @5.68
+
+# Outward tilt of entire TMR6.
+LET @axis = @7.53 - @5.58
+ROTATE @axis &tilt6a %6.55 %6.28 %6.58
+LET @axis = @2.40 - @5.58
+ROTATE @axis &tilt6b %6.55 %6.28 %6.58
+LET @axis = @3.38 - @6.48
+ROTATE @axis &tilt6c %6.49 %6.28 %6.49
+
+# Slight tilt of TMR5 away from TMR3.
+LET @axis = @2.40 - @5.58
+ROTATE @axis &tilt5 %5.32 %5.32 %5.68
+
+# Move TMR5 closer to TMR7.
+LET &magn = @axis
+LET &mult = &mov5 / &magn
+LET @axis *= &mult
+LET @axis.x += &mov5x
+LET @axis.z -= &mov5z
+MOVEREL %5.32 %5.68 @axis
+
+# Tilt entire TMR5 toward TMR7.
+PTALIGN @5.68 @7.53 @5.32 @axis &angle
+ROTATE @axis &tilt57 %5.32 %5.32 %5.68
+
+# Slight bend of TMR5.
+PTALIGN @5.68 @7.53 @5.50 @axis &angle
+ROTATE @axis &bend5 %5.50 %5.32 %5.68
+
+# TMR7 to make room for TMR6.
+PTALIGN @45.54 @7.31 @7.49 @axis &angle
+ROTATE @axis &room76b %7.49 %7.31 %7.49
+PTALIGN @5.36 @7.31 @7.49 @axis &angle
+ROTATE @axis &room76a %7.49 %7.31 %7.49
+
+# TMR7 away from TMR2.
+PTALIGN @2.43 @7.53 @7.49 @axis &angle
+ROTATE @axis &room7E %7.49 %7.49 %7.56
+
+# TMR7 away from EXR2.
+PTALIGN @45.57 @7.31 @7.49 @axis &angle
+ROTATE @axis &room72 %7.49 %7.31 %7.49
+
+# TMR1 tilt.
+LET @axis = @2.38 - @1.57
+ROTATE @axis &tilt1 %1.58 %1.28 %1.58
+
+
+BRIDGE %45.51 %6.55
+MEASURE %45.51 $ExO1 %6.55 "OH" &d6
+MEASURE %45.51 $ExO2 %6.55 "OH" &d6a
+IF &d6a < &d6 THEN LET &d6 = &d6a
+ECHO "New 45-6 Distance: " &d6
+ATOMTO %7.53 "OH" @5.58
+BRIDGE %5.58 %7.53
+MEASURE %5.58 "OH" %7.53 "OH" &d57
+ECHO "New 5-7 Distance: " &d57
+
+
+DELETE 1 %1.35
+DELETE %1.59 %2.37
+DELETE %2.68 %3.20
+DELETE %3.57 %4.38
+DELETE %4.66 %45.45
+DELETE %45.57 %5.31
+DELETE %5.69 %6.27
+DELETE %6.60 %7.30
+DELETE %7.57 9999
+
+MINC
+MINC
+ECHO "Total clashes: " ~
+INTC ~
+ECHO ", worst between " ~
+WORST ~
+ECHO "."
+
+SAVE "pdbs/OR5/OR5P3.bound.pdb"

src/classes/atom.cpp

@@ -2225,7 +2225,7 @@ void Atom::save_pdb_line(FILE* pf, unsigned int atomno)
     fprintf(pf, "%c%c%c %c", aa3let[0] & 0x5f, aa3let[1] & 0x5f, aa3let[2] & 0x5f, pdbchain);
 
     sprintf(numbuf, "%d", residue);
-    fprintf(pf, "%4s    ", numbuf);
+    fprintf(pf, "%4s     ", numbuf);
 
     if (!location.x) location.x = 0;
     sprintf(numbuf, "%4.3f", location.x);
@@ -2239,7 +2239,7 @@ void Atom::save_pdb_line(FILE* pf, unsigned int atomno)
     sprintf(numbuf, "%4.3f", location.z);
     fprintf(pf, "%7s ", numbuf);
 
-    fprintf(pf, "  1.00001.00           %s%c\n", get_elem_sym(), fabs(charge) > hydrophilicity_cutoff ? (charge > 0 ? '+' : '-') : ' ' );
+    fprintf(pf, " 1.00  1.00           %s%c\n", get_elem_sym(), fabs(charge) > hydrophilicity_cutoff ? (charge > 0 ? '+' : '-') : ' ' );
 }
 
 void Atom::stream_pdb_line(ostream& os, unsigned int atomno)
@@ -2278,7 +2278,7 @@ void Atom::stream_pdb_line(ostream& os, unsigned int atomno)
     if (!location.z) location.z = 0;
     os << fixed << setprecision(3) << setw(8) << setfill(' ') << location.z;
 
-    os << "  1.00001.00           " << get_elem_sym();
+    os << "  1.00  1.00           " << get_elem_sym();
     os << (fabs(charge) > hydrophilicity_cutoff ? (charge > 0 ? '+' : '-') : ' ');
     os << endl;
 }

src/classes/protein.cpp

@@ -506,23 +506,17 @@ float Protein::get_internal_clashes(int sr, int er, bool repack, int repack_iter
                     if (r > rr) continue;
                 }
 
-                Atom* ia = residues[i]->get_nearest_atom(residues[j]->get_CA_location());
-                Atom* ja = residues[j]->get_nearest_atom(ia->get_location());
-                r = ia->distance_to(ja);
-                if (r > ia->get_vdW_radius() + ja->get_vdW_radius()) continue;
-
                 if (abs(resno2-resno) <= 1) continue;
-                if (resno2 >= sr && resno2 <= er) continue;
 
                 float f = residues[i]->get_intermol_clashes(residues[j]);
                 result += f;
+                // if (f > 1e4) cout << *residues[i] << " ! " << *residues[j] << ": " << f << endl;
                 SCoord dpos = residues[j]->get_CA_location().subtract(residues[i]->get_CA_location());
                 dpos.r = f;
                 clashttl = clashttl.add(dpos);
-                float limit = unconnected_residue_mindist + sqrt(residues[i]->get_reach()) + sqrt(residues[j]->get_reach());
-                f = fmax(limit - dpos.r, 0);
                 if (f > clash_worst)
                 {
+                    // cout << *residues[i] << " ! " << *residues[j] << " are the worst: " << f << endl;
                     clash_worst = f;
                     stop1 = residues[i];
                     stop2 = residues[j];

src/classes/scoring.cpp

@@ -211,6 +211,22 @@ DockResult::DockResult(Protein* protein, Molecule* ligand, Point size, int* addl
     {
         if (!reaches_spheroid[i]) continue;
         if (!protein->aa_ptr_in_range(reaches_spheroid[i])) continue;
+
+        Atom* rsa = reaches_spheroid[i]->get_nearest_atom(ligand->get_barycenter());
+        if (!rsa) continue;
+        if (rsa->is_backbone) continue;
+        Atom* la = ligand->get_nearest_atom(rsa->get_location());
+        if (!la) continue;
+        rsa = ligand->get_nearest_atom(la->get_location());
+        if (!rsa) continue;
+        la = ligand->get_nearest_atom(rsa->get_location());
+        if (!la) continue;
+        rsa = ligand->get_nearest_atom(la->get_location());
+        if (!rsa) continue;
+        if (rsa->is_backbone) continue;
+        float r = la->distance_to(rsa);
+        if (r > _INTERA_R_CUTOFF) continue;
+
         reaches_spheroid[i]->clear_atom_binding_energies();
         int resno = reaches_spheroid[i]->get_residue_no();
 
@@ -232,7 +248,7 @@ DockResult::DockResult(Protein* protein, Molecule* ligand, Point size, int* addl
             res_clash_dir[resno] = res_clash_dir[resno].add(clashdir);
         }
 
-        lb -= fmax(reaches_spheroid[i]->total_eclipses() - reaches_spheroid[i]->initial_eclipses, 0);
+        // lb -= (reaches_spheroid[i]->total_eclipses() - reaches_spheroid[i]->initial_eclipses);
 
         #if _dbg_51e2_ionic
         if (resno == 262)