Bug fixe: only heavy atoms when atom list set to all + upp_bound not fixed...

Bug fixe: only heavy atoms when atom list set to all + upp_bound not fixed when ambig_restraint turn on

ariaec/conf/aria_ec.ini

@@ -34,8 +34,9 @@ hb_dplus:                                   0.5
 ; --------------------------- Contacts parameters ---------------------------- #
 ; evfold_weight                 : use EVFold weight -> 10/i (i:contact rank)
 ; neighborhood_contact          : add neighbors when writing restraint
-; pair_list                     : all, min [min]; use all atom pairs or
-;                                 minimized list (CA-CA, CB-CB, SC-SC)
+; pair_list                     : all, heavy, min [min]; use all atom pairs,
+;                                 heavy atms or minimized list
+;                                 (CA-CA, CB-CB, SC-SC)
 ; ambiguous_distance_restraint  : use Ambiguous Distance Restraints
 ; ambiguous_distance_type       : min, deff [min]; compute Deff from target
 ;                                 distance iif ADR
@@ -53,6 +54,7 @@ groupby_method:                             min
 restraint_distance:                         2.5
 lower_bound:                                1.0
 def_upper_bound:                            4.0
+; Specific tresholds only for unambig restraints
 ca_upper_bound:                             7.0
 cb_upper_bound:                             7.0
 ; ---------------------------- Filter parameters ----------------------------- #
@@ -106,10 +108,10 @@ md_parameters_random_seed:                  89764443
 md_parameters_steps_high:                   10000
 md_parameters_steps_cool1:                  5000
 md_parameters_steps_cool2:                  4000
-iteration_n_structures:                     1000
+iteration_n_structures:                     100
 iteration_n_best_structures:                15
 iteration_n_kept_structures:                0
-max_contribution:                           20
+max_contribution:                           1000
 merging_method:                             standard
 violation_analysis_violation_threshold:     0.5
 water_refinement_solvent:                   water

ariaec/econverter.py

@@ -454,16 +454,7 @@ class AriaEcXMLConverter(AriaXMLConverter):
         conv_table = ConversionTable.ConversionTable().table['AMINO_ACID'][
             'iupac']
 
-        # Dict giving all atms foreach aa
-        aa_atm = dict((AminoAcid.AminoAcid(aa)[0],
-                       filter(contactmap.heavy_reg.match, v.keys())) for
-                      aa, v in conv_table.items())
-
-        # Dict giving atoms product for each res-res pair
-        aa_atm_pair = {(aa1, aa2): list(itertools.product(atms1, atms2))
-                       for aa1, atms1 in aa_atm.items()
-                       for aa2, atms2 in aa_atm.items()
-                       }
+        # TODO: atm_pair method
         if pair_flag == "min":
             def min_atms(aa1, aa2, atms):
                 # Function to minimize atom pair list between aa1 & aa2
@@ -472,6 +463,51 @@ class AriaEcXMLConverter(AriaXMLConverter):
                         ('CA', 'CA'),
                         ('CB', 'CB'),
                         self.settings.scsc_min[aa1][aa2])]
+            aa_atm = dict((AminoAcid.AminoAcid(aa)[0],
+                           filter(contactmap.heavy_reg.match, atms.keys()))
+                          for aa, atms in conv_table.items())
+            aa_atm_pair = \
+                {
+                    (aa1, aa2): min_atms(aa1, aa2, list(itertools.product(atms1,
+                                                                          atms2)))
+                    for aa1, atms1 in aa_atm.items()
+                    for aa2, atms2 in aa_atm.items()
+                }
+        elif pair_flag == "heavy":
+            # Dict giving all atms foreach aa
+            aa_heav_atm = dict((AminoAcid.AminoAcid(aa)[0],
+                                filter(contactmap.heavy_reg.match, atms.keys()))
+                               for aa, atms in conv_table.items())
+
+            # Dict giving atoms product for each res-res pair
+            aa_atm_pair = {(aa1, aa2): list(itertools.product(atms1, atms2))
+                                for aa1, atms1 in aa_heav_atm.items()
+                                for aa2, atms2 in aa_heav_atm.items()
+                                }
+        elif pair_flag == "all":
+            # Dict giving all atms foreach aa
+            aa_atm = dict((AminoAcid.AminoAcid(aa)[0],
+                           filter(contactmap.all_reg.match, atms.keys()))
+                          for aa, atms in conv_table.items())
+
+            # Dict giving atoms product for each res-res pair
+            aa_atm_pair = {(aa1, aa2): list(itertools.product(atms1, atms2))
+                           for aa1, atms1 in aa_atm.items()
+                           for aa2, atms2 in aa_atm.items()
+                           }
+        else:
+            logger.error("Wrong pair_list option. Pair_list set to min")
+
+            def min_atms(aa1, aa2, atms):
+                # Function to minimize atom pair list between aa1 & aa2
+                return [
+                    atmpair for atmpair in atms if atmpair in (
+                        ('CA', 'CA'),
+                        ('CB', 'CB'),
+                        self.settings.scsc_min[aa1][aa2])]
+            aa_atm = dict((AminoAcid.AminoAcid(aa)[0],
+                           filter(contactmap.heavy_reg.match, atms.keys()))
+                          for aa, atms in conv_table.items())
             aa_atm_pair = \
                 {
                     (aa1, aa2): min_atms(aa1, aa2, list(itertools.product(atms1,
@@ -510,7 +546,7 @@ class AriaEcXMLConverter(AriaXMLConverter):
             if adr_flag:
                 rest_id += 1
                 # TODO: Autre dist_list if target != ResMap
-                # dist_list containt all atm-atm dist foren_factorach res-res (+ neigh)
+                # dist_list containt all atm-atm dist foreach res-res (+ neigh)
                 #  if target is a  ResAtmMap or only one dist foreach res-res
                 # !!!!! Utiliser stack pour eviter les distances en doubles ????
                 if isinstance(targetdist, ResAtmMap) or len(
@@ -573,9 +609,9 @@ class AriaEcXMLConverter(AriaXMLConverter):
                         if restraint_dict.get(rest_id, None) is None:
                             # If restraint doesn't exist
                             # TODO: set upper and lower bound if sdv not None
-                            if atm_pair == ("CA", "CA"):
+                            if atm_pair == ("CA", "CA") and not adr_flag:
                                 upp_bound = self.settings.setup.config["ca_upper_bound"]
-                            elif atm_pair == ("CB", "CB"):
+                            elif atm_pair == ("CB", "CB") and not adr_flag:
                                 upp_bound = self.settings.setup.config["cb_upper_bound"]
                             else:
                                 upp_bound = self.settings.setup.config[
@@ -698,7 +734,7 @@ class AriaEcXMLConverter(AriaXMLConverter):
 
         aria_project_dict = {}
 
-        # Check if ratio cool1 / cool2 is 5:4
+        # Check if ratio cool1 / cool2 is 5:4 ... (relevant ?)
         steps_1 = float(self.settings.setup.config['md_parameters_steps_cool1'])
         steps_2 = float(self.settings.setup.config['md_parameters_steps_cool2'])
         if steps_1 / steps_2 != float(5) / int(4):

ariaec/protmap.py

@@ -511,8 +511,10 @@ class ResAtmMap(ProteinMap):
             CG    1.391024  1.379875  1.478053  2.505414  0.000000  2.790403
 
     """
+    # TODO: matrix should be not only for heavy atoms
     # Matrix only for heavy atoms.
     heavy_reg = re.compile(r"[CNOS][ABGDEZH][0-9]?")
+    all_reg = re.compile(r"^((?!cns|dyana).*)$")
     # TODO: Autre methodes de dist
     distance_method = 'euclidean'