       1 package Atom;
       2 #
       3 # $RCSfile: Atom.pm,v $
       4 # $Date: 2015/02/28 20:47:02 $
       5 # $Revision: 1.62 $
       6 #
       7 # Author: Manish Sud <msud@san.rr.com>
       8 #
       9 # Copyright (C) 2015 Manish Sud. All rights reserved.
      10 #
      11 # This file is part of MayaChemTools.
      12 #
      13 # MayaChemTools is free software; you can redistribute it and/or modify it under
      14 # the terms of the GNU Lesser General Public License as published by the Free
      15 # Software Foundation; either version 3 of the License, or (at your option) any
      16 # later version.
      17 #
      18 # MayaChemTools is distributed in the hope that it will be useful, but without
      19 # any warranty; without even the implied warranty of merchantability of fitness
      20 # for a particular purpose.  See the GNU Lesser General Public License for more
      21 # details.
      22 #
      23 # You should have received a copy of the GNU Lesser General Public License
      24 # along with MayaChemTools; if not, see <http://www.gnu.org/licenses/> or
      25 # write to the Free Software Foundation Inc., 59 Temple Place, Suite 330,
      26 # Boston, MA, 02111-1307, USA.
      27 #
      28
      29 use strict;
      30 use Carp;
      31 use Exporter;
      32 use Storable ();
      33 use Scalar::Util ();
      34 use ObjectProperty;
      35 use PeriodicTable;
      36 use Vector;
      37 use MathUtil;
      38 use Text::ParseWords;
      39 use TextUtil;
      40 use FileUtil;
      41
      42 use vars qw(@ISA @EXPORT @EXPORT_OK %EXPORT_TAGS);
      43
      44 @ISA = qw(ObjectProperty Exporter);
      45 @EXPORT = qw();
      46 @EXPORT_OK = qw();
      47
      48 %EXPORT_TAGS = (all  => [@EXPORT, @EXPORT_OK]);
      49
      50 # Setup class variables...
      51 my($ClassName, $ObjectID, %MDLValenceModelDataMap, %DaylightValenceModelDataMap);
      52 _InitializeClass();
      53
      54 # Overload Perl functions...
      55 use overload '""' => 'StringifyAtom';
      56
      57 # Class constructor...
      58 sub new {
      59   my($Class, %NamesAndValues) = @_;
      60
      61   # Initialize object...
      62   my $This = {};
      63   bless $This, ref($Class) || $Class;
      64   $This->_InitializeAtom();
      65
      66   $This->_InitializeAtomProperties(%NamesAndValues);
      67
      68   return $This;
      69 }
      70
      71 # Initialize object data...
      72 sub _InitializeAtom {
      73   my($This) = @_;
      74   my($ObjectID) = _GetNewObjectID();
      75
      76   # All other property names and values along with all Set/Get<PropertyName> methods
      77   # are implemented on-demand using ObjectProperty class.
      78   $This->{ID} = $ObjectID;
      79   $This->{Name} = "Atom ${ObjectID}";
      80   $This->{AtomSymbol} = '';
      81   $This->{AtomicNumber} = 0;
      82   $This->{XYZ} = Vector::ZeroVector;
      83 }
      84
      85 # Initialize atom properties...
      86 sub _InitializeAtomProperties {
      87   my($This, %NamesAndValues) = @_;
      88
      89   my($Name, $Value, $MethodName);
      90   while (($Name, $Value) = each  %NamesAndValues) {
      91     $MethodName = "Set${Name}";
      92     $This->$MethodName($Value);
      93   }
      94   if (!exists $NamesAndValues{'AtomSymbol'}) {
      95     carp "Warning: ${ClassName}->new: Atom object instantiated without setting atom symbol...";
      96   }
      97
      98   return $This;
      99 }
     100
     101 # Initialize class ...
     102 sub _InitializeClass {
     103   #Class name...
     104   $ClassName = __PACKAGE__;
     105
     106   # ID to keep track of objects...
     107   $ObjectID = 0;
     108
     109   # Load atom class data...
     110   _LoadAtomClassData();
     111 }
     112
     113 # Setup an explicit SetID method to block setting of ID by AUTOLOAD function...
     114 sub SetID {
     115   my($This, $Value) = @_;
     116
     117   carp "Warning: ${ClassName}->SetID: Object ID can't be changed: it's used for internal tracking...";
     118
     119   return $This;
     120 }
     121
     122 # Setup an explicit SetMolecule method to block setting of ID by AUTOLOAD function...
     123 sub SetMolecule {
     124   my($This, $Value) = @_;
     125
     126   carp "Warning: ${ClassName}->SetMolecule: Molecule property can't be changed: it's used for internal tracking...";
     127
     128   return $This;
     129 }
     130
     131 # Assign atom to  molecule...
     132 sub _SetMolecule {
     133   my($This, $Molecule) = @_;
     134
     135   $This->{Molecule} = $Molecule;
     136
     137   # Weaken the reference to disable increment of reference count; otherwise,
     138   # it it becomes a circular reference and destruction of Molecule object doesn't
     139   # get initiated which in turn disables destruction of atom object.
     140   #
     141   Scalar::Util::weaken($This->{Molecule});
     142
     143   return $This;
     144 }
     145
     146 # Setup atom symbol and atomic number for the element...
     147 #
     148 # Possible atom symbol values:
     149 #    . An element symbol or some other type of atom: L - Atom list; LP - Lone pair; R# - R group;
     150 #       A, Q, * - unknown atom; or something else?
     151 #
     152 # Default mass number corresponds to the most abundant natural isotope unless it's explicity
     153 # set using "MassNumber" property.
     154 #
     155 sub SetAtomSymbol {
     156   my($This, $AtomSymbol) = @_;
     157   my($AtomicNumber);
     158
     159   $This->{AtomSymbol} = $AtomSymbol;
     160
     161   $AtomicNumber = PeriodicTable::GetElementAtomicNumber($AtomSymbol);
     162   $This->{AtomicNumber} = (defined $AtomicNumber) ? $AtomicNumber : 0;
     163
     164   return $This;
     165 }
     166
     167 # Setup atom symbol and atomic number for the element...
     168 sub SetAtomicNumber {
     169   my($This, $AtomicNumber) = @_;
     170   my($AtomSymbol);
     171
     172   $AtomSymbol = PeriodicTable::GetElementAtomSymbol($AtomicNumber);
     173   if (!defined $AtomSymbol) {
     174     carp "Warning: ${ClassName}->SetAtomicNumber: Didn't set atomic number: Invalid atomic number, $AtomicNumber, specified...";
     175     return;
     176   }
     177   $This->{AtomicNumber} = $AtomicNumber;
     178   $This->{AtomSymbol} = $AtomSymbol;
     179
     180   return $This;
     181 }
     182
     183 # Set atom as stereo center...
     184 #
     185 sub SetStereoCenter {
     186   my($This, $StereoCenter) = @_;
     187
     188   $This->SetProperty('StereoCenter', $StereoCenter);
     189
     190   return $This;
     191 }
     192
     193 # Is it a stereo center?
     194 #
     195 sub IsStereoCenter {
     196   my($This) = @_;
     197   my($StereoCenter);
     198
     199   $StereoCenter = $This->GetProperty('StereoCenter');
     200
     201   return (defined($StereoCenter) && $StereoCenter) ? 1 : 0;
     202 }
     203
     204 # Set atom stereochemistry.
     205 #
     206 # Supported values are: R, S.
     207 #
     208 # Notes:
     209 #
     210 # . After the ligands around a central stereocenter has been ranked using CIP priority scheme and
     211 # the lowest ranked ligand lies behind the center atom, then R and S values correspond to:
     212 #
     213 # R: Clockwise arrangement of remaining ligands around the central atom going from highest to lowest ranked ligand
     214 # S: CounterClockwise arrangement of remaining ligands around the central atom going from highest to lowest ranked ligand
     215 #
     216 # . Assignment of any other arbitray values besides R and S is also allowed; however, a warning is printed.
     217 #
     218 sub SetStereochemistry {
     219   my($This, $Stereochemistry) = @_;
     220
     221   if ($Stereochemistry !~ /^(R|S)$/i) {
     222     carp "Warning: ${ClassName}->SetStereochemistry: Assigning non-supported Stereochemistry value of $Stereochemistry. Supported values: R, S...";
     223   }
     224
     225   $This->SetProperty('StereoCenter', 1);
     226   $This->SetProperty('Stereochemistry', $Stereochemistry);
     227
     228   return $This;
     229 }
     230
     231 # Setup mass number for atom...
     232 sub SetMassNumber {
     233   my($This, $MassNumber) = @_;
     234   my($AtomicNumber, $AtomSymbol);
     235
     236   $AtomicNumber = $This->{AtomicNumber};
     237   $AtomSymbol = $This->{AtomSymbol};
     238   if (!$AtomicNumber) {
     239     carp "Warning: ${ClassName}->SetMassNumber: Didn't set mass number: Non standard atom with atomic number, $AtomicNumber, and atomic symbol, $AtomSymbol...";
     240     return;
     241   }
     242   if (!PeriodicTable::IsElementNaturalIsotopeMassNumber($AtomicNumber, $MassNumber)) {
     243     carp "Warning: ${ClassName}->SetMassNumber: Unknown mass number, $MassNumber, specified for atom with atomic number, $AtomicNumber, and atomic symbol, $AtomSymbol. Don't forget to Set ExactMass property explicitly; otherwise, GetExactMass method would return mass of most abundant isotope...";
     244   }
     245   $This->SetProperty('MassNumber', $MassNumber);
     246
     247   return $This;
     248 }
     249
     250 # Get mass number...
     251 #
     252 sub GetMassNumber {
     253   my($This) = @_;
     254
     255   # Is mass number explicity set?
     256   if ($This->HasProperty('MassNumber')) {
     257     return $This->GetProperty('MassNumber');
     258   }
     259
     260   # Is it an element symbol?
     261   my($AtomicNumber) = $This->{AtomicNumber};
     262   if (!$AtomicNumber) {
     263     return 0;
     264   }
     265
     266   # Return most abundant mass number...
     267   return PeriodicTable::GetElementMostAbundantNaturalIsotopeMassNumber($AtomicNumber);
     268 }
     269
     270 # Get atomic weight:
     271 #   . Explicitly set by the caller
     272 #   . Using atomic number
     273 #
     274 sub GetAtomicWeight {
     275   my($This) = @_;
     276
     277   # Is atomic weight explicity set?
     278   if ($This->HasProperty('AtomicWeight')) {
     279     return $This->GetProperty('AtomicWeight');
     280   }
     281
     282   # Is it an element symbol?
     283   my($AtomicNumber) = $This->{AtomicNumber};
     284   if (!$AtomicNumber) {
     285     return 0;
     286   }
     287
     288   # Return its atomic weight...
     289   return PeriodicTable::GetElementAtomicWeight($AtomicNumber);
     290 }
     291
     292 # Get exact mass weight:
     293 #   . Explicitly set by the caller
     294 #   . Using atomic number and mass number explicity set by the caller
     295 #   . Using atomic number and most abundant isotope
     296 #
     297 sub GetExactMass {
     298   my($This) = @_;
     299
     300   # Is exact mass explicity set?
     301   if ($This->HasProperty('ExactMass')) {
     302     return $This->GetProperty('ExactMass');
     303   }
     304
     305   # Is it an element symbol?
     306   my($AtomicNumber) = $This->{AtomicNumber};
     307   if (!$AtomicNumber) {
     308     return 0;
     309   }
     310
     311   # Is mass number explicitly set?
     312   if ($This->HasProperty('MassNumber')) {
     313     my($MassNumber) = $This->GetProperty('MassNumber');
     314     if (PeriodicTable::IsElementNaturalIsotopeMassNumber($AtomicNumber, $MassNumber)) {
     315       return PeriodicTable::GetElementNaturalIsotopeMass($AtomicNumber, $MassNumber);
     316     }
     317   }
     318
     319   # Return most abundant isotope mass...
     320   return PeriodicTable::GetElementMostAbundantNaturalIsotopeMass($AtomicNumber);
     321 }
     322
     323 # Get formal charge:
     324 #   . Explicitly set by the caller
     325 #   . Or return zero insetad of undef
     326 #
     327 sub GetFormalCharge {
     328   my($This) = @_;
     329   my($FormalCharge);
     330
     331   $FormalCharge = 0;
     332   if ($This->HasProperty('FormalCharge')) {
     333     $FormalCharge = $This->GetProperty('FormalCharge');
     334   }
     335
     336   return defined($FormalCharge) ? $FormalCharge : 0;
     337 }
     338
     339 # Get spin multiplicity:
     340 #   . Explicitly set by the caller
     341 #   . From FreeRadicalElectrons value explicitly set by the caller
     342 #   . Or return zero insetad of undef
     343 #
     344 sub GetSpinMultiplicity {
     345   my($This) = @_;
     346   my($SpinMultiplicity);
     347
     348   $SpinMultiplicity = 0;
     349   if ($This->HasProperty('SpinMultiplicity')) {
     350     $SpinMultiplicity = $This->GetProperty('SpinMultiplicity');
     351     return defined($SpinMultiplicity) ? $SpinMultiplicity : 0;
     352   }
     353
     354   if ($This->HasProperty('FreeRadicalElectrons')) {
     355     my($FreeRadicalElectrons);
     356     $FreeRadicalElectrons = $This->GetProperty('FreeRadicalElectrons');
     357
     358     SPINMULTIPLICITY: {
     359       if ($FreeRadicalElectrons == 1) { $SpinMultiplicity = 2; last SPINMULTIPLICITY;}
     360       if ($FreeRadicalElectrons == 2) { $SpinMultiplicity = 1; last SPINMULTIPLICITY;}
     361       carp "Warning: ${ClassName}->GetSpinMultiplicity: It's not possible to determine spin multiplicity from the specified free radical electrons value, $FreeRadicalElectrons. It has been set to 0...";
     362       $SpinMultiplicity = 0;
     363     }
     364   }
     365
     366   return $SpinMultiplicity;
     367 }
     368
     369 # Get number of free radical electrons:
     370 #   . Explicitly set by the caller
     371 #   . From SpinMultiplicity value explicitly set by the caller
     372 #   . Or return zero insetad of undef
     373 #
     374 # Notes:
     375 #  . For atoms with explicit assignment of SpinMultiplicity property values corresponding to
     376 #    Singlet (two unpaired electrons corresponding to one spin state), Doublet (free radical; an unpaired
     377 #    electron corresponding to two spin states), and Triplet (two unparied electrons corresponding to
     378 #    three spin states; divalent carbon atoms (carbenes)), FreeRadicalElectrons are calculated as follows:
     379 #
     380 #       SpinMultiplicity: Doublet(2); FreeRadicalElectrons: 1
     381 #       SpinMultiplicity: Singlet(1)/Triplet(3); FreeRadicalElectrons: 2
     382 #
     383 sub GetFreeRadicalElectrons {
     384   my($This) = @_;
     385   my($FreeRadicalElectrons);
     386
     387   $FreeRadicalElectrons = 0;
     388
     389   if ($This->HasProperty('FreeRadicalElectrons')) {
     390     $FreeRadicalElectrons = $This->GetProperty('FreeRadicalElectrons');
     391     return defined($FreeRadicalElectrons) ? $FreeRadicalElectrons : 0;
     392   }
     393
     394   if ($This->HasProperty('SpinMultiplicity')) {
     395     my($SpinMultiplicity);
     396     $SpinMultiplicity = $This->GetProperty('SpinMultiplicity');
     397
     398     SPINMULTIPLICITY: {
     399       if ($SpinMultiplicity == 1) { $FreeRadicalElectrons = 2; last SPINMULTIPLICITY;}
     400       if ($SpinMultiplicity == 2) { $FreeRadicalElectrons = 1; last SPINMULTIPLICITY;}
     401       if ($SpinMultiplicity == 3) { $FreeRadicalElectrons = 2; last SPINMULTIPLICITY;}
     402       carp "Warning: ${ClassName}->GetFreeRadicalElectrons: It's not possible to determine free radical electrons from the specified spin multiplicity value, $FreeRadicalElectrons. It has been set to 0...";
     403       $FreeRadicalElectrons = 0;
     404     }
     405   }
     406
     407   return $FreeRadicalElectrons;
     408 }
     409
     410 # Set atom coordinates using:
     411 # . An array reference with three values
     412 # . An array containg three values
     413 # . A 3D vector
     414 #
     415 sub SetXYZ {
     416   my($This, @Values) = @_;
     417
     418   if (!@Values) {
     419     carp "Warning: ${ClassName}->SetXYZ: No values specified...";
     420     return;
     421   }
     422
     423   $This->{XYZ}->SetXYZ(@Values);
     424   return $This;
     425 }
     426
     427 # Set X value...
     428 sub SetX {
     429   my($This, $Value) = @_;
     430
     431   if (!defined $Value) {
     432     carp "Warning: ${ClassName}->SetX: Undefined X value...";
     433     return;
     434   }
     435   $This->{XYZ}->SetX($Value);
     436   return $This;
     437 }
     438
     439 # Set Y value...
     440 sub SetY {
     441   my($This, $Value) = @_;
     442
     443   if (!defined $Value) {
     444     carp "Warning: ${ClassName}->SetY: Undefined Y value...";
     445     return;
     446   }
     447   $This->{XYZ}->SetY($Value);
     448   return $This;
     449 }
     450
     451 # Set Z value...
     452 sub SetZ {
     453   my($This, $Value) = @_;
     454
     455   if (!defined $Value) {
     456     carp "Warning: ${ClassName}->SetZ: Undefined Z value...";
     457     return;
     458   }
     459   $This->{XYZ}->SetZ($Value);
     460   return $This;
     461 }
     462
     463 # Return XYZ as:
     464 # . Reference to an array
     465 # . An array
     466 #
     467 sub GetXYZ {
     468   my($This) = @_;
     469
     470   return $This->{XYZ}->GetXYZ();
     471 }
     472
     473 # Return XYZ as a vector object...
     474 #
     475 sub GetXYZVector {
     476   my($This) = @_;
     477
     478   return $This->{XYZ};
     479 }
     480
     481 # Get X value...
     482 sub GetX {
     483   my($This) = @_;
     484
     485   return $This->{XYZ}->GetX();
     486 }
     487
     488 # Get Y value...
     489 sub GetY {
     490   my($This) = @_;
     491
     492   return $This->{XYZ}->GetY();
     493 }
     494
     495 # Get Z value...
     496 sub GetZ {
     497   my($This) = @_;
     498
     499   return $This->{XYZ}->GetZ();
     500 }
     501
     502 # Delete atom...
     503 sub DeleteAtom {
     504   my($This) = @_;
     505
     506   # Is this atom in a molecule?
     507   if (!$This->HasProperty('Molecule')) {
     508     # Nothing to do...
     509     return $This;
     510   }
     511   my($Molecule) = $This->GetProperty('Molecule');
     512
     513   return $Molecule->_DeleteAtom($This);
     514 }
     515
     516 # Get atom neighbor objects as array. In scalar conetxt, return number of neighbors...
     517 sub GetNeighbors {
     518   my($This, @ExcludeNeighbors) = @_;
     519
     520   # Is this atom in a molecule?
     521   if (!$This->HasProperty('Molecule')) {
     522     return undef;
     523   }
     524   my($Molecule) = $This->GetProperty('Molecule');
     525
     526   if (@ExcludeNeighbors) {
     527     return $This->_GetAtomNeighbors(@ExcludeNeighbors);
     528   }
     529   else {
     530     return $This->_GetAtomNeighbors();
     531   }
     532 }
     533
     534 # Get atom neighbor objects as array. In scalar conetxt, return number of neighbors...
     535 sub _GetAtomNeighbors {
     536   my($This, @ExcludeNeighbors) = @_;
     537   my($Molecule) = $This->GetProperty('Molecule');
     538
     539   if (!@ExcludeNeighbors) {
     540     return $Molecule->_GetAtomNeighbors($This);
     541   }
     542
     543   # Setup a map for neigbhors to exclude...
     544   my($ExcludeNeighbor, $ExcludeNeighborID, %ExcludeNeighborsIDsMap);
     545
     546   %ExcludeNeighborsIDsMap = ();
     547   for $ExcludeNeighbor (@ExcludeNeighbors) {
     548     $ExcludeNeighborID = $ExcludeNeighbor->GetID();
     549     $ExcludeNeighborsIDsMap{$ExcludeNeighborID} = $ExcludeNeighborID;
     550   }
     551
     552   # Generate a filtered neighbors list...
     553   my($Neighbor, $NeighborID, @FilteredAtomNeighbors);
     554   @FilteredAtomNeighbors = ();
     555   NEIGHBOR: for $Neighbor ($Molecule->_GetAtomNeighbors($This)) {
     556       $NeighborID = $Neighbor->GetID();
     557       if (exists $ExcludeNeighborsIDsMap{$NeighborID}) {
     558         next NEIGHBOR;
     559       }
     560     push @FilteredAtomNeighbors, $Neighbor;
     561   }
     562
     563   return wantarray ? @FilteredAtomNeighbors : scalar @FilteredAtomNeighbors;
     564 }
     565
     566 # Get specific atom neighbor objects as array. In scalar conetxt, return number of neighbors.
     567 #
     568 # Notes:
     569 #   . AtomSpecification correspond to any valid AtomicInvariant based atomic specifications
     570 #     as implemented in DoesAtomNeighborhoodMatch method.
     571 #   . Multiple atom specifications can be used in a string delimited by comma.
     572 #
     573 sub GetNeighborsUsingAtomSpecification {
     574   my($This, $AtomSpecification, @ExcludeNeighbors) = @_;
     575   my(@AtomNeighbors);
     576
     577   @AtomNeighbors = ();
     578   @AtomNeighbors = $This->GetNeighbors(@ExcludeNeighbors);
     579
     580   # Does atom has any neighbors and do they need to be filtered?
     581   if (!(@AtomNeighbors && defined($AtomSpecification) && $AtomSpecification)) {
     582     return wantarray ? @AtomNeighbors : scalar @AtomNeighbors;
     583   }
     584
     585   # Filter neighbors using atom specification...
     586   my($AtomNeighbor, @FilteredAtomNeighbors);
     587
     588   @FilteredAtomNeighbors = ();
     589   NEIGHBOR: for $AtomNeighbor (@AtomNeighbors) {
     590     if (!$AtomNeighbor->_DoesAtomSpecificationMatch($AtomSpecification)) {
     591       next NEIGHBOR;
     592     }
     593     push @FilteredAtomNeighbors, $AtomNeighbor;
     594   }
     595
     596   return wantarray ? @FilteredAtomNeighbors : scalar @FilteredAtomNeighbors;
     597 }
     598
     599
     600 # Get non-hydrogen atom neighbor objects as array. In scalar context, return number of neighbors...
     601 sub GetHeavyAtomNeighbors {
     602   my($This) = @_;
     603
     604   return $This->GetNonHydrogenAtomNeighbors();
     605 }
     606
     607 # Get non-hydrogen atom neighbor objects as array. In scalar context, return number of neighbors...
     608 sub GetNonHydrogenAtomNeighbors {
     609   my($This) = @_;
     610
     611   # Is this atom in a molecule?
     612   if (!$This->HasProperty('Molecule')) {
     613     return undef;
     614   }
     615   my($NonHydrogenAtomsOnly, $HydrogenAtomsOnly) = (1, 0);
     616
     617   return $This->_GetFilteredAtomNeighbors($NonHydrogenAtomsOnly, $HydrogenAtomsOnly);
     618 }
     619
     620 # Get hydrogen atom neighbor objects as array. In scalar context, return numbe of neighbors...
     621 sub GetHydrogenAtomNeighbors {
     622   my($This) = @_;
     623
     624   # Is this atom in a molecule?
     625   if (!$This->HasProperty('Molecule')) {
     626     return undef;
     627   }
     628   my($NonHydrogenAtomsOnly, $HydrogenAtomsOnly) = (0, 1);
     629
     630   return $This->_GetFilteredAtomNeighbors($NonHydrogenAtomsOnly, $HydrogenAtomsOnly);
     631 }
     632
     633 # Get non-hydrogen neighbor of hydrogen atom...
     634 #
     635 sub GetNonHydrogenNeighborOfHydrogenAtom {
     636   my($This) = @_;
     637
     638   # Is it Hydrogen?
     639   if (!$This->IsHydrogen()) {
     640     return undef;
     641   }
     642   my(@Neighbors);
     643
     644   @Neighbors = $This->GetNonHydrogenAtomNeighbors();
     645
     646   return (@Neighbors == 1) ? $Neighbors[0] : undef;
     647 }
     648
     649 # Get filtered atom atom neighbors
     650 sub _GetFilteredAtomNeighbors {
     651   my($This, $NonHydrogenAtomsOnly, $HydrogenAtomsOnly) = @_;
     652
     653   # Check flags...
     654   if (!defined $NonHydrogenAtomsOnly) {
     655     $NonHydrogenAtomsOnly = 0;
     656   }
     657   if (!defined $HydrogenAtomsOnly) {
     658     $HydrogenAtomsOnly = 0;
     659   }
     660   my($Neighbor, @FilteredAtomNeighbors);
     661
     662   @FilteredAtomNeighbors = ();
     663   NEIGHBOR: for $Neighbor ($This->GetNeighbors()) {
     664     if ($NonHydrogenAtomsOnly && $Neighbor->IsHydrogen()) {
     665       next NEIGHBOR;
     666     }
     667     if ($HydrogenAtomsOnly && (!$Neighbor->IsHydrogen())) {
     668       next NEIGHBOR;
     669     }
     670     push @FilteredAtomNeighbors, $Neighbor;
     671   }
     672
     673   return wantarray ? @FilteredAtomNeighbors : scalar @FilteredAtomNeighbors;
     674 }
     675
     676 # Get number of neighbors...
     677 #
     678 sub GetNumOfNeighbors {
     679   my($This) = @_;
     680   my($NumOfNeighbors);
     681
     682   $NumOfNeighbors = $This->GetNeighbors();
     683
     684   return (defined $NumOfNeighbors) ? $NumOfNeighbors : undef;
     685 }
     686
     687 # Get number of neighbors which are non-hydrogen atoms...
     688 sub GetNumOfHeavyAtomNeighbors {
     689   my($This) = @_;
     690
     691   return $This->GetNumOfNonHydrogenAtomNeighbors();
     692 }
     693
     694 # Get number of neighbors which are non-hydrogen atoms...
     695 sub GetNumOfNonHydrogenAtomNeighbors {
     696   my($This) = @_;
     697   my($NumOfNeighbors);
     698
     699   $NumOfNeighbors = $This->GetNonHydrogenAtomNeighbors();
     700
     701   return (defined $NumOfNeighbors) ? $NumOfNeighbors : undef;
     702 }
     703
     704 # Get number of neighbors which are hydrogen atoms...
     705 sub GetNumOfHydrogenAtomNeighbors {
     706   my($This) = @_;
     707   my($NumOfNeighbors);
     708
     709   $NumOfNeighbors = $This->GetHydrogenAtomNeighbors();
     710
     711   return (defined $NumOfNeighbors) ? $NumOfNeighbors : undef;
     712 }
     713
     714 # Get bond objects as array. In scalar context, return number of bonds...
     715 sub GetBonds {
     716   my($This) = @_;
     717
     718   # Is this atom in a molecule?
     719   if (!$This->HasProperty('Molecule')) {
     720     return undef;
     721   }
     722   my($Molecule) = $This->GetProperty('Molecule');
     723
     724   return $Molecule->_GetAtomBonds($This);
     725 }
     726
     727 # Get bond to specified atom...
     728 sub GetBondToAtom {
     729   my($This, $Other) = @_;
     730
     731   # Is this atom in a molecule?
     732   if (!$This->HasProperty('Molecule')) {
     733     return undef;
     734   }
     735   my($Molecule) = $This->GetProperty('Molecule');
     736
     737   return $Molecule->_GetBondToAtom($This, $Other);
     738 }
     739
     740 # It it bonded to a specified atom?
     741 sub IsBondedToAtom {
     742   my($This, $Other) = @_;
     743
     744   # Is this atom in a molecule?
     745   if (!$This->HasProperty('Molecule')) {
     746     return undef;
     747   }
     748   my($Molecule) = $This->GetProperty('Molecule');
     749
     750   return $Molecule->_IsBondedToAtom($This, $Other);
     751 }
     752
     753 # Get bond objects to non-hydrogen atoms as array. In scalar context, return number of bonds...
     754 sub GetBondsToHeavyAtoms {
     755   my($This) = @_;
     756
     757   return $This->GetBondsToNonHydrogenAtoms();
     758 }
     759
     760 # Get bond objects to non-hydrogen atoms as array. In scalar context, return number of bonds...
     761 sub GetBondsToNonHydrogenAtoms {
     762   my($This) = @_;
     763
     764   # Is this atom in a molecule?
     765   if (!$This->HasProperty('Molecule')) {
     766     return undef;
     767   }
     768   my($BondsToNonHydrogenAtomsOnly, $BondsToHydrogenAtomsOnly) = (1, 0);
     769
     770   return $This->_GetFilteredBonds($BondsToNonHydrogenAtomsOnly, $BondsToHydrogenAtomsOnly);
     771 }
     772
     773 # Get bond objects to hydrogen atoms as array. In scalar context, return number of bonds...
     774 sub GetBondsToHydrogenAtoms {
     775   my($This) = @_;
     776
     777   # Is this atom in a molecule?
     778   if (!$This->HasProperty('Molecule')) {
     779     return undef;
     780   }
     781   my($BondsToNonHydrogenAtomsOnly, $BondsToHydrogenAtomsOnly) = (0, 1);
     782
     783   return $This->_GetFilteredBonds($BondsToNonHydrogenAtomsOnly, $BondsToHydrogenAtomsOnly);
     784 }
     785
     786 # Get filtered bonds...
     787 sub _GetFilteredBonds {
     788   my($This, $BondsToNonHydrogenAtomsOnly, $BondsToHydrogenAtomsOnly) = @_;
     789
     790   # Check flags...
     791   if (!defined $BondsToNonHydrogenAtomsOnly) {
     792     $BondsToNonHydrogenAtomsOnly = 0;
     793   }
     794   if (!defined $BondsToHydrogenAtomsOnly) {
     795     $BondsToHydrogenAtomsOnly = 0;
     796   }
     797
     798   my($Bond, $BondedAtom, @FilteredBonds);
     799
     800   @FilteredBonds = ();
     801   BOND: for $Bond ($This->GetBonds()) {
     802     $BondedAtom = $Bond->GetBondedAtom($This);
     803     if ($BondsToNonHydrogenAtomsOnly && $BondedAtom->IsHydrogen()) {
     804       next BOND;
     805     }
     806     if ($BondsToHydrogenAtomsOnly && (!$BondedAtom->IsHydrogen())) {
     807       next BOND;
     808     }
     809     push @FilteredBonds, $Bond;
     810   }
     811
     812   return wantarray ? @FilteredBonds : (scalar @FilteredBonds);
     813 }
     814
     815 # Get number of bonds...
     816 #
     817 sub GetNumOfBonds {
     818   my($This) = @_;
     819   my($NumOfBonds);
     820
     821   $NumOfBonds = $This->GetBonds();
     822
     823   return (defined $NumOfBonds) ? ($NumOfBonds) : undef;
     824 }
     825
     826 # Get number of bonds to non-hydrogen atoms...
     827 sub GetNumOfBondsToHeavyAtoms {
     828   my($This) = @_;
     829
     830   return $This->GetNumOfBondsToNonHydrogenAtoms();
     831 }
     832
     833 # Get number of bonds to non-hydrogen atoms...
     834 sub GetNumOfBondsToNonHydrogenAtoms {
     835   my($This) = @_;
     836   my($NumOfBonds);
     837
     838   $NumOfBonds = $This->GetBondsToNonHydrogenAtoms();
     839
     840   return (defined $NumOfBonds) ? ($NumOfBonds) : undef;
     841 }
     842
     843 # Get number of single bonds to heavy atoms...
     844 sub GetNumOfSingleBondsToHeavyAtoms {
     845   my($This) = @_;
     846
     847   return $This->GetNumOfSingleBondsToNonHydrogenAtoms();
     848 }
     849
     850 # Get number of single bonds to non-hydrogen atoms...
     851 sub GetNumOfSingleBondsToNonHydrogenAtoms {
     852   my($This) = @_;
     853
     854   # Is this atom in a molecule?
     855   if (!$This->HasProperty('Molecule')) {
     856     return undef;
     857   }
     858   return $This->_GetNumOfBondsWithSpecifiedBondOrderToNonHydrogenAtoms(1);
     859 }
     860
     861 # Get number of double bonds to heavy atoms...
     862 sub GetNumOfDoubleBondsToHeavyAtoms {
     863   my($This) = @_;
     864
     865   return $This->GetNumOfDoubleBondsToNonHydrogenAtoms();
     866 }
     867
     868 # Get number of double bonds to non-hydrogen atoms...
     869 sub GetNumOfDoubleBondsToNonHydrogenAtoms {
     870   my($This) = @_;
     871
     872   # Is this atom in a molecule?
     873   if (!$This->HasProperty('Molecule')) {
     874     return undef;
     875   }
     876   return $This->_GetNumOfBondsWithSpecifiedBondOrderToNonHydrogenAtoms(2);
     877 }
     878
     879 # Get number of triple bonds to heavy atoms...
     880 sub GetNumOfTripleBondsToHeavyAtoms {
     881   my($This) = @_;
     882
     883   return $This->GetNumOfTripleBondsToNonHydrogenAtoms();
     884 }
     885
     886 # Get number of triple bonds to non-hydrogen atoms...
     887 sub GetNumOfTripleBondsToNonHydrogenAtoms {
     888   my($This) = @_;
     889
     890   # Is this atom in a molecule?
     891   if (!$This->HasProperty('Molecule')) {
     892     return undef;
     893   }
     894   return $This->_GetNumOfBondsWithSpecifiedBondOrderToNonHydrogenAtoms(3);
     895 }
     896
     897 # Get number of bonds of specified bond order to non-hydrogen atoms...
     898 sub _GetNumOfBondsWithSpecifiedBondOrderToNonHydrogenAtoms {
     899   my($This, $SpecifiedBondOrder) = @_;
     900   my($NumOfBonds, $Bond, $BondOrder, @Bonds);
     901
     902   $NumOfBonds = 0;
     903   @Bonds = $This->GetBondsToNonHydrogenAtoms();
     904   for $Bond (@Bonds) {
     905     $BondOrder = $Bond->GetBondOrder();
     906     if ($SpecifiedBondOrder == $BondOrder) {
     907       $NumOfBonds++;
     908     }
     909   }
     910   return $NumOfBonds;
     911 }
     912
     913 # Get number of aromatic bonds to heavy atoms...
     914 sub GetNumOfAromaticBondsToHeavyAtoms {
     915   my($This) = @_;
     916
     917   return $This->GetNumOfAromaticBondsToNonHydrogenAtoms();
     918 }
     919
     920 # Get number of aromatic bonds to non-hydrogen atoms...
     921 sub GetNumOfAromaticBondsToNonHydrogenAtoms {
     922   my($This) = @_;
     923   my($NumOfBonds, $Bond, @Bonds);
     924
     925   # Is this atom in a molecule?
     926   if (!$This->HasProperty('Molecule')) {
     927     return undef;
     928   }
     929
     930   $NumOfBonds = 0;
     931   @Bonds = $This->GetBondsToNonHydrogenAtoms();
     932   for $Bond (@Bonds) {
     933     if ($Bond->IsAromatic()) { $NumOfBonds++; }
     934   }
     935   return $NumOfBonds;
     936 }
     937
     938 # Get number of different bond types to non-hydrogen atoms...
     939 #
     940 sub GetNumOfBondTypesToHeavyAtoms {
     941   my($This, $CountAromaticBonds) = @_;
     942
     943   return $This->GetNumOfBondTypesToNonHydrogenAtoms($CountAromaticBonds);
     944 }
     945
     946 # Get number of single, double, triple, and aromatic bonds from an atom to all other
     947 # non-hydrogen atoms. Value of CountAtomaticBonds parameter controls whether
     948 # number of aromatic bonds is returned; default is not to count aromatic bonds. During
     949 # counting of aromatic bonds, the bond marked aromatic is not included in the count
     950 # of other bond types.
     951 #
     952 sub GetNumOfBondTypesToNonHydrogenAtoms {
     953   my($This, $CountAromaticBonds) = @_;
     954   my($NumOfSingleBonds, $NumOfDoubleBonds, $NumOfTripleBonds, $NumOfAromaticBonds, $None, $Bond, @Bonds);
     955
     956   $CountAromaticBonds = defined($CountAromaticBonds) ? $CountAromaticBonds : 0;
     957
     958   ($NumOfSingleBonds, $NumOfDoubleBonds, $NumOfTripleBonds) = ('0') x 3;
     959   $NumOfAromaticBonds = $CountAromaticBonds ? 0 : undef;
     960
     961   # Is this atom in a molecule?
     962   if (!$This->HasProperty('Molecule')) {
     963     return ($NumOfSingleBonds, $NumOfDoubleBonds, $NumOfTripleBonds, $NumOfAromaticBonds);
     964   }
     965
     966   @Bonds = $This->GetBondsToNonHydrogenAtoms();
     967
     968   for $Bond (@Bonds) {
     969     BONDTYPE: {
     970       if ($CountAromaticBonds) {
     971         if ($Bond->IsAromatic()) { $NumOfAromaticBonds++; last BONDTYPE; }
     972       }
     973       if ($Bond->IsSingle()) { $NumOfSingleBonds++; last BONDTYPE; }
     974       if ($Bond->IsDouble()) { $NumOfDoubleBonds++; last BONDTYPE; }
     975       if ($Bond->IsTriple()) { $NumOfTripleBonds++; last BONDTYPE; }
     976       $None = 1;
     977     }
     978   }
     979   return ($NumOfSingleBonds, $NumOfDoubleBonds, $NumOfTripleBonds, $NumOfAromaticBonds);
     980 }
     981
     982 # Get number of sigma and pi bonds to heavy atoms...
     983 #
     984 sub GetNumOfSigmaAndPiBondsToHeavyAtoms {
     985   my($This) = @_;
     986
     987   return $This->GetNumOfSigmaAndPiBondsToNonHydrogenAtoms();
     988 }
     989
     990 # Get number of sigma and pi bonds from an atom to all other non-hydrogen atoms.
     991 # Sigma and pi bonds are counted using the following methodology: a single bond
     992 # correspond to one sigma bond; a double bond contributes one to sigma bond count
     993 # and one to pi bond count; a triple bond contributes one to sigma bond count and
     994 # two to pi bond count.
     995 #
     996 sub GetNumOfSigmaAndPiBondsToNonHydrogenAtoms {
     997   my($This) = @_;
     998   my($NumOfSingleBonds, $NumOfDoubleBonds, $NumOfTripleBonds, $NumOfSigmaBonds, $NumOfPiBonds);
     999
    1000   ($NumOfSingleBonds, $NumOfDoubleBonds, $NumOfTripleBonds) = $This->GetNumOfBondTypesToNonHydrogenAtoms();
    1001
    1002   $NumOfSigmaBonds = $NumOfSingleBonds + $NumOfDoubleBonds + $NumOfTripleBonds;
    1003   $NumOfPiBonds = $NumOfDoubleBonds + 2*$NumOfTripleBonds;
    1004
    1005   return ($NumOfSigmaBonds, $NumOfPiBonds);
    1006 }
    1007
    1008 # Get information related to atoms for all heavy atoms attached to an atom..
    1009 #
    1010 sub GetHeavyAtomNeighborsAtomInformation {
    1011   my($This) = @_;
    1012
    1013   return $This->GetNonHydrogenAtomNeighborsAtomInformation();
    1014 }
    1015
    1016 # Get information related to atoms for all non-hydrogen atoms attached to an atom..
    1017 #
    1018 # The following values are returned:
    1019 #   . Number of non-hydrogen atom neighbors
    1020 #   . A reference to an array containing atom objects correpsonding to non-hydrogen
    1021 #     atom neighbors
    1022 #   . Number of different types of non-hydrogen atom neighbors
    1023 #   . A reference to a hash containing atom symbol as key with value corresponding
    1024 #     to its count for non-hydrogen atom neighbors
    1025 #
    1026 sub GetNonHydrogenAtomNeighborsAtomInformation {
    1027   my($This) = @_;
    1028
    1029   # Is this atom in a molecule?
    1030   if (!$This->HasProperty('Molecule')) {
    1031     return (undef, undef, undef, undef);
    1032   }
    1033   my($AtomSymbol, $AtomNeighbor, $NumOfAtomNeighbors, $NumOfAtomNeighborsType, @AtomNeighbors, %AtomNeighborsTypeMap);
    1034
    1035   $NumOfAtomNeighbors = 0; @AtomNeighbors = ();
    1036   $NumOfAtomNeighborsType = 0; %AtomNeighborsTypeMap = ();
    1037
    1038   @AtomNeighbors = $This->GetNonHydrogenAtomNeighbors();
    1039   $NumOfAtomNeighbors = scalar @AtomNeighbors;
    1040
    1041   for $AtomNeighbor (@AtomNeighbors) {
    1042     $AtomSymbol = $AtomNeighbor->{AtomSymbol};
    1043     if (exists $AtomNeighborsTypeMap{$AtomSymbol}) {
    1044       $AtomNeighborsTypeMap{$AtomSymbol} += 1;
    1045     }
    1046     else {
    1047       $AtomNeighborsTypeMap{$AtomSymbol} = 1;
    1048       $NumOfAtomNeighborsType++;
    1049     }
    1050   }
    1051
    1052   return ($NumOfAtomNeighbors, \@AtomNeighbors, $NumOfAtomNeighborsType, \%AtomNeighborsTypeMap);
    1053 }
    1054
    1055 # Get information related to bonds for all heavy atoms attached to an atom..
    1056 #
    1057 sub GetHeavyAtomNeighborsBondformation {
    1058   my($This) = @_;
    1059
    1060   return $This->GetNonHydrogenAtomNeighborsBondInformation();
    1061 }
    1062
    1063 # Get information related to bonds for all non-hydrogen atoms attached to an atom..
    1064 #
    1065 # The following values are returned:
    1066 #   . Number of bonds to non-hydrogen atom neighbors
    1067 #   . A reference to an array containing bond objects correpsonding to non-hydrogen
    1068 #     atom neighbors
    1069 #   . A reference to a hash containing bond type as key with value corresponding
    1070 #     to its count for non-hydrogen atom neighbors. Bond types are: Single, Double or Triple
    1071 #   . A reference to a hash containing atom symbol as key pointing to bond type as second
    1072 #     key with values correponding to count of bond types for atom symbol for non-hydrogen
    1073 #     atom neighbors
    1074 #   . A reference to a hash containing atom symbol as key pointing to bond type as second
    1075 #     key with values correponding to atom objects array involved in corresponding bond type for
    1076 #     atom symbol for non-hydrogen atom neighbors
    1077 #
    1078 sub GetNonHydrogenAtomNeighborsBondInformation {
    1079   my($This) = @_;
    1080
    1081   # Is this atom in a molecule?
    1082   if (!$This->HasProperty('Molecule')) {
    1083     return (undef, undef, undef, undef, undef);
    1084   }
    1085   my($BondedAtom, $BondedAtomSymbol, $BondType, $None, $Bond, $NumOfBonds, @Bonds, %BondTypeCountMap, %AtomsBondTypesCountMap, %AtomsBondTypeAtomsMap);
    1086
    1087   $NumOfBonds = 0; @Bonds = ();
    1088   %BondTypeCountMap = ();
    1089   %AtomsBondTypesCountMap = (); %AtomsBondTypeAtomsMap = ();
    1090
    1091   $BondTypeCountMap{Single} = 0;
    1092   $BondTypeCountMap{Double} = 0;
    1093   $BondTypeCountMap{Triple} = 0;
    1094
    1095   @Bonds = $This->GetBondsToNonHydrogenAtoms();
    1096   $NumOfBonds = scalar @Bonds;
    1097
    1098   BOND: for $Bond (@Bonds) {
    1099     $BondType = $Bond->IsSingle() ? "Single" : ($Bond->IsDouble() ? "Double" : ($Bond->IsTriple() ? "Triple" : ""));
    1100     if (!$BondType) {
    1101       next BOND;
    1102     }
    1103
    1104     # Track bond types...
    1105     if (exists $BondTypeCountMap{$BondType}) {
    1106       $BondTypeCountMap{$BondType} += 1;
    1107     }
    1108     else {
    1109       $BondTypeCountMap{$BondType} = 1;
    1110     }
    1111
    1112     $BondedAtom = $Bond->GetBondedAtom($This);
    1113     $BondedAtomSymbol = $BondedAtom->{AtomSymbol};
    1114
    1115     # Track bond types count for atom types involved in specific bond types...
    1116     if (!exists $AtomsBondTypesCountMap{$BondedAtomSymbol}) {
    1117       %{$AtomsBondTypesCountMap{$BondedAtomSymbol}} = ();
    1118     }
    1119     if (exists $AtomsBondTypesCountMap{$BondedAtomSymbol}{$BondType}) {
    1120       $AtomsBondTypesCountMap{$BondedAtomSymbol}{$BondType} += 1;
    1121     }
    1122     else {
    1123       $AtomsBondTypesCountMap{$BondedAtomSymbol}{$BondType} = 1;
    1124     }
    1125
    1126     # Track atoms involved in specific bond types for specific atom types...
    1127     if (!exists $AtomsBondTypeAtomsMap{$BondedAtomSymbol}) {
    1128       %{$AtomsBondTypeAtomsMap{$BondedAtomSymbol}} = ();
    1129     }
    1130     if (!exists $AtomsBondTypeAtomsMap{$BondedAtomSymbol}{$BondType}) {
    1131       @{$AtomsBondTypeAtomsMap{$BondedAtomSymbol}{$BondType}} = ();
    1132     }
    1133     push @{$AtomsBondTypeAtomsMap{$BondedAtomSymbol}{$BondType}}, $BondedAtom;
    1134   }
    1135
    1136   return ($NumOfBonds, \@Bonds, \%BondTypeCountMap, \%AtomsBondTypesCountMap, \%AtomsBondTypeAtomsMap);
    1137 }
    1138
    1139 # Get number of bonds to hydrogen atoms...
    1140 sub GetNumOfBondsToHydrogenAtoms {
    1141   my($This) = @_;
    1142   my($NumOfBonds);
    1143
    1144   $NumOfBonds = $This->GetBondsToHydrogenAtoms();
    1145
    1146   return (defined $NumOfBonds) ? ($NumOfBonds) : undef;
    1147 }
    1148
    1149 # Get sum of bond orders to all bonded atoms...
    1150 #
    1151 sub GetSumOfBondOrders {
    1152   my($This) = @_;
    1153
    1154   # Is this atom in a molecule?
    1155   if (!$This->HasProperty('Molecule')) {
    1156     return undef;
    1157   }
    1158
    1159   return $This->_GetSumOfBondOrders();
    1160 }
    1161
    1162 # Get sum of bond orders to non-hydrogen atoms only...
    1163 #
    1164 sub GetSumOfBondOrdersToHeavyAtoms {
    1165   my($This) = @_;
    1166
    1167   return $This->GetSumOfBondOrdersToNonHydrogenAtoms();
    1168 }
    1169
    1170 # Get sum of bond orders to non-hydrogen atoms only...
    1171 #
    1172 sub GetSumOfBondOrdersToNonHydrogenAtoms {
    1173   my($This) = @_;
    1174
    1175   # Is this atom in a molecule?
    1176   if (!$This->HasProperty('Molecule')) {
    1177     return undef;
    1178   }
    1179   my($ToNonHydrogenAtomsOnly, $ToHydrogenAtomsOnly) = (1, 0);
    1180
    1181   return $This->_GetSumOfBondOrders($ToNonHydrogenAtomsOnly, $ToHydrogenAtomsOnly);
    1182 }
    1183
    1184 # Get sum of bond orders to hydrogen atoms only...
    1185 #
    1186 sub GetSumOfBondOrdersToHydrogenAtoms {
    1187   my($This) = @_;
    1188
    1189   # Is this atom in a molecule?
    1190   if (!$This->HasProperty('Molecule')) {
    1191     return undef;
    1192   }
    1193   my($ToNonHydrogenAtomsOnly, $ToHydrogenAtomsOnly) = (0, 1);
    1194
    1195   return $This->_GetSumOfBondOrders($ToNonHydrogenAtomsOnly, $ToHydrogenAtomsOnly);
    1196 }
    1197
    1198 # Get sum of bond orders to all bonded atoms,  non-hydrogen or hydrogen bonded atoms...
    1199 #
    1200 sub _GetSumOfBondOrders {
    1201   my($This, $ToNonHydrogenAtomsOnly, $ToHydrogenAtomsOnly) = @_;
    1202
    1203   # Check flags...
    1204   if (!defined $ToNonHydrogenAtomsOnly) {
    1205     $ToNonHydrogenAtomsOnly = 0;
    1206   }
    1207   if (!defined $ToHydrogenAtomsOnly) {
    1208     $ToHydrogenAtomsOnly = 0;
    1209   }
    1210   my($Bond, $SumOfBondOrders, @Bonds);
    1211   @Bonds = ();
    1212
    1213   if ($ToNonHydrogenAtomsOnly) {
    1214     @Bonds = $This->GetBondsToNonHydrogenAtoms();
    1215   }
    1216   elsif ($ToHydrogenAtomsOnly) {
    1217     @Bonds = $This->GetBondsToHydrogenAtoms();
    1218   }
    1219   else {
    1220     # All bonds...
    1221     @Bonds = $This->GetBonds();
    1222   }
    1223
    1224   $SumOfBondOrders = 0;
    1225   for $Bond (@Bonds) {
    1226     $SumOfBondOrders += $Bond->GetBondOrder();
    1227   }
    1228
    1229   if ($SumOfBondOrders =~ /\./) {
    1230     #
    1231     # Change any fractional bond order to next largest integer...
    1232     #
    1233     # As long as aromatic bond orders in a ring are correctly using using 4n + 2 Huckel rule
    1234     # (BondOrder: 1.5) or explicity set as Kekule  bonds (alternate single/double),
    1235     # SumOfBondOrders should add up to an integer.
    1236     #
    1237     $SumOfBondOrders = ceil($SumOfBondOrders);
    1238   }
    1239
    1240   return $SumOfBondOrders;
    1241 }
    1242
    1243 # Get largest bond order to any bonded atoms...
    1244 #
    1245 sub GetLargestBondOrder {
    1246   my($This) = @_;
    1247
    1248   # Is this atom in a molecule?
    1249   if (!$This->HasProperty('Molecule')) {
    1250     return undef;
    1251   }
    1252
    1253   return $This->_GetLargestBondOrder();
    1254 }
    1255
    1256 # Get largest bond order to bonded non-hydrogen atoms...
    1257 #
    1258 sub GetLargestBondOrderToHeavyAtoms {
    1259   my($This) = @_;
    1260
    1261   return $This->GetLargestBondOrderToNonHydrogenAtoms();
    1262 }
    1263
    1264 # Get largest bond order to bonded non-hydrogen atoms...
    1265 #
    1266 sub GetLargestBondOrderToNonHydrogenAtoms {
    1267   my($This) = @_;
    1268
    1269   # Is this atom in a molecule?
    1270   if (!$This->HasProperty('Molecule')) {
    1271     return undef;
    1272   }
    1273
    1274   my($ToNonHydrogenAtomsOnly) = (1);
    1275
    1276   return $This->_GetLargestBondOrder($ToNonHydrogenAtomsOnly);
    1277 }
    1278
    1279 # Get largest bond order to all bonded atoms, non-hydrogen or hydrogen bonded atoms...
    1280 #
    1281 sub _GetLargestBondOrder {
    1282   my($This, $ToNonHydrogenAtomsOnly, $ToHydrogenAtomsOnly) = @_;
    1283
    1284   # Check flags...
    1285   if (!defined $ToNonHydrogenAtomsOnly) {
    1286     $ToNonHydrogenAtomsOnly = 0;
    1287   }
    1288   if (!defined $ToHydrogenAtomsOnly) {
    1289     $ToHydrogenAtomsOnly = 0;
    1290   }
    1291   my($Bond, $LargestBondOrder, $BondOrder, @Bonds);
    1292   @Bonds = ();
    1293
    1294   if ($ToNonHydrogenAtomsOnly) {
    1295     @Bonds = $This->GetBondsToNonHydrogenAtoms();
    1296   }
    1297   elsif ($ToHydrogenAtomsOnly) {
    1298     @Bonds = $This->GetBondsToHydrogenAtoms();
    1299   }
    1300   else {
    1301     # All bonds...
    1302     @Bonds = $This->GetBonds();
    1303   }
    1304
    1305   $LargestBondOrder = 0;
    1306   for $Bond (@Bonds) {
    1307     $BondOrder = $Bond->GetBondOrder();
    1308     if ($BondOrder > $LargestBondOrder) {
    1309       $LargestBondOrder = $BondOrder;
    1310     }
    1311   }
    1312
    1313   return $LargestBondOrder;
    1314 }
    1315
    1316 # Get number of implicit hydrogen for atom...
    1317 #
    1318 sub GetImplicitHydrogens {
    1319   my($This) = @_;
    1320
    1321   return $This->GetNumOfImplicitHydrogens();
    1322 }
    1323
    1324 # Get number of implicit hydrogen for atom...
    1325 #
    1326 sub GetNumOfImplicitHydrogens {
    1327   my($This) = @_;
    1328
    1329   # Is this atom in a molecule?
    1330   if (!$This->HasProperty('Molecule')) {
    1331     return undef;
    1332   }
    1333
    1334   # Is ImplicitHydrogens property explicitly set?
    1335   if ($This->HasProperty('ImplicitHydrogens')) {
    1336     return $This->GetProperty('ImplicitHydrogens');
    1337   }
    1338
    1339   # Is it an element symbol?
    1340   if (!$This->{AtomicNumber}) {
    1341     return 0;
    1342   }
    1343
    1344   my($ImplicitHydrogens, $PotentialTotalValence, $SumOfBondOrders);
    1345
    1346   $ImplicitHydrogens = 0;
    1347   $SumOfBondOrders = $This->GetSumOfBondOrders();
    1348   $PotentialTotalValence = $This->GetPotentialTotalCommonValence();
    1349
    1350   if (defined($PotentialTotalValence) && defined($SumOfBondOrders)) {
    1351     # Subtract sum of bond orders to non-hydrogen and hydrogen atom neighbors...
    1352     $ImplicitHydrogens = $PotentialTotalValence - $SumOfBondOrders;
    1353   }
    1354
    1355   return $ImplicitHydrogens > 0 ? $ImplicitHydrogens : 0;
    1356 }
    1357
    1358 # Get number of bonds available to form additional bonds with heavy atoms, excluding
    1359 # any implicit bonds to hydrogens set using ImplicitHydrogens property.
    1360 #
    1361 # It's different from number of implicit or missing hydrogens, both of which are equivalent.
    1362 #
    1363 # For example, in a SMILES string, [nH] ring atom corresponds to an aromatic nitrogen.
    1364 # Although the hydrogen specified for n is treated internally as implicit hydrogen and shows
    1365 # up in missing hydrogen count, it's not available to participate in double bonds to additional
    1366 # heavy atoms.
    1367 #
    1368 sub GetNumOfBondsAvailableForHeavyAtoms {
    1369   my($This) = @_;
    1370
    1371   return $This->GetNumOfBondsAvailableForNonHydrogenAtoms();
    1372 }
    1373
    1374 # It's another name for GetNumOfBondsAvailableForHeavyAtoms
    1375 #
    1376 sub GetNumOfBondsAvailableForNonHydrogenAtoms {
    1377   my($This) = @_;
    1378   my($NumOfAvailableBonds, $PotentialTotalValence, $SumOfBondOrders);
    1379
    1380   $NumOfAvailableBonds = 0;
    1381
    1382   $SumOfBondOrders = $This->GetSumOfBondOrders();
    1383   $PotentialTotalValence = $This->GetPotentialTotalCommonValence();
    1384
    1385   if (defined($PotentialTotalValence) && defined($SumOfBondOrders)) {
    1386     # Subtract sum of bond orders to non-hydrogen and hydrogen atom neighbors...
    1387     $NumOfAvailableBonds = $PotentialTotalValence - $SumOfBondOrders;
    1388   }
    1389
    1390   if ($This->HasProperty('ImplicitHydrogens')) {
    1391     $NumOfAvailableBonds -= $This->GetProperty('ImplicitHydrogens');
    1392   }
    1393
    1394   return $NumOfAvailableBonds > 0 ? $NumOfAvailableBonds : 0;
    1395 }
    1396
    1397 # Disable setting of explicit hydrogens property...
    1398 sub SetExplicitHydrogens {
    1399   my($This, $Value) = @_;
    1400
    1401   carp "Warning: ${ClassName}->SetExplicitHydrogens: Setting of explicit hydrogens is not supported...";
    1402
    1403   return $This;
    1404 }
    1405
    1406 # Get number of explicit hydrogens for atom...
    1407 #
    1408 sub GetExplicitHydrogens {
    1409   my($This) = @_;
    1410
    1411   return $This->GetNumOfExplicitHydrogens();
    1412 }
    1413
    1414 # Get number of explicit hydrogens for atom...
    1415 #
    1416 sub GetNumOfExplicitHydrogens {
    1417   my($This) = @_;
    1418   my($HydrogenAtomNbrs);
    1419
    1420   # Is this atom in a molecule?
    1421   if (!$This->HasProperty('Molecule')) {
    1422     return undef;
    1423   }
    1424
    1425   $HydrogenAtomNbrs = $This->GetNumOfHydrogenAtomNeighbors();
    1426
    1427   return defined $HydrogenAtomNbrs ? $HydrogenAtomNbrs : 0;
    1428 }
    1429
    1430 # Get num of missing hydrogens...
    1431 #
    1432 sub GetMissingHydrogens {
    1433   my($This) = @_;
    1434
    1435   return $This->GetNumOfMissingHydrogens();
    1436 }
    1437
    1438 # Get num of missing hydrogens...
    1439 #
    1440 sub GetNumOfMissingHydrogens {
    1441   my($This) = @_;
    1442
    1443   # Is this atom in a molecule?
    1444   if (!$This->HasProperty('Molecule')) {
    1445     return undef;
    1446   }
    1447
    1448   return $This->GetNumOfImplicitHydrogens();
    1449 }
    1450
    1451 # Get total number of hydrogens...
    1452 #
    1453 sub GetHydrogens {
    1454   my($This) = @_;
    1455
    1456   return $This->GetNumOfHydrogens();
    1457 }
    1458
    1459 # Get total number of hydrogens...
    1460 #
    1461 sub GetNumOfHydrogens {
    1462   my($This) = @_;
    1463
    1464   # Is this atom in a molecule?
    1465   if (!$This->HasProperty('Molecule')) {
    1466     return undef;
    1467   }
    1468
    1469   return $This->GetNumOfImplicitHydrogens() + $This->GetNumOfExplicitHydrogens();
    1470 }
    1471
    1472 # Valence corresponds to the number of electrons used by an atom in bonding:
    1473 #
    1474 #   Valence = ValenceElectrons - ValenceFreeElectrons = BondingElectrons
    1475 #
    1476 # Single, double, triple bonds with bond orders of 1, 2 and 3 correspond to contribution of
    1477 # 1, 2, and 3 bonding electrons. So Valence can be computed using:
    1478 #
    1479 #   Valence = SumOfBondOrders + NumOfMissingHydrogens + FormalCharge
    1480 #
    1481 # where positive and negative values of FormalCharge increase and decrease the number
    1482 # of bonding electrons respectively.
    1483 #
    1484 # The current release of MayaChemTools supports the following three valence models, which
    1485 # are used during calculation of implicit hydrogens: MDLValenceModel, DaylightValenceModel,
    1486 # InternalValenceModel or MayaChemToolsValenceModel.
    1487 #
    1488 # Notes:
    1489 #  . This doesn't always corresponds to explicit valence.
    1490 #  . Missing hydrogens are included in the valence.
    1491 #  . For neutral molecules, valence and sum of bond order are equal.
    1492 #  . For molecules containing only single bonds, SumOfBondOrders and NumOfBonds are equal.
    1493 #  . Free radical electrons lead to the decrease in valence. For atoms with explicit assignment
    1494 #    of SpinMultiplicity property values corresponding to Singlet (two unparied electrons
    1495 #    corresponding to one spin state), Doublet (free radical; an unpaired electron corresponding
    1496 #    to two spin states), and Triplet (two unparied electrons corresponding to three spin states;
    1497 #    divalent carbon atoms (carbenes)), FreeRadicalElectrons are calculated as follows:
    1498 #
    1499 #       SpinMultiplicity: Doublet(2); FreeRadicalElectrons: 1 (one valence electron not available for bonding)
    1500 #       SpinMultiplicity: Singlet(1)/Triplet(3); FreeRadicalElectrons: 2 (two valence electrons not available for bonding)
    1501 #
    1502 sub GetValence {
    1503   my($This) = @_;
    1504
    1505   # Is this atom in a molecule?
    1506   if (!$This->HasProperty('Molecule')) {
    1507     return undef;
    1508   }
    1509
    1510   # Is Valence property explicitly set?
    1511   if ($This->HasProperty('Valence')) {
    1512     return $This->GetProperty('Valence');
    1513   }
    1514   my($Valence);
    1515
    1516   $Valence = $This->GetSumOfBondOrders() + $This->GetNumOfMissingHydrogens() + $This->GetFormalCharge();
    1517
    1518   return $Valence > 0 ? $Valence : 0;
    1519 }
    1520
    1521 # Get free non-bodning valence electrons left on atom after taking into account
    1522 # sum of bond orders, missing hydrogens and formal charged on the atom. Free
    1523 # radical electrons are included in the valence free electrons count by default.
    1524 #
    1525 # Valence corresponds to number of electrons used by atom in bonding:
    1526 #
    1527 #   Valence = ValenceElectrons - ValenceFreeElectrons
    1528 #
    1529 # Additionally, valence can also be calculated by:
    1530 #
    1531 #   Valence = SumOfBondOrders + NumOfMissingHydrogens + FormalCharge
    1532 #
    1533 # Valence and SumOfBondOrders are equal for neutral molecules.
    1534 #
    1535 # From two formulas for Valence described above, non-bonding free electrons
    1536 # left can be computed by:
    1537 #
    1538 #  ValenceFreeElectrons = ValenceElectrons - Valence
    1539 #                       = ValenceElectrons - SumOfBondOrders -
    1540 #                         NumOfMissingHydrogens - FormalCharge
    1541 #
    1542 # . Notes:
    1543 #    . Missing hydrogens are excluded from the valence free electrons.
    1544 #    . Any free radical electrons are considered part of the valence free electrons
    1545 #      by default.
    1546 #
    1547 # Examples:
    1548 #
    1549 # o NH3: ValenceFreeElectrons = 5 - 3 = 5 - 3 - 0 - 0 = 2
    1550 # o NH2: ValenceFreeElectrons = 5 - 3 = 5 - 2 - 1 - 0 = 2
    1551 # o NH4+; ValenceFreeElectrons = 5 - 5 = 5 - 4 - 0 - 1 = 0
    1552 # o NH3+; ValenceFreeElectrons = 5 - 5 = 5 - 3 - 1 - 1 = 0
    1553 # o C(=O)O- : ValenceFreeElectrons on O- = 6 - 0 = 6 - 1 - 0 - (-1) = 6
    1554 # o C(=O)O- : ValenceFreeElectrons on =O = 6 - 2 = 6 - 2 - 0 - 0 = 4
    1555 #
    1556 #
    1557 sub GetValenceFreeElectrons {
    1558   my($This, $ExcludeFreeRadicalElectrons) = @_;
    1559
    1560   # Is this atom in a molecule?
    1561   if (!$This->HasProperty('Molecule')) {
    1562     return undef;
    1563   }
    1564
    1565   # Is ValenceFreeElectrons property explicitly set?
    1566   if ($This->HasProperty('ValenceFreeElectrons')) {
    1567     return $This->GetProperty('ValenceFreeElectrons');
    1568   }
    1569
    1570   if (!$This->{AtomicNumber}) {
    1571     return 0;
    1572   }
    1573
    1574   my($ValenceFreeElectrons);
    1575
    1576   $ValenceFreeElectrons = $This->GetValenceElectrons() - $This->GetValence();
    1577   if ($ExcludeFreeRadicalElectrons) {
    1578     $ValenceFreeElectrons -= $This->GetFreeRadicalElectrons();
    1579   }
    1580
    1581   return $ValenceFreeElectrons > 0 ? $ValenceFreeElectrons : 0;
    1582 }
    1583
    1584 # Get potential total common valence for calculating the number of implicit hydrogens
    1585 # using the specified common valence model or default internal model for a molecule...
    1586 #
    1587 sub GetPotentialTotalCommonValence {
    1588   my($This) = @_;
    1589
    1590   # Is this atom in a molecule?
    1591   if (!$This->HasProperty('Molecule')) {
    1592     return undef;
    1593   }
    1594   my($PotentialTotalValence, $ValenceModel);
    1595
    1596   $PotentialTotalValence = 0;
    1597   $ValenceModel = $This->GetProperty('Molecule')->GetValenceModel();
    1598
    1599   VALENCEMODEL: {
    1600     if ($ValenceModel =~ /^MDLValenceModel$/i) {
    1601       $PotentialTotalValence = $This->_GetPotentialTotalCommonValenceUsingMDLValenceModel();
    1602       last VALENCEMODEL;
    1603     }
    1604     if ($ValenceModel =~ /^DaylightValenceModel$/i) {
    1605        $PotentialTotalValence = $This->_GetPotentialTotalCommonValenceUsingDaylightValenceModel();
    1606       last VALENCEMODEL;
    1607     }
    1608     if ($ValenceModel !~ /^(InternalValenceModel|MayaChemToolsValenceModel)$/i) {
    1609       carp "Warning: ${ClassName}->GetPotentialTotalCommonValence: The current release of MayaChemTools doesn't support the specified valence model $ValenceModel. Supported valence models: MDLValenceModel, DaylightValenceModel, InternalValenceModel. Using internal valence model...";
    1610     }
    1611     # Use internal valence model as the default valence model...
    1612     $PotentialTotalValence = $This->_GetPotentialTotalCommonValenceUsingInternalValenceModel();
    1613   }
    1614
    1615   return $PotentialTotalValence;
    1616 }
    1617
    1618 # Get potential total common valence using data for MDL valence model available in file,
    1619 # lib/data/MDLValenceModelData.csv, distributed with the package...
    1620 #
    1621 sub _GetPotentialTotalCommonValenceUsingMDLValenceModel {
    1622   my($This) = @_;
    1623
    1624   return $This->_GetPotentialTotalCommonValenceUsingValenceModelData(\%MDLValenceModelDataMap);
    1625
    1626 }
    1627
    1628 # Get potential total common valence using data for Daylight valence model available in file,
    1629 # lib/data/DaylightValenceModelData.csv, distributed with the release...
    1630 #
    1631 sub _GetPotentialTotalCommonValenceUsingDaylightValenceModel {
    1632   my($This) = @_;
    1633
    1634   return $This->_GetPotentialTotalCommonValenceUsingValenceModelData(\%DaylightValenceModelDataMap);
    1635 }
    1636
    1637 # Get potential total common valence using data for a specific valence model...
    1638 #
    1639 sub _GetPotentialTotalCommonValenceUsingValenceModelData {
    1640   my($This, $ValenceModelDataRef) = @_;
    1641   my($AtomicNumber, $FormalCharge);
    1642
    1643   $AtomicNumber = $This->{AtomicNumber};
    1644   if (!$AtomicNumber) {
    1645     return 0;
    1646   }
    1647
    1648   $FormalCharge = $This->GetFormalCharge();
    1649
    1650   # Is any valence model data available for atomic number and formal charge?
    1651   if (!exists $ValenceModelDataRef->{$AtomicNumber}) {
    1652     return 0;
    1653   }
    1654   if (!exists $ValenceModelDataRef->{$AtomicNumber}{$FormalCharge}) {
    1655     return 0;
    1656   }
    1657
    1658   my($PotentialTotalValence, $SumOfBondOrders, $CurrentEffectiveValence, $AvailableCommonValence);
    1659
    1660   $SumOfBondOrders = $This->GetSumOfBondOrders();
    1661   if (!defined $SumOfBondOrders) {
    1662     $SumOfBondOrders = 0;
    1663   }
    1664   $CurrentEffectiveValence = $SumOfBondOrders + $This->GetFreeRadicalElectrons();
    1665
    1666   $PotentialTotalValence = 0;
    1667   VALENCE: for $AvailableCommonValence (@{$ValenceModelDataRef->{$AtomicNumber}{$FormalCharge}{CommonValences}}) {
    1668       if ($CurrentEffectiveValence <= $AvailableCommonValence) {
    1669         $PotentialTotalValence = $AvailableCommonValence;
    1670         last VALENCE;
    1671       }
    1672   }
    1673
    1674   return $PotentialTotalValence;
    1675 }
    1676
    1677 #
    1678 # For elements with one one common valence, potential total common valence used
    1679 # during the calculation for number of implicit hydrogens during InternalValenceMode
    1680 # corresponds to:
    1681 #
    1682 #   CommonValence + FormalCharge - FreeRadicalElectrons
    1683 #
    1684 # For elements with multiple common valences, each common valence is used to
    1685 # calculate total potential common valence as shown above, and the first total potential
    1686 # common valence gerater than the sum of bond orderes is selected as the final total
    1687 # common valence.
    1688 #
    1689 # Group numbers > 14 - Group numbers 15 (N), 16 (O), 17 (F), 18 (He)
    1690 #
    1691 # Formal charge sign is not adjusted. Positive and negative values result in the
    1692 # increase and decrease of valence.
    1693 #
    1694 # Group 14 containing C, Si, Ge, Sn, Pb...
    1695 #
    1696 # Formal charge sign is reversed for positive values. Both positive and negative
    1697 # values result in the decrease of valence.
    1698 #
    1699 # Group 13 containing B, Al, Ga, In, Tl...
    1700 #
    1701 # Formal charge sign is always reversed. Positive and negative values result in the
    1702 # decrease and increase of valence.
    1703 #
    1704 # Groups 1 (H) through 12 (Zn)...
    1705 #
    1706 # Formal charge sign is reversed for positive values. Both positive and negative
    1707 # values result in the decrease of valence.
    1708 #
    1709 # Lanthanides and actinides...
    1710 #
    1711 # Formal charge sign is reversed for positive values. Both positive and negative
    1712 # values result in the decrease of valence.
    1713 #
    1714 # Notes:
    1715 #  . CommonValence and HighestCommonValence available from PeriodicTable module
    1716 #    are equivalent to most common and highest sum of bond orders for an element. For
    1717 #    neutral atoms involved only in single bonds, it corresponds to highest number of
    1718 #    allowed bonds for the atom.
    1719 #  . FormalCharge sign is reversed for electropositive elements with positive formal charge
    1720 #    during common valence calculations. Electropositive elements, metals and transition elements,
    1721 #    have usually plus formal charge and it leads to decrease in common valence; the negative
    1722 #    formal charge should result in the decrease of common valence.
    1723 #  . For carbon, both plus/minus formal charge cause decrease in common valence
    1724 #  . For elements on the right of carbon in periodic table, electronegative elements, plus formal
    1725 #    charge causes common valence to increase and minus formal charge cause it to decrease.
    1726 #
    1727 sub _GetPotentialTotalCommonValenceUsingInternalValenceModel {
    1728   my($This) = @_;
    1729   my($AtomicNumber, $CommonValences);
    1730
    1731   $AtomicNumber = $This->{AtomicNumber};
    1732   if (!$AtomicNumber) {
    1733     return 0;
    1734   }
    1735
    1736   $CommonValences = PeriodicTable::GetElementCommonValences($AtomicNumber);
    1737   if (!$CommonValences) {
    1738     return 0;
    1739   }
    1740
    1741   my($PotentialTotalValence, $AdjustedFormalCharge, $FreeRadicalElectrons, $SumOfBondOrders, $AvailableCommonValence, @AvailableCommonValences);
    1742
    1743   $AdjustedFormalCharge = $This->_GetFormalChargeAdjustedForInternalValenceModel();
    1744   $FreeRadicalElectrons = $This->GetFreeRadicalElectrons();
    1745
    1746   $SumOfBondOrders = $This->GetSumOfBondOrders();
    1747   if (!defined $SumOfBondOrders) {
    1748     $SumOfBondOrders = 0;
    1749   }
    1750
    1751   @AvailableCommonValences = split /\,/, $CommonValences;
    1752
    1753   if (@AvailableCommonValences == 1) {
    1754     # Calculate potential total valence using the only available common valence...
    1755     $PotentialTotalValence = $AvailableCommonValences[0] + $AdjustedFormalCharge - $FreeRadicalElectrons;
    1756   }
    1757   else {
    1758     # Calculate potential total valence using common valence from a list of available valences
    1759     # that makes it higher than sum of bond orders or using the highest common valence...
    1760     VALENCE: for $AvailableCommonValence (@AvailableCommonValences) {
    1761       $PotentialTotalValence = $AvailableCommonValence + $AdjustedFormalCharge - $FreeRadicalElectrons;
    1762
    1763       if ($PotentialTotalValence < 0 || $PotentialTotalValence >= $SumOfBondOrders) {
    1764         last VALENCE;
    1765       }
    1766     }
    1767   }
    1768
    1769   return $PotentialTotalValence > 0 ? $PotentialTotalValence : 0;
    1770 }
    1771
    1772 # Adjust sign of the formal charge for potential total common valence calculation
    1773 # used during internal valence model to figure out number of implicit hydrogens.
    1774 #
    1775 sub _GetFormalChargeAdjustedForInternalValenceModel {
    1776   my($This) = @_;
    1777   my($FormalCharge, $GroupNumber, $SwitchSign);
    1778
    1779   $FormalCharge = $This->GetFormalCharge();
    1780   if ($FormalCharge == 0) {
    1781     return 0;
    1782   }
    1783
    1784   $GroupNumber = $This->GetGroupNumber();
    1785   if (!defined $GroupNumber) {
    1786     return $FormalCharge;
    1787   }
    1788
    1789   # Group numbers > 14 - Group numbers 15 (N), 16 (O), 17 (F), 18 (He)
    1790   #
    1791   # Formal charge sign is not adjusted. Positive and negative values result in the
    1792   # increase and decrease of valence.
    1793   #
    1794   # Group 14 containing C, Si, Ge, Sn, Pb...
    1795   #
    1796   # Formal charge sign is reversed for positive values. Both positive and negative
    1797   # values result in the decrease of valence.
    1798   #
    1799   # Group 13 containing B, Al, Ga, In, Tl...
    1800   #
    1801   # Formal charge sign is always reversed. Positive and negative values result in the
    1802   # decrease and increase of valence.
    1803   #
    1804   # Groups 1 (H) through 12 (Zn)...
    1805   #
    1806   # Formal charge sign is reversed for positive values. Both positive and negative
    1807   # values result in the decrease of valence.
    1808   #
    1809   # Lanthanides and actinides...
    1810   #
    1811   # Formal charge sign is reversed for positive values. Both positive and negative
    1812   # values result in the decrease of valence.
    1813   #
    1814
    1815   $SwitchSign = 0;
    1816   if (length $GroupNumber) {
    1817     GROUPNUMBER: {
    1818       if ($GroupNumber > 14) {
    1819         # Groups on the right side of C group in the periodic table...
    1820         $SwitchSign = 0;
    1821         last GROUPNUMBER;
    1822       }
    1823       if ($GroupNumber == 14) {
    1824         # Group containing C, Si, Ge, Sn, Pb...
    1825         $SwitchSign = ($FormalCharge > 0) ? 1 : 0;
    1826         last GROUPNUMBER;
    1827       }
    1828       if ($GroupNumber == 13) {
    1829         # Group containing B, Al, Ga, In, Tl...
    1830         $SwitchSign = 1;
    1831         last GROUPNUMBER;
    1832       }
    1833       # Groups 1 (H) through 12 (Zn)...
    1834       if ($GroupNumber >=1 && $GroupNumber <= 12) {
    1835         # Groups 1 (H) through 12 (Zn)...
    1836         $SwitchSign = ($FormalCharge > 0) ? 1 : 0;
    1837         last GROUPNUMBER;
    1838       }
    1839     }
    1840   }
    1841   else {
    1842     # Lanthanides and actinides...
    1843     $SwitchSign = ($FormalCharge > 0) ? 1 : 0;
    1844   }
    1845
    1846   if ($SwitchSign) {
    1847     $FormalCharge *= -1.0;
    1848   }
    1849
    1850   return $FormalCharge;
    1851 }
    1852
    1853 # Get lowest common valence...
    1854 sub GetLowestCommonValence {
    1855   my($This) = @_;
    1856
    1857   # Is LowestCommonValence property explicitly set?
    1858   if ($This->HasProperty('LowestCommonValence')) {
    1859     return $This->GetProperty('LowestCommonValence');
    1860   }
    1861   my($AtomicNumber, $LowestCommonValence);
    1862
    1863   $AtomicNumber = $This->{AtomicNumber};
    1864   if (!$AtomicNumber) {
    1865     return 0;
    1866   }
    1867   # Any need to differentiate between internal and other valence models...
    1868
    1869   # LowestCommonValence is not set for all elements...
    1870   $LowestCommonValence = PeriodicTable::GetElementLowestCommonValence($AtomicNumber);
    1871   if (!$LowestCommonValence) {
    1872     $LowestCommonValence = undef;
    1873   }
    1874
    1875   return $LowestCommonValence;
    1876 }
    1877
    1878 # Get highest common valence...
    1879 sub GetHighestCommonValence {
    1880   my($This) = @_;
    1881
    1882   # Is HighestCommonValence property explicitly set?
    1883   if ($This->HasProperty('HighestCommonValence')) {
    1884     return $This->GetProperty('HighestCommonValence');
    1885   }
    1886   my($AtomicNumber, $HighestCommonValence);
    1887
    1888   $AtomicNumber = $This->{AtomicNumber};
    1889   if (!$AtomicNumber) {
    1890     return 0;
    1891   }
    1892
    1893   # Any need to differentiate between internal and other valence models...
    1894
    1895   # HighestCommonValence is not set for all elements...
    1896   $HighestCommonValence = PeriodicTable::GetElementHighestCommonValence($AtomicNumber);
    1897   if (!$HighestCommonValence) {
    1898     $HighestCommonValence = undef;
    1899   }
    1900
    1901   return $HighestCommonValence;
    1902 }
    1903
    1904 # Get valence electrons...
    1905 sub GetValenceElectrons {
    1906   my($This) = @_;
    1907
    1908   # Is ValenceElectrons property explicitly set?
    1909   if ($This->HasProperty('ValenceElectrons')) {
    1910     return $This->GetProperty('ValenceElectrons');
    1911   }
    1912   my($AtomicNumber, $ValenceElectrons);
    1913
    1914   $AtomicNumber = $This->{AtomicNumber};
    1915   if (!$AtomicNumber) {
    1916     return 0;
    1917   }
    1918
    1919   $ValenceElectrons = PeriodicTable::GetElementValenceElectrons($AtomicNumber);
    1920
    1921   return $ValenceElectrons;
    1922 }
    1923
    1924 # Add hydrogens to specified atom in molecule and return number of hydrogens added:
    1925 #
    1926 #   o HydrogensToAdd = ImplicitHydrogenCount - ExplicitHydrogenCount
    1927 #
    1928 #   o XYZ are set to ZeroVector
    1929 #
    1930 sub AddHydrogens {
    1931   my($This, $HydrogenPositionsWarning) = @_;
    1932
    1933   # Is this atom in a molecule?
    1934   if (!$This->HasProperty('Molecule')) {
    1935     return undef;
    1936   }
    1937   if (!defined $HydrogenPositionsWarning) {
    1938     $HydrogenPositionsWarning = 1;
    1939   }
    1940   if ($HydrogenPositionsWarning) {
    1941     carp "Warning: ${ClassName}->AddHydrogens: The current release of MayaChemTools doesn't assign any hydrogen positions...";
    1942   }
    1943
    1944   # Is it an element symbol?
    1945   if (!$This->{AtomicNumber}) {
    1946     return 0;
    1947   }
    1948
    1949   my($Molecule, $HydrogensAdded, $HydrogensToAdd);
    1950
    1951   $Molecule = $This->GetProperty('Molecule');
    1952   $HydrogensAdded = 0;
    1953   $HydrogensToAdd = $This->GetNumOfMissingHydrogens();
    1954   if ($HydrogensToAdd <= 0) {
    1955     return $HydrogensAdded;
    1956   }
    1957
    1958   my($Count, $Hydrogen);
    1959
    1960   for $Count (1 .. $HydrogensToAdd) {
    1961     $HydrogensAdded++;
    1962
    1963     $Hydrogen = $Molecule->NewAtom('AtomSymbol' => 'H', 'XYZ' => [0, 0, 0]);
    1964     $Molecule->NewBond('Atoms' => [$This, $Hydrogen], 'BondOrder' => 1);
    1965   }
    1966
    1967   return $HydrogensAdded;
    1968 }
    1969
    1970 # Delete hydrogens attached to atom in molecule and return total number of hydrogens deleted...
    1971 sub DeleteHydrogens {
    1972   my($This) = @_;
    1973
    1974   # Is this atom in a molecule?
    1975   if (!$This->HasProperty('Molecule')) {
    1976     return undef;
    1977   }
    1978
    1979   # Is it an element symbol?
    1980   if (!$This->{AtomicNumber}) {
    1981     return 0;
    1982   }
    1983
    1984   my($Molecule, $Neighbor, $HydrogensDeleted, @Neighbors);
    1985
    1986   $Molecule = $This->GetProperty('Molecule');
    1987   $HydrogensDeleted = 0;
    1988   @Neighbors = $This->GetNeighbors();
    1989
    1990   NEIGHBOR: for $Neighbor (@Neighbors) {
    1991     if (!$Neighbor->IsHydrogen()) {
    1992       next NEIGHBOR;
    1993     }
    1994     $Molecule->_DeleteAtom($Neighbor);
    1995     $HydrogensDeleted++;
    1996   }
    1997
    1998   return $HydrogensDeleted;
    1999 }
    2000
    2001 # Copy atom and all its associated data...
    2002 sub Copy {
    2003   my($This) = @_;
    2004   my($Atom);
    2005
    2006   $Atom = Storable::dclone($This);
    2007
    2008   return $Atom;
    2009 }
    2010
    2011 # Get atomic invariant value...
    2012 #
    2013 sub GetAtomicInvariantValue {
    2014   my($This, $AtomicInvariant) = @_;
    2015   my($Value);
    2016
    2017   $Value = "";
    2018
    2019   ATOMICVARIANT: {
    2020     if ($AtomicInvariant =~ /^(AS|AtomSymbol|ElementSymbol)$/i) {
    2021       $Value = $This->GetAtomSymbol();
    2022       last ATOMICVARIANT;
    2023     }
    2024     if ($AtomicInvariant =~ /^(X|NumOfNonHydrogenAtomNeighbors|NumOfHeavyAtomNeighbors)$/i) {
    2025       $Value = $This->GetNumOfNonHydrogenAtomNeighbors();
    2026       last ATOMICVARIANT;
    2027     }
    2028     if ($AtomicInvariant =~ /^(BO|SumOfBondOrdersToNonHydrogenAtoms|SumOfBondOrdersToHeavyAtoms)$/i) {
    2029       $Value = $This->GetSumOfBondOrdersToNonHydrogenAtoms();
    2030       last ATOMICVARIANT;
    2031     }
    2032     if ($AtomicInvariant =~ /^(LBO|LargestBondOrderToNonHydrogenAtoms|LargestBondOrderToHeavyAtoms)$/i) {
    2033       $Value = $This->GetLargestBondOrderToNonHydrogenAtoms();
    2034       last ATOMICVARIANT;
    2035     }
    2036     if ($AtomicInvariant =~ /^(H|NumOfImplicitAndExplicitHydrogens)$/i) {
    2037       $Value = $This->GetNumOfHydrogens();
    2038       last ATOMICVARIANT;
    2039     }
    2040     if ($AtomicInvariant =~ /^(SB|NumOfSingleBondsToNonHydrogenAtoms|NumOfSingleBondsToHeavyAtoms)$/i) {
    2041       $Value = $This->GetNumOfSingleBondsToNonHydrogenAtoms();
    2042       last ATOMICVARIANT;
    2043     }
    2044     if ($AtomicInvariant =~ /^(DB|NumOfDoubleBondsToNonHydrogenAtoms|NumOfDoubleBondsToHeavyAtoms)$/i) {
    2045       $Value = $This->GetNumOfDoubleBondsToNonHydrogenAtoms();
    2046       last ATOMICVARIANT;
    2047     }
    2048     if ($AtomicInvariant =~ /^(TB|NumOfTripleBondsToNonHydrogenAtoms|NumOfTripleBondsToHeavyAtoms)$/i) {
    2049       $Value = $This->GetNumOfTripleBondsToNonHydrogenAtoms();
    2050       last ATOMICVARIANT;
    2051     }
    2052     if ($AtomicInvariant =~ /^(AB|NumOfAromaticBondsToNonHydrogenAtoms|NumOfAromaticBondsToHeavyAtoms)$/i) {
    2053       $Value = $This->GetNumOfAromaticBondsToNonHydrogenAtoms();
    2054       last ATOMICVARIANT;
    2055     }
    2056     if ($AtomicInvariant =~ /^(FC|FormalCharge)$/i) {
    2057       $Value = $This->GetFormalCharge();
    2058       $Value = defined $Value ? $Value : 0;
    2059       last ATOMICVARIANT;
    2060     }
    2061     if ($AtomicInvariant =~ /^(T|TotalNumOfAtomNeighbors)$/i) {
    2062       $Value = $This->GetNumOfNonHydrogenAtomNeighbors() + $This->GetNumOfHydrogens();
    2063       last ATOMICVARIANT;
    2064     }
    2065     if ($AtomicInvariant =~ /^(TSB|TotalNumOfSingleBonds)$/i) {
    2066       $Value = $This->GetNumOfSingleBondsToNonHydrogenAtoms() + $This->GetNumOfHydrogens();
    2067       last ATOMICVARIANT;
    2068     }
    2069     if ($AtomicInvariant =~ /^(Ar|Aromatic)$/i) {
    2070       $Value = $This->IsAromatic() ? 1 : 0;
    2071       last ATOMICVARIANT;
    2072     }
    2073     if ($AtomicInvariant =~ /^(RA|RingAtom)$/i) {
    2074       $Value = $This->IsInRing() ? 1 : 0;
    2075       last ATOMICVARIANT;
    2076     }
    2077     if ($AtomicInvariant =~ /^(Str|Stereochemistry)$/i) {
    2078       $Value = $This->GetStereochemistry();
    2079       $Value= (defined($Value) && ($Value =~ /^(R|S)$/i)) ? $Value : '';
    2080       last ATOMICVARIANT;
    2081     }
    2082     if ($AtomicInvariant =~ /^(AN|AtomicNumber)$/i) {
    2083       $Value = $This->GetAtomicNumber();
    2084       last ATOMICVARIANT;
    2085     }
    2086     if ($AtomicInvariant =~ /^(AM|AtomicMass)$/i) {
    2087       $Value = round($This->GetExactMass(), 4) + 0;
    2088       last ATOMICVARIANT;
    2089     }
    2090     if ($AtomicInvariant =~ /^(MN|MassNumber)$/i) {
    2091       $Value = $This->GetMassNumber();
    2092       last ATOMICVARIANT;
    2093     }
    2094     if ($AtomicInvariant =~ /^(SM|SpinMultiplicity)$/i) {
    2095       $Value = $This->GetSpinMultiplicity();
    2096       $Value = defined $Value ? $Value : '';
    2097       last ATOMICVARIANT;
    2098     }
    2099     $Value = "";
    2100     carp "Warning: ${ClassName}->GetAtomicInvariantValue: Unknown atomic invariant $AtomicInvariant...";
    2101   }
    2102
    2103   return $Value;
    2104 }
    2105
    2106 # Get period number of the atom..
    2107 #
    2108 sub GetPeriodNumber {
    2109   my($This) = @_;
    2110
    2111   # Is PeriodNumber property explicitly set?
    2112   if ($This->HasProperty('PeriodNumber')) {
    2113     return $This->GetProperty('PeriodNumber');
    2114   }
    2115   my($AtomicNumber, $PeriodNumber);
    2116
    2117   $AtomicNumber = $This->{AtomicNumber};
    2118   if (!$AtomicNumber) {
    2119     return 0;
    2120   }
    2121
    2122   $PeriodNumber = PeriodicTable::GetElementPeriodNumber($AtomicNumber);
    2123
    2124   return $PeriodNumber;
    2125 }
    2126
    2127 # Get group number of the atom..
    2128 #
    2129 sub GetGroupNumber {
    2130   my($This) = @_;
    2131
    2132   # Is GroupNumber property explicitly set?
    2133   if ($This->HasProperty('GroupNumber')) {
    2134     return $This->GetProperty('GroupNumber');
    2135   }
    2136   my($AtomicNumber, $GroupNumber);
    2137
    2138   $AtomicNumber = $This->{AtomicNumber};
    2139   if (!$AtomicNumber) {
    2140     return 0;
    2141   }
    2142
    2143   $GroupNumber = PeriodicTable::GetElementGroupNumber($AtomicNumber);
    2144
    2145   return $GroupNumber;
    2146 }
    2147
    2148 # Is it a specified topological pharmacophore atom type?
    2149 #
    2150 sub IsTopologicalPharmacophoreType {
    2151   my($This, $Type) = @_;
    2152
    2153   return $This->_IsFunctionalClassType($Type);
    2154 }
    2155
    2156 # Is it a specified functional class atom type?
    2157 #
    2158 sub IsFunctionalClassType {
    2159   my($This, $Type) = @_;
    2160
    2161   return $This->_IsFunctionalClassType($Type);
    2162 }
    2163
    2164 # Is it a specified functional/topological pharmacophore atom type?
    2165 #
    2166 sub _IsFunctionalClassType {
    2167   my($This, $Type) = @_;
    2168   my($Value);
    2169
    2170   $Value = 0;
    2171
    2172   TYPE: {
    2173     if ($Type =~ /^(HBD|HydrogenBondDonor)$/i) {
    2174       $Value = $This->IsHydrogenBondDonor();
    2175       last TYPE;
    2176     }
    2177     if ($Type =~ /^(HBA|HydrogenBondAcceptor)$/i) {
    2178       $Value = $This->IsHydrogenBondAcceptor();
    2179       last TYPE;
    2180     }
    2181     if ($Type =~ /^(PI|PositivelyIonizable)$/i) {
    2182       $Value = $This->IsPositivelyIonizable();
    2183       last TYPE;
    2184     }
    2185     if ($Type =~ /^(NI|NegativelyIonizable)$/i) {
    2186       $Value = $This->IsNegativelyIonizable();
    2187       last TYPE;
    2188     }
    2189     if ($Type =~ /^(H|Hydrophobic)$/i) {
    2190       $Value = $This->IsHydrophobic();
    2191       last TYPE;
    2192     }
    2193     if ($Type =~ /^(Ar|Aromatic)$/i) {
    2194       $Value = $This->IsAromatic();
    2195       last TYPE;
    2196     }
    2197     if ($Type =~ /^(Hal|Halogen)$/i) {
    2198       $Value = $This->IsHalogen();
    2199       last TYPE;
    2200     }
    2201     if ($Type =~ /^(RA|RingAtom)$/i) {
    2202       $Value = $This->IsInRing();
    2203       last TYPE;
    2204     }
    2205     if ($Type =~ /^(CA|ChainAtom)$/i) {
    2206       $Value = $This->IsNotInRing();
    2207       last TYPE;
    2208     }
    2209     $Value = 0;
    2210     carp "Warning: ${ClassName}->_IsType: Unknown functional/pharmacohore type $Type...";
    2211   }
    2212   return $Value;
    2213 }
    2214
    2215 # Is it a Hydrogen atom?
    2216 sub IsHydrogen {
    2217   my($This) = @_;
    2218
    2219   return ($This->{AtomicNumber} == 1) ? 1 : 0;
    2220 }
    2221
    2222 # Is it a Carbon atom?
    2223 sub IsCarbon {
    2224   my($This) = @_;
    2225
    2226   return ($This->{AtomicNumber} == 6) ? 1 : 0;
    2227 }
    2228
    2229 # Is it a Nitrogen atom?
    2230 sub IsNitrogen {
    2231   my($This) = @_;
    2232
    2233   return ($This->{AtomicNumber} == 7) ? 1 : 0;
    2234 }
    2235
    2236 # Is it a Oxygen atom?
    2237 sub IsOxygen {
    2238   my($This) = @_;
    2239
    2240   return ($This->{AtomicNumber} == 8) ? 1 : 0;
    2241 }
    2242
    2243 # Is it a Fluorine atom?
    2244 sub IsFluorine {
    2245   my($This) = @_;
    2246
    2247   return ($This->{AtomicNumber} == 9) ? 1 : 0;
    2248 }
    2249
    2250 # Is it a Silicon atom?
    2251 sub IsSilicon {
    2252   my($This) = @_;
    2253
    2254   return ($This->{AtomicNumber} == 14) ? 1 : 0;
    2255 }
    2256
    2257 # Is it a Phosphorus atom?
    2258 sub IsPhosphorus {
    2259   my($This) = @_;
    2260
    2261   return ($This->{AtomicNumber} == 15) ? 1 : 0;
    2262 }
    2263
    2264 # Is it a Sulphur atom?
    2265 sub IsSulphur {
    2266   my($This) = @_;
    2267
    2268   return $This->IsSulfur();
    2269 }
    2270
    2271 # Is it a Sulfur atom?
    2272 sub IsSulfur {
    2273   my($This) = @_;
    2274
    2275   return ($This->{AtomicNumber} == 16) ? 1 : 0;
    2276 }
    2277
    2278 # Is it a Chlorine atom?
    2279 sub IsChlorine {
    2280   my($This) = @_;
    2281
    2282   return ($This->{AtomicNumber} == 17) ? 1 : 0;
    2283 }
    2284
    2285 # Is it a Arsenic atom?
    2286 sub IsArsenic {
    2287   my($This) = @_;
    2288
    2289   return ($This->{AtomicNumber} == 33) ? 1 : 0;
    2290 }
    2291
    2292 # Is it a Selenium atom?
    2293 sub IsSelenium {
    2294   my($This) = @_;
    2295
    2296   return ($This->{AtomicNumber} == 34) ? 1 : 0;
    2297 }
    2298
    2299 # Is it a Bromine atom?
    2300 sub IsBromine {
    2301   my($This) = @_;
    2302
    2303   return ($This->{AtomicNumber} == 35) ? 1 : 0;
    2304 }
    2305
    2306 # Is it a Tellurium atom?
    2307 sub IsTellurium {
    2308   my($This) = @_;
    2309
    2310   return ($This->{AtomicNumber} == 52) ? 1 : 0;
    2311 }
    2312
    2313 # Is it a Iodine atom?
    2314 sub IsIodine {
    2315   my($This) = @_;
    2316
    2317   return ($This->{AtomicNumber} == 53) ? 1 : 0;
    2318 }
    2319
    2320 # Is it a hetro atom? (N, O, F, P, S, Cl, Br, I)
    2321 sub IsHeteroAtom {
    2322   my($This) = @_;
    2323
    2324   return ($This->{AtomicNumber} =~ /^(7|8|9|15|16|17|35|53)$/) ? 1 : 0;
    2325 }
    2326
    2327 # Is it a halogen atom? (F, Cl, Br, I)
    2328 sub IsHalogen {
    2329   my($This) = @_;
    2330
    2331   return ($This->{AtomicNumber} =~ /^(9|17|35|53)$/) ? 1 : 0;
    2332 }
    2333
    2334 # Is it classified as metallic?
    2335 sub IsMetallic {
    2336   my($This) = @_;
    2337   my($Classification);
    2338
    2339   $Classification = PeriodicTable::GetElementClassification($This->{AtomicNumber});
    2340
    2341   return ($Classification =~ /^Metallic$/i) ? 1 : 0;
    2342 }
    2343
    2344 # Is it a non carbon or hydrogen atom? (C, H)
    2345 sub IsNonCarbonOrHydrogen {
    2346   my($This) = @_;
    2347
    2348   return ($This->{AtomicNumber} =~ /^(1|6)$/) ? 0 : 1;
    2349 }
    2350
    2351 # Is it a polar atom? ( N, O,  P, S)
    2352 sub IsPolarAtom {
    2353   my($This) = @_;
    2354
    2355   return ($This->{AtomicNumber} =~ /^(7|8|15|16)$/) ? 1 : 0;
    2356 }
    2357
    2358 # Is it an isotope?
    2359 sub IsIsotope {
    2360   my($This) = @_;
    2361
    2362   my($AtomicNumber) = $This->{AtomicNumber};
    2363   if (!$AtomicNumber) {
    2364     return 0;
    2365   }
    2366
    2367   if (!$This->HasProperty('MassNumber')) {
    2368     return 0;
    2369   }
    2370   my($MassNumber, $MostAbundantMassNumber);
    2371
    2372   $MassNumber = $This->GetProperty('MassNumber');
    2373   $MostAbundantMassNumber = PeriodicTable::GetElementMostAbundantNaturalIsotopeMassNumber($AtomicNumber);
    2374
    2375   return ($MassNumber == $MostAbundantMassNumber) ? 0 : 1;
    2376 }
    2377
    2378 # Is it a terminal atom?
    2379 sub IsTerminal {
    2380   my($This) = @_;
    2381
    2382   # Is this atom in a molecule?
    2383   if (!$This->HasProperty('Molecule')) {
    2384     return undef;
    2385   }
    2386
    2387   return ($This->GetNumOfNonHydrogenAtomNeighbors() <= 1) ? 1 : 0
    2388
    2389 }
    2390
    2391 # Is aromatic property set for the atom?
    2392 sub IsAromatic {
    2393   my($This) = @_;
    2394   my($Aromatic);
    2395
    2396   $Aromatic = $This->GetAromatic();
    2397
    2398   return (defined($Aromatic) && $Aromatic) ? 1 : 0;
    2399 }
    2400
    2401 # Is this a hydrogen atom and attached to one of these atoms: N, O, P, S
    2402 sub IsPolarHydrogen {
    2403   my($This) = @_;
    2404
    2405   if (!$This->IsHydrogen()) {
    2406     return 0;
    2407   }
    2408
    2409   my(@Bonds);
    2410   @Bonds = $This->GetBonds();
    2411   if (@Bonds > 1) {
    2412     return 0;
    2413   }
    2414
    2415   my($Bond, $BondedAtom);
    2416   ($Bond) = @Bonds;
    2417   $BondedAtom = $Bond->GetBondedAtom($This);
    2418
    2419   return $BondedAtom->IsPolarAtom() ? 1 : 0;
    2420 }
    2421
    2422 # Is it a hydrogen bond donor atom?
    2423 #
    2424 sub IsHBondDonor {
    2425   my($This, $HydrogenBondsType) = @_;
    2426
    2427   return $This->IsHydrogenBondDonor($HydrogenBondsType);
    2428 }
    2429
    2430 # The currrent release of MayaChemTools supports identification of two types of
    2431 # hydrogen bond donor and acceptor atoms with these names:
    2432 #
    2433 # HBondsType1 or HydrogenBondsType1
    2434 # HBondsType2 or HydrogenBondsType2
    2435 #
    2436 # The names of these hydrogen bond types are rather arbitrary. However, their
    2437 # definitions have specific meaning and are as follows:
    2438 #
    2439 # HydrogenBondsType1 [ Ref 60-61, Ref 65-66 ]:
    2440 #   . Donor: NH, NH2, NH3, OH - Any N and O with available H
    2441 #   . Acceptor: N[!H], O - Any N without available H and any O
    2442 #
    2443 # HydrogenBondsType2 [ Ref 91 ]:
    2444 #   . Donor: NH, NH2, NH3, OH - Any N and O with availabe H
    2445 #   . Acceptor: N, O - Any N and O
    2446 #
    2447 # Note:
    2448 #   . HydrogenBondsType2 definition corresponds to Rule of 5.
    2449 #
    2450
    2451 # Is it a hydrogen bond donor atom?
    2452 #
    2453 # The currrent release of MayaChemTools supports identification of two types of
    2454 sub IsHydrogenBondDonor {
    2455   my($This, $HydrogenBondsType) = @_;
    2456   my($Status);
    2457
    2458   $HydrogenBondsType = defined $HydrogenBondsType ? $HydrogenBondsType : 'HBondsType1';
    2459   $Status = 0;
    2460
    2461   HYDROGENBONDSTYPE: {
    2462
    2463       if ($HydrogenBondsType =~ /^(HBondsType1|HydrogenBondsType1)$/i) {
    2464         $Status = $This->_IsHydrogenBondDonorOfType1();
    2465         last HYDROGENBONDSTYPE;
    2466       }
    2467
    2468       if ($HydrogenBondsType =~ /^(HBondsType2|HydrogenBondsType2)$/i) {
    2469         $Status = $This->_IsHydrogenBondDonorOfType2();
    2470         last HYDROGENBONDSTYPE;
    2471       }
    2472
    2473       $Status = 0;
    2474       carp "Warning: ${ClassName}->IsHydrogenBondDonor: The current release of MayaChemTools doesn't support specified value, $HydrogenBondsType, for HydrogenBondsType. Valid values: HBondsType1, HydrogenBondsType1, HBondsType2 HydrogenBondsType2 ...";
    2475   }
    2476
    2477   return $Status;
    2478 }
    2479
    2480 # Is it a MayaChemTools HBondType1 hydrogen bond donor atom?
    2481 #
    2482 sub _IsHydrogenBondDonorOfType1 {
    2483   my($This) = @_;
    2484
    2485   return $This->_IsHydrogenBondDonorOfType1OrType2();
    2486 }
    2487
    2488 # Is it a MayaChemTools HBondType2 hydrogen bond donor atom?
    2489 #
    2490 sub _IsHydrogenBondDonorOfType2 {
    2491   my($This) = @_;
    2492
    2493   return $This->_IsHydrogenBondDonorOfType1OrType2();
    2494 }
    2495
    2496 # Is it a hydrogen bond donor atom of MayaChemTools Type1 or Type2?
    2497 #
    2498 # HydrogenBondDonor definition [ Ref 60-61, Ref 65-66, Ref 91 ]: NH, NH2, OH
    2499 #
    2500 # In other words:
    2501 #   . NH, NH2 - Nitrogen atom with available hydrogen
    2502 #   . OH - Oxygen atom with avilable hydrogen
    2503 #
    2504 sub _IsHydrogenBondDonorOfType1OrType2 {
    2505   my($This) = @_;
    2506
    2507   # Is this atom in a molecule?
    2508   if (!$This->HasProperty('Molecule')) {
    2509     return 0;
    2510   }
    2511
    2512   # Is it N or O?
    2513   if ($This->{AtomicNumber} !~ /^(7|8)$/) {
    2514     return 0;
    2515   }
    2516
    2517   # Any explicitly attached hydrogens?
    2518   if ($This->GetExplicitHydrogens()) {
    2519     return 1;
    2520   }
    2521
    2522   # Any missing hydrogens?
    2523   return $This->GetNumOfMissingHydrogens() ? 1 : 0;
    2524 }
    2525
    2526 # Is it a hydrogen bond acceptor atom?
    2527 #
    2528 sub IsHBondAcceptor {
    2529   my($This, $HydrogenBondsType) = @_;
    2530
    2531   return $This->IsHydrogenBondAcceptor($HydrogenBondsType);
    2532 }
    2533
    2534 # Is it a hydrogen bond acceptor atom?
    2535 #
    2536 sub IsHydrogenBondAcceptor {
    2537   my($This, $HydrogenBondsType) = @_;
    2538   my($Status);
    2539
    2540   $HydrogenBondsType = defined $HydrogenBondsType ? $HydrogenBondsType : 'HBondsType1';
    2541   $Status = 0;
    2542
    2543   HYDROGENBONDSTYPE: {
    2544
    2545       if ($HydrogenBondsType =~ /^(HBondsType1|HydrogenBondsType1)$/i) {
    2546         $Status = $This->_IsHydrogenBondAcceptorOfType1();
    2547         last HYDROGENBONDSTYPE;
    2548       }
    2549
    2550       if ($HydrogenBondsType =~ /^(HBondsType2|HydrogenBondsType2)$/i) {
    2551         $Status = $This->_IsHydrogenBondAcceptorOfType2();
    2552         last HYDROGENBONDSTYPE;
    2553       }
    2554
    2555       $Status = 0;
    2556       carp "Warning: ${ClassName}->IsHydrogenBondAcceptor: The current release of MayaChemTools doesn't support specified value, $HydrogenBondsType, for HydrogenBondsType. Valid values: HBondsType1, HydrogenBondsType1, HBondsType2 HydrogenBondsType2 ...";
    2557   }
    2558
    2559   return $Status;
    2560 }
    2561
    2562 # Is it a MayaChemTools HBondType1 hydrogen bond acceptor atom?
    2563 #
    2564 # HydrogenBondAcceptor definition [ Ref 60-61, Ref 65-66 ]: N[!H], O
    2565 #
    2566 # In other words:
    2567 #   . N[!H] - Nitrogen atom with no hydrogen
    2568 #   . O - Oxygen atom
    2569 #
    2570 sub _IsHydrogenBondAcceptorOfType1 {
    2571   my($This) = @_;
    2572
    2573   # Is this atom in a molecule?
    2574   if (!$This->HasProperty('Molecule')) {
    2575     return 0;
    2576   }
    2577
    2578   # Is it N or O?
    2579   if ($This->{AtomicNumber} !~ /^(7|8)$/) {
    2580     return 0;
    2581   }
    2582
    2583   # Is it O?
    2584   if ($This->{AtomicNumber} == 8 ) {
    2585     return 1;
    2586   }
    2587
    2588   # Any explicitly attached hydrogens?
    2589   if ($This->GetExplicitHydrogens()) {
    2590     return 0;
    2591   }
    2592
    2593   # Any missing hydrogens?
    2594   return $This->GetNumOfMissingHydrogens() ? 0 : 1;
    2595 }
    2596
    2597 # Is it a MayaChemTools HBondType2 hydrogen bond acceptor atom?
    2598 #
    2599 # HydrogenBondAcceptor definition [ Ref 91 ]: N, O
    2600 #
    2601 # In other words:
    2602 #   . Any Nitrogen or Oxygen atom
    2603 #
    2604 # Note:
    2605 #   . HydrogenBondsType2 definition corresponds to Rule of 5.
    2606 #
    2607 sub _IsHydrogenBondAcceptorOfType2 {
    2608   my($This) = @_;
    2609
    2610   # Is this atom in a molecule?
    2611   if (!$This->HasProperty('Molecule')) {
    2612     return 0;
    2613   }
    2614
    2615   return ($This->{AtomicNumber} =~ /^(7|8)$/) ? 1 : 0;
    2616 }
    2617
    2618 # Is it a positively ionizable atom?
    2619 #
    2620 # PositivelyIonizable defintion [ Ref 60-61, Ref 65-66 ]: +, NH2
    2621 #
    2622 # In other words:
    2623 #   . Any atom with positve formal charge
    2624 #   . NH2 - Nitogen atom in amino group
    2625 #
    2626 sub IsPositivelyIonizable {
    2627   my($This) = @_;
    2628   my($FormalCharge);
    2629
    2630   # Is this atom in a molecule?
    2631   if (!$This->HasProperty('Molecule')) {
    2632     return 0;
    2633   }
    2634
    2635   # Any explicit positive formal charge?
    2636   $FormalCharge = $This->GetFormalCharge();
    2637   if (defined($FormalCharge) && $FormalCharge > 0) {
    2638     return 1;
    2639   }
    2640
    2641   # Is it  N?
    2642   if ($This->{AtomicNumber} != 7 ) {
    2643     return 0;
    2644   }
    2645
    2646   return ($This->GetNumOfHydrogens() == 2) ? 1 : 0;
    2647 }
    2648
    2649 # Is it a negatively ionizable atom?
    2650 #
    2651 # NegativelyIonizable definition [ Ref 60-61, Ref 65-66 ]: -, C(=O)OH, S(=O)OH, P(=O)OH
    2652 #
    2653 # In other words:
    2654 #   . Any atom with negative formal charge
    2655 #   . Carbon atom in C(=O)OH group
    2656 #   . Phosphorous in P(=O)OH group
    2657 #   . Sulfur atom in S(=O)OH group
    2658 #
    2659 sub IsNegativelyIonizable {
    2660   my($This) = @_;
    2661   my($FormalCharge);
    2662
    2663   # Is this atom in a molecule?
    2664   if (!$This->HasProperty('Molecule')) {
    2665     return 0;
    2666   }
    2667
    2668   # Any explicit negative formal charge?
    2669   $FormalCharge = $This->GetFormalCharge();
    2670   if (defined($FormalCharge) && $FormalCharge < 0) {
    2671     return 1;
    2672   }
    2673
    2674   # Is it C, P or S?
    2675   if ($This->{AtomicNumber} !~ /^(6|15|16)$/ ) {
    2676     return 0;
    2677   }
    2678
    2679   # Collect oxygens connected to C, P or S with single or double bonds and not connected to
    2680   # any other heavy atom...
    2681   my($Neighbor, $NeighborOxygenBondOrder, $NumOfNeighborOxygensWithSingleBonds, $NumOfNeighborOxygensWithDoubleBonds);
    2682
    2683   $NumOfNeighborOxygensWithSingleBonds = 0; $NumOfNeighborOxygensWithDoubleBonds = 0;
    2684
    2685   NEIGHBOR: for $Neighbor ($This->GetNeighbors()) {
    2686     # Is it an oxygen?
    2687     if ($Neighbor->{AtomicNumber} != 8) {
    2688       next NEIGHBOR;
    2689     }
    2690     # Is oxygent connected to only heavy atom?
    2691     if ($Neighbor->GetNumOfHeavyAtomNeighbors() != 1) {
    2692       next NEIGHBOR;
    2693     }
    2694     $NeighborOxygenBondOrder = $This->GetBondToAtom($Neighbor)->GetBondOrder();
    2695
    2696     if ($NeighborOxygenBondOrder == 2) {
    2697       $NumOfNeighborOxygensWithDoubleBonds++;
    2698     }
    2699     elsif ($NeighborOxygenBondOrder == 1) {
    2700       $NumOfNeighborOxygensWithSingleBonds++;
    2701     }
    2702   }
    2703   return ($NumOfNeighborOxygensWithDoubleBonds >= 1 && $NumOfNeighborOxygensWithSingleBonds >= 1) ? 1 : 0;
    2704 }
    2705
    2706 # Is it a liphophilic atom?
    2707 #
    2708 # Lipophilic definition [ Ref 60-61, Ref 65-66 ]: C(C)(C)(C)(C), Cl, Br, I, S(C)(C)
    2709 #
    2710 # In other words:
    2711 #   . C(C)(C)(C)(C) - Carbon atom connected to only other carbons
    2712 #   . Chlorine, Bromine or Iodine atom
    2713 #   . S(C)(C) - Sulfur connected to two carbons
    2714 #
    2715 sub IsLipophilic {
    2716   my($This) = @_;
    2717
    2718   # Is this atom in a molecule?
    2719   if (!$This->HasProperty('Molecule')) {
    2720     return 0;
    2721   }
    2722
    2723   # Is it Cl, Br, I?
    2724   if ($This->{AtomicNumber} =~ /^(17|35|53)$/) {
    2725     return 1;
    2726   }
    2727
    2728   # Is it C, S?
    2729   if ($This->{AtomicNumber} !~ /^(6|16)$/) {
    2730     return 0;
    2731   }
    2732
    2733   # Are all heavy atom neighbors Carbons?
    2734   my($HeavyAtomNeighbor, @HeavyAtomNeighbors);
    2735   @HeavyAtomNeighbors = ();
    2736   @HeavyAtomNeighbors = $This->GetHeavyAtomNeighbors();
    2737
    2738   for $HeavyAtomNeighbor (@HeavyAtomNeighbors) {
    2739     if ($HeavyAtomNeighbor->{AtomicNumber} != 6) {
    2740       return 0;
    2741     }
    2742   }
    2743
    2744   # Does sulfur has two carbon neighbors?
    2745   if ($This->{AtomicNumber} == 16) {
    2746     if (@HeavyAtomNeighbors != 2) {
    2747       return 0;
    2748     }
    2749   }
    2750   return 1;
    2751 }
    2752
    2753 # Is it hydrophobic?
    2754 #
    2755 sub IsHydrophobic {
    2756   my($This) = @_;
    2757
    2758   return $This->IsLipophilic();
    2759 }
    2760
    2761 # Is it a Nitrogen atom in Guadinium group?
    2762 #
    2763 sub IsGuadiniumNitrogen {
    2764   my($This) = @_;
    2765
    2766   # Is it Nitrogen?
    2767   if (!$This->IsNitrogen()) {
    2768     return 0;
    2769   }
    2770
    2771   # Is it connected to a Guadinium Carbon?
    2772   my($AtomNeighbor);
    2773
    2774   for $AtomNeighbor ($This->GetNonHydrogenAtomNeighbors()) {
    2775     if ($AtomNeighbor->IsGuadiniumCarbon()) {
    2776       return 1;
    2777     }
    2778   }
    2779
    2780   return 0;
    2781 }
    2782
    2783 # Is it a Carbon atom in Guadinium group?
    2784 #
    2785 # Guadinium group definition:
    2786 #
    2787 #   R2N-C(=NR)-(NR2) or R2N-C(=NR2+)-(NR2)
    2788 #
    2789 #   where:
    2790 #      . R = Hydrogens or group of atoms attached through Carbon
    2791 #      . Only one of the three Nitrogens has a double bond to Carbon and has optional
    2792 #        formal charge allowing it to be neutral or charged state
    2793 #
    2794 sub IsGuadiniumCarbon {
    2795   my($This) = @_;
    2796
    2797   # Is it Carbon?
    2798   if (!$This->IsCarbon()) {
    2799     return 0;
    2800   }
    2801
    2802   # Match atom neighborhood...
    2803   my($CentralAtomSpec, @NbrAtomSpecsRef, @NbrBondSpecsRef, @NbrOfNbrAtomSpecsRef);
    2804
    2805   $CentralAtomSpec = 'C.X3.BO4';
    2806   @NbrAtomSpecsRef = ('N.FC0', 'N.FC0', 'N.FC0,N.FC+1');
    2807   @NbrBondSpecsRef = ('-', '-', '=');
    2808   @NbrOfNbrAtomSpecsRef = ('C,H', 'C,H', 'C,H');
    2809
    2810   if ($This->DoesAtomNeighborhoodMatch($CentralAtomSpec, \@NbrAtomSpecsRef, \@NbrBondSpecsRef, \@NbrOfNbrAtomSpecsRef)) {
    2811     return 1;
    2812   }
    2813
    2814   return 0;
    2815 }
    2816
    2817 # Is it a Nitrogen atom in Amide group?
    2818 #
    2819 sub IsAmideNitrogen {
    2820   my($This) = @_;
    2821
    2822   # Is it Nitrogen?
    2823   if (!$This->IsNitrogen()) {
    2824     return 0;
    2825   }
    2826
    2827   # Is it connected to a Amide Carbon?
    2828   my($AtomNeighbor);
    2829
    2830   for $AtomNeighbor ($This->GetNonHydrogenAtomNeighbors()) {
    2831     if ($AtomNeighbor->IsAmideCarbon()) {
    2832       return 1;
    2833     }
    2834   }
    2835
    2836   return 0;
    2837 }
    2838
    2839 # Is it a Carbon atom in Amide group?
    2840 #
    2841 # Amide group definition: R-C(=O)-N(-R')-R''
    2842 #
    2843 #   where:
    2844 #      . R = Hydrogen or groups of atoms attached through Carbon
    2845 #      . R' = Hydrogens or groups of atoms attached through Carbon or hetro atoms
    2846 #      . R'' = Hydrogens or groups of atoms attached through Carbon or hetro atoms
    2847 #
    2848 sub IsAmideCarbon {
    2849   my($This) = @_;
    2850
    2851   # Is this atom in a molecule?
    2852   if (!$This->HasProperty('Molecule')) {
    2853     return 0;
    2854   }
    2855
    2856   # Is it Carbon?
    2857   if (!$This->IsCarbon()) {
    2858     return 0;
    2859   }
    2860
    2861   # Match atom neighborhood...
    2862   my($CentralAtomSpec, @NbrAtomSpecsRef, @NbrBondSpecsRef, @NbrOfNbrAtomSpecsRef);
    2863
    2864   $CentralAtomSpec = 'C.X3.BO4,C.X2.BO3';
    2865   @NbrAtomSpecsRef = ('C,H', 'O', 'N');
    2866   @NbrBondSpecsRef = ('-', '=', '-');
    2867   @NbrOfNbrAtomSpecsRef = ('C,H', 'C', 'C,H,N,O,S,P,F,Cl,Br,I');
    2868
    2869   if ($This->DoesAtomNeighborhoodMatch($CentralAtomSpec, \@NbrAtomSpecsRef, \@NbrBondSpecsRef, \@NbrOfNbrAtomSpecsRef)) {
    2870     return 1;
    2871   }
    2872
    2873   return 0;
    2874 }
    2875
    2876 # Is it a Oxygen atom in Carboxylate group?
    2877 #
    2878 sub IsCarboxylateOxygen {
    2879   my($This) = @_;
    2880
    2881   return $This->_MatchCarboxylateAndOrCarboxylOxygen('Carboxylate');
    2882 }
    2883
    2884 # Is it a Carbon atom in Carboxylate group?
    2885 #
    2886 # Carboxyl group definition: R-C(=O)-O-
    2887 #
    2888 sub IsCarboxylateCarbon {
    2889   my($This) = @_;
    2890
    2891   return $This->_MatchCarboxylateAndOrCarboxylCarbon('Carboxylate');
    2892 }
    2893
    2894 # Is it a Oxygen atom in Carboxyl group?
    2895 #
    2896 sub IsCarboxylOxygen {
    2897   my($This) = @_;
    2898
    2899   return $This->_MatchCarboxylateAndOrCarboxylOxygen('Carboxyl');
    2900 }
    2901
    2902 # Is it a Carbon atom in Carboxyl group?
    2903 #
    2904 # Carboxyl group definition: R-C(=O)-OH
    2905 #
    2906 sub IsCarboxylCarbon {
    2907   my($This) = @_;
    2908
    2909   return $This->_MatchCarboxylateAndOrCarboxylCarbon('Carboxyl');
    2910 }
    2911
    2912 # Match Carboxylate and/or Carboxyl oxygen...
    2913 #
    2914 sub _MatchCarboxylateAndOrCarboxylOxygen {
    2915   my($This, $Mode) = @_;
    2916
    2917   # Is it Oxygen?
    2918   if (!$This->IsOxygen()) {
    2919     return 0;
    2920   }
    2921
    2922   # Is it connected to a Carboxylate Carbon?
    2923   my($AtomNeighbor);
    2924
    2925   for $AtomNeighbor ($This->GetNonHydrogenAtomNeighbors()) {
    2926     if ($AtomNeighbor->_MatchCarboxylateAndOrCarboxylCarbon($Mode)) {
    2927       return 1;
    2928     }
    2929   }
    2930
    2931   return 0;
    2932 }
    2933
    2934 # Match Carboxylate and Carboxyl Carbon
    2935 #
    2936 # Carboxylate group definition: R-C(=O)-O-
    2937 # Carboxyl group definition: R-C(=O)-OH
    2938 #
    2939 #   where:
    2940 #      . R = Hydrogens or groups of atoms attached through Carbon
    2941 #
    2942 sub _MatchCarboxylateAndOrCarboxylCarbon {
    2943   my($This, $Mode) = @_;
    2944
    2945   # Is this atom in a molecule?
    2946   if (!$This->HasProperty('Molecule')) {
    2947     return 0;
    2948   }
    2949
    2950   # Is it Carbon?
    2951   if (!$This->IsCarbon()) {
    2952     return 0;
    2953   }
    2954
    2955   # Match atom neighborhood...
    2956   my($CentralAtomSpec, @NbrAtomSpecsRef, @NbrBondSpecsRef, @NbrOfNbrAtomSpecsRef);
    2957
    2958   $CentralAtomSpec = 'C.X3.BO4,C.X2.BO3';
    2959   MODE: {
    2960     if ($Mode =~ /^Carboxylate$/i) {
    2961       @NbrAtomSpecsRef = ('C,H', 'O', 'O.X1.FC-1');
    2962       last MODE;
    2963     }
    2964     if ($Mode =~ /^Carboxyl$/i) {
    2965       @NbrAtomSpecsRef = ('C,H', 'O', 'O.X1.FC0');
    2966       last MODE;
    2967     }
    2968     if ($Mode =~ /^CarboxylateOrCarboxyl$/i) {
    2969       @NbrAtomSpecsRef = ('C,H', 'O', 'O.X1.FC-1,O.X1.FC0');
    2970       last MODE;
    2971     }
    2972     carp "Warning: ${ClassName}->_MatchCarboxylateAndCarboxylCarbon.: Unknown mode $Mode...";
    2973     return 0;
    2974   }
    2975   @NbrBondSpecsRef = ('-', '=', '-');
    2976   @NbrOfNbrAtomSpecsRef = ('C,H', 'C', 'C');
    2977
    2978   if ($This->DoesAtomNeighborhoodMatch($CentralAtomSpec, \@NbrAtomSpecsRef, \@NbrBondSpecsRef, \@NbrOfNbrAtomSpecsRef)) {
    2979     return 1;
    2980   }
    2981
    2982   return 0;
    2983 }
    2984
    2985 # Is it a Oxygen atom in Phosphate group?
    2986 #
    2987 sub IsPhosphateOxygen {
    2988   my($This) = @_;
    2989
    2990   # Is it Oxygen?
    2991   if (!$This->IsOxygen()) {
    2992     return 0;
    2993   }
    2994
    2995   # Is it connected to a Phosphate Phosphorus?
    2996   my($AtomNeighbor);
    2997
    2998   for $AtomNeighbor ($This->GetNonHydrogenAtomNeighbors()) {
    2999     if ($AtomNeighbor->IsPhosphatePhosphorus()) {
    3000       return 1;
    3001     }
    3002   }
    3003
    3004   return 0;
    3005 }
    3006
    3007 # Is it a Phosphorus atom in Phosphate group?
    3008 #
    3009 # Phosphate group definition: AO-(O=)P(-OA)-OA
    3010 #
    3011 #   where:
    3012 #      . A = Any Groups of atoms including hydrogens
    3013 #
    3014 sub IsPhosphatePhosphorus {
    3015   my($This) = @_;
    3016
    3017   # Is this atom in a molecule?
    3018   if (!$This->HasProperty('Molecule')) {
    3019     return 0;
    3020   }
    3021
    3022   # Is it Phosphorus?
    3023   if (!$This->IsPhosphorus()) {
    3024     return 0;
    3025   }
    3026
    3027   # Match atom neighborhood...
    3028   my($CentralAtomSpec, @NbrAtomSpecsRef, @NbrBondSpecsRef, @NbrOfNbrAtomSpecsRef);
    3029
    3030   $CentralAtomSpec = 'P.X4.BO5';
    3031   @NbrAtomSpecsRef = ('O', 'O', 'O', 'O');
    3032   @NbrBondSpecsRef = ('-', '=', '-', '-');
    3033   @NbrOfNbrAtomSpecsRef = (undef, undef, undef, undef);
    3034
    3035   if ($This->DoesAtomNeighborhoodMatch($CentralAtomSpec, \@NbrAtomSpecsRef, \@NbrBondSpecsRef, \@NbrOfNbrAtomSpecsRef)) {
    3036     return 1;
    3037   }
    3038
    3039   return 0;
    3040 }
    3041
    3042
    3043 # Match central atom and its neighborhood using specified atom and bonds specifications...
    3044 #
    3045 # Let:
    3046 #   AS = Atom symbol corresponding to element symbol, atomic number (#n) or any
    3047 #        atom (A)
    3048 #
    3049 #   X<n>   = Number of non-hydrogen atom neighbors or heavy atoms attached to atom
    3050 #   T<n>   = Total number of atom neighbors including implcit and explicit hydrogens
    3051 #   BO<n> = Sum of bond orders to non-hydrogen atom neighbors or heavy atoms attached to atom
    3052 #   LBO<n> = Largest bond order of non-hydrogen atom neighbors or heavy atoms attached to atom
    3053 #   SB<n> = Number of single bonds to non-hydrogen atom neighbors or heavy atoms attached to atom
    3054 #   TSB<n> = Total number of single bonds to atom neighbors including implcit and explicit hydrogens
    3055 #   DB<n> = Number of double bonds to non-hydrogen atom neighbors or heavy atoms attached to atom
    3056 #   TB<n> = Number of triple bonds to non-hydrogen atom neighbors or heavy atoms attached to atom
    3057 #   H<n>   = Number of implicit and explicit hydrogens for atom
    3058 #   Ar     = Aromatic annotation indicating whether atom is aromatic
    3059 #   RA or RA<n>  = Ring atom annotation indicating whether atom is a ring
    3060 #   TR<n>  = Total number of rings containing atom
    3061 #   FC<+n/-n> = Formal charge assigned to atom
    3062 #   MN<n> = Mass number indicating isotope other than most abundant isotope
    3063 #   SM<n> = Spin multiplicity of atom. Possible values: 1 (singlet), 2 (doublet) or 3 (triplet)
    3064 #
    3065 # Then:
    3066 #
    3067 #   Atom specification corresponds to:
    3068 #
    3069 #     AS.X<n>.T<n>.BO<n>.LBO<n>.<SB><n>.TSB<n>.<DB><n>.<TB><n>.H<n>.Ar.RA<n>.TR<n>FC<+n/-n>.MN<n>.SM<n>
    3070 #
    3071 # Except for AS which is a required atomic invariant in atom specification, all other atomic invariants are
    3072 # optional. For an atom specification to match an atom, the values of all specified atomic invariants must
    3073 # match. Exclamation in from of atomic invariant can be used to negate its effect during the match.
    3074 #
    3075 # A comma delimited atom specification string is used to match any one of the specifed atom specification.
    3076 #
    3077 # Notes:
    3078 #   . During atom specification match to an atom, the first atomic invariant is always assumed to
    3079 #     atom symbol.
    3080 #   . Atom match specfication is based on AtomicInvariantAtomTypes implemented in
    3081 #     AotmTypes::AtomicInvariantAtomType.pm module
    3082 #
    3083 # Examples:
    3084 #     . ('N', 'N', 'N')
    3085 #     . ('N.FC0', 'N.FC0', 'N,N.FC+1.H1')
    3086 #     . ('N.H2', 'N.H2', 'N.H1')
    3087 #     . ('C,N', '!N', '!H')
    3088 #     . ('C,N', 'N.Ar', 'N.R5')
    3089 #
    3090 # Let:
    3091 #   -|1|s|Single = Single bond
    3092 #   =|2|d|Double = Double bond
    3093 #   #|3|t|Triple  = Triple bond
    3094 #   :|1.5|a|Ar|Aromatic = Aromatic bond
    3095 #
    3096 #   @|RB|Ring = Ring bond
    3097 #   ~|*|Any = Any bond
    3098 #
    3099 # Then:
    3100 #
    3101 #   Bond specification corresponds to:
    3102 #
    3103 #     -.:
    3104 #     =.@
    3105 #     Double.Aromatic
    3106 #
    3107 # For a bond specification to match bond between two atoms, the values of all specified bond symbols must
    3108 # match. Exclamation in from of bond symbol can be used to negate its effect during the match.
    3109 #
    3110 # A comma delimited bond specification string is used to match any one of the specifed atom specification.
    3111 #
    3112 sub DoesAtomNeighborhoodMatch {
    3113   my($CentralAtom, $CentralAtomSpec, $NbrAtomSpecsRef, $NbrBondSpecsRef, $NbrOfNbrAtomSpecsRef) = @_;
    3114   my($NumOfNbrAtomSpecs, $NumOfNbrBondSpecs, $NumOfNbrOfNbrAtomSpecs);
    3115
    3116   # Is this atom in a molecule?
    3117   if (!$CentralAtom->HasProperty('Molecule')) {
    3118     return 0;
    3119   }
    3120
    3121   $NumOfNbrAtomSpecs = defined $NbrAtomSpecsRef ? scalar @{$NbrAtomSpecsRef} : 0;
    3122   $NumOfNbrBondSpecs = defined $NbrBondSpecsRef ? scalar @{$NbrBondSpecsRef} : 0;
    3123   $NumOfNbrOfNbrAtomSpecs = defined $NbrOfNbrAtomSpecsRef ? scalar @{$NbrOfNbrAtomSpecsRef} : 0;
    3124
    3125   # Validate number of specifications...
    3126   if ($NumOfNbrBondSpecs && ($NumOfNbrAtomSpecs != $NumOfNbrBondSpecs)) {
    3127     carp "Warning: ${ClassName}->DoesAtomNeighborhoodMatch: Number of specified central atom, $NumOfNbrAtomSpecs, and bond, $NumOfNbrBondSpecs, specifications must be same; No neighborhood match performed ...";
    3128     return 0;
    3129   }
    3130
    3131   if ($NumOfNbrOfNbrAtomSpecs && ($NumOfNbrOfNbrAtomSpecs != $NumOfNbrAtomSpecs)) {
    3132     carp "Warning: ${ClassName}->DoesAtomNeighborhoodMatch: Number of specified central atom, $NumOfNbrAtomSpecs, and neighbor of neighbor atoms specifications, $NumOfNbrOfNbrAtomSpecs, must be same; No neighborhood match performed ...";
    3133     return 0;
    3134   }
    3135
    3136   # Sort atom and bond specifications in terms of being most specific to least specific..
    3137   ($NbrAtomSpecsRef, $NbrBondSpecsRef, $NbrOfNbrAtomSpecsRef) = $CentralAtom->_SortSpecificationsForAtomNeighborhoodMatch($NbrAtomSpecsRef, $NbrBondSpecsRef, $NbrOfNbrAtomSpecsRef);
    3138
    3139   # Does central atom specification match?
    3140   if (!$CentralAtom->_DoesAtomSpecificationMatch($CentralAtomSpec)) {
    3141     return 0;
    3142   }
    3143
    3144   # No neighbors to match...
    3145   if (!$NumOfNbrAtomSpecs) {
    3146     return 1;
    3147   }
    3148
    3149   # Match neighbors...
    3150   my($NbrSpecsMatched, $NbrSpecCount, $NbrSpecMatchCount, %NbrSpecAlreadyMatchedMap);
    3151
    3152   $NbrSpecCount = $NumOfNbrAtomSpecs;
    3153   $NbrSpecMatchCount = 0;
    3154
    3155   %NbrSpecAlreadyMatchedMap = ();
    3156   ($NbrSpecsMatched, $NbrSpecMatchCount) = $CentralAtom->_MatchAtomNeighborhoodUsingAtomBondSpecs($NbrSpecCount, $NbrSpecMatchCount, \%NbrSpecAlreadyMatchedMap, $NbrAtomSpecsRef, $NbrBondSpecsRef, $NbrOfNbrAtomSpecsRef);
    3157
    3158   if ($NbrSpecsMatched) {
    3159     # It's match...
    3160     return 1;
    3161   }
    3162
    3163   # Match central atom's missing hydrogens with any unmatched atom
    3164   # and bond specifications...
    3165   #
    3166   ($NbrSpecsMatched, $NbrSpecMatchCount) = $CentralAtom->_MatchAtomNeighborhoodUsingMissingHydrogens($NbrSpecCount, $NbrSpecMatchCount, \%NbrSpecAlreadyMatchedMap, $NbrAtomSpecsRef, $NbrBondSpecsRef, $NbrOfNbrAtomSpecsRef);
    3167
    3168   if ($NbrSpecsMatched) {
    3169     # It's match...
    3170     return 1;
    3171   }
    3172
    3173   # No match...
    3174   return 0;
    3175 }
    3176
    3177 # Match central atom neighborhood atom and bond specifications...
    3178 #
    3179 sub _MatchAtomNeighborhoodUsingAtomBondSpecs {
    3180   my($CentralAtom, $NbrSpecCount, $NbrSpecMatchCount, $NbrSpecAlreadyMatchedRef, $NbrAtomSpecsRef, $NbrBondSpecsRef, $NbrOfNbrAtomSpecsRef) = @_;
    3181   my($Index, $NbrAtom, $NbrAtomSpec, $NbrBondSpec, $NbrOfNbrAtom, $NbrOfNbrAtomSpec, $MatchNbrOfNbrAtomSpecs, $NbrSpecsMatched);
    3182
    3183   $MatchNbrOfNbrAtomSpecs = (defined $NbrOfNbrAtomSpecsRef && scalar @{$NbrOfNbrAtomSpecsRef}) ? 1 : 0;
    3184
    3185   $NbrSpecsMatched = 0;
    3186
    3187   # Match central atom's  immediate neighbors atom and bond specifications...
    3188   NBRATOM:  for $NbrAtom ($CentralAtom->GetNeighbors()) {
    3189     NBRATOMSPEC: for $Index (0 .. ($NbrSpecCount - 1)) {
    3190       if (exists $NbrSpecAlreadyMatchedRef->{$Index}) {
    3191         next NBRATOMSPEC;
    3192       }
    3193       $NbrAtomSpec = $NbrAtomSpecsRef->[$Index];
    3194       $NbrBondSpec = $NbrBondSpecsRef->[$Index];
    3195
    3196       $NbrOfNbrAtomSpec = $MatchNbrOfNbrAtomSpecs ? $NbrOfNbrAtomSpecsRef->[$Index] : undef;
    3197
    3198       # Match neighbor atom specification...
    3199       if (!$NbrAtom->_DoesAtomSpecificationMatch($NbrAtomSpec)) {
    3200         next NBRATOMSPEC;
    3201       }
    3202
    3203       # Match central atom to neighbor atom bond specification...
    3204       if (!$CentralAtom->_DoesBondSpecificationMatch($NbrAtom, $NbrBondSpec)) {
    3205         next NBRATOMSPEC;
    3206       }
    3207
    3208       # Match any neighbor of neighbor atom specifications...
    3209       if (defined $NbrOfNbrAtomSpec) {
    3210         # Go over the neighbors of central atom skipping the central atom...
    3211         for $NbrOfNbrAtom ($NbrAtom->GetNeighbors($CentralAtom)) {
    3212           if (!$NbrOfNbrAtom->_DoesAtomSpecificationMatch($NbrOfNbrAtomSpec)) {
    3213             next NBRATOMSPEC;
    3214           }
    3215         }
    3216       }
    3217
    3218       # It's a match for a neighbor atom specification...
    3219       $NbrSpecAlreadyMatchedRef->{$Index} = $Index;
    3220       $NbrSpecMatchCount++;
    3221
    3222       if ($NbrSpecMatchCount == $NbrSpecCount) {
    3223         # It's match...
    3224         $NbrSpecsMatched = 1;
    3225         last NBRATOM;
    3226       }
    3227       # Match next neighbor atom...
    3228       next NBRATOM;
    3229     }
    3230   }
    3231   return ($NbrSpecsMatched, $NbrSpecMatchCount);
    3232 }
    3233
    3234 # Match central atom's missing hydrogens with any unmatched atom and bond
    3235 # specifications...
    3236 #
    3237 sub _MatchAtomNeighborhoodUsingMissingHydrogens {
    3238   my($CentralAtom, $NbrSpecCount, $NbrSpecMatchCount, $NbrSpecAlreadyMatchedRef, $NbrAtomSpecsRef, $NbrBondSpecsRef, $NbrOfNbrAtomSpecsRef) = @_;
    3239   my($Index, $NbrAtom, $NbrAtomSpec, $NbrBondSpec, $NumOfMissingHydrogens, $MissingHydrogensIndex, $NbrSpecsMatched, $AtomSpecMatched, $AtomSpec, $AtomSymbol);
    3240
    3241   $NbrSpecsMatched = 0;
    3242
    3243   $NumOfMissingHydrogens = $CentralAtom->GetNumOfMissingHydrogens();
    3244   if (($NbrSpecCount - $NbrSpecMatchCount) > $NumOfMissingHydrogens) {
    3245     # It won't match...
    3246     return ($NbrSpecsMatched, $NbrSpecMatchCount);
    3247   }
    3248
    3249   MISSINGHYDROGENNBR: for $MissingHydrogensIndex (0 .. ($NumOfMissingHydrogens - 1)) {
    3250     NBRATOMSPEC: for $Index (0 .. ($NbrSpecCount - 1)) {
    3251       if (exists $NbrSpecAlreadyMatchedRef->{$Index}) {
    3252         next NBRATOMSPEC;
    3253       }
    3254       $NbrAtomSpec = $NbrAtomSpecsRef->[$Index];
    3255       $NbrBondSpec = $NbrBondSpecsRef->[$Index];
    3256
    3257       $NbrAtomSpec =~ s/ //g;
    3258
    3259       # Match neighbor atom specification hydrogen atom symbol...
    3260       $AtomSpecMatched = 0;
    3261       ATOMSPEC: for $AtomSpec (split /\,/, $NbrAtomSpec) {
    3262         ($AtomSymbol) = split /\./, $AtomSpec;
    3263         if ($AtomSymbol =~ /^(H|A|\*)$/i) {
    3264           $AtomSpecMatched = 1;
    3265           last ATOMSPEC;
    3266         }
    3267       }
    3268       if (!$AtomSpecMatched) {
    3269         next NBRATOMSPEC;
    3270       }
    3271
    3272       # Match neighbor atom bond specification to singal bond...
    3273       if (defined $NbrBondSpec) {
    3274         $NbrBondSpec =~ s/ //g;
    3275         if ($NbrBondSpec !~ /^(-|1|s|Single|\~|\*|Any)/i) {
    3276           next NBRATOMSPEC;
    3277         }
    3278       }
    3279
    3280       # It's a match for a neighbor atom specification...
    3281       $NbrSpecAlreadyMatchedRef->{$Index} = $Index;
    3282       $NbrSpecMatchCount++;
    3283
    3284       if ($NbrSpecMatchCount == $NbrSpecCount) {
    3285         # It's match...
    3286         $NbrSpecsMatched = 1;
    3287         last MISSINGHYDROGENNBR;
    3288       }
    3289       # Match next missing hydrogen neighbor...
    3290       next MISSINGHYDROGENNBR;
    3291     }
    3292   }
    3293
    3294   return ($NbrSpecsMatched, $NbrSpecMatchCount);
    3295 }
    3296
    3297 # Sort atom and bond specifications base on neighbor atom specifications going
    3298 # from most to least specific atom specifications.
    3299 #
    3300 # Atom specifications are sorted at the following two levels:
    3301 #
    3302 #   o By atom specification count with in each specification going from most specific
    3303 #      to least specific, where count is determined by the number of "," in each
    3304 #      specification. Wild card containing specifications are considered least specific
    3305 #      and end up at the end of the sorted list.
    3306 #   o By atomic invariant count with in each sorted list going from most specific to
    3307 #      least specific, where count is determined by the number of "." in each atom
    3308 #      specification.
    3309 #
    3310 #A single atom specification,
    3311 # without any commas in atom specification, is is considered most specific...
    3312 #
    3313 sub _SortSpecificationsForAtomNeighborhoodMatch {
    3314   my($This, $NbrAtomSpecsRef, $NbrBondSpecsRef, $NbrOfNbrAtomSpecsRef) = @_;
    3315   my($Index, $NeedToSort, $NumOfNbrAtomSpecs, $NbrAtomSpecCount,  $NbrAtomSpecAtomicInvarintCount, $NbrAtomSpecToMatch, $NbrAtomSpec, $FirstAtomicInvariant, $WildCardInNbrAtomSpec, @NbrAtomSpecs, @NbrAtomSpecAtomicInvariants, @SortedNbrAtomSpecs, @SortedNbrBondSpecs, @SortedNbrOfNbrAtomSpecs, %NbrAtomSpecDataMap);
    3316
    3317   $NumOfNbrAtomSpecs = defined $NbrAtomSpecsRef ? scalar @{$NbrAtomSpecsRef} : 0;
    3318
    3319   # Figure out whether sorting is necessary...
    3320   $NeedToSort = 0;
    3321   if ($NumOfNbrAtomSpecs > 1) {
    3322     ATOMSPEC: for $NbrAtomSpecToMatch (@{$NbrAtomSpecsRef}) {
    3323       if ($NbrAtomSpecToMatch =~ /(,|\.|A|\*)/i) {
    3324         $NeedToSort = 1;
    3325         last ATOMSPEC;
    3326       }
    3327     }
    3328   }
    3329   if (!$NeedToSort) {
    3330     # Nothing to do...
    3331     return ($NbrAtomSpecsRef, $NbrBondSpecsRef, $NbrOfNbrAtomSpecsRef);
    3332   }
    3333
    3334   %NbrAtomSpecDataMap = ();
    3335
    3336   for $Index (0 .. ($NumOfNbrAtomSpecs - 1)) {
    3337     $NbrAtomSpecToMatch = $NbrAtomSpecsRef->[$Index];
    3338     $NbrAtomSpecToMatch =~ s/ //g;
    3339
    3340     @NbrAtomSpecs = split /\,/, $NbrAtomSpecToMatch;
    3341     $NbrAtomSpecCount = scalar @NbrAtomSpecs;
    3342
    3343     # Does neighbor specification contains a wild card in atom symbol specification?
    3344     #
    3345     if ($NbrAtomSpecToMatch =~ /(A|\*)/i) {
    3346       $WildCardInNbrAtomSpec = 0;
    3347       NBRATOMSPEC: for $NbrAtomSpec (@NbrAtomSpecs) {
    3348         ($FirstAtomicInvariant) = split /\./, $NbrAtomSpec;
    3349         if ($FirstAtomicInvariant =~ /^!/) {
    3350           $FirstAtomicInvariant =~ s/^!//;
    3351         }
    3352         $WildCardInNbrAtomSpec = $This->_IsWildCardAtomSymbolAtomicInvariant($FirstAtomicInvariant);
    3353         if ($WildCardInNbrAtomSpec) {
    3354           last NBRATOMSPEC;
    3355         }
    3356       }
    3357       if ($WildCardInNbrAtomSpec) {
    3358         # Set NbrAtomSpecCount arbitrarily high to make the spec containing wild
    3359         # card last on the sorted list while maintaining its original order in the list...
    3360         $NbrAtomSpecCount = 999;
    3361       }
    3362     }
    3363
    3364     if (!exists $NbrAtomSpecDataMap{$NbrAtomSpecCount}) {
    3365       %{$NbrAtomSpecDataMap{$NbrAtomSpecCount}} = ();
    3366     }
    3367
    3368     # Use first NbrAtomSpec available in @NbrAtomSpecs to determine atomic invariant count
    3369     # with in each NbrAtomSpecToMatch, as @NbrAtomSpecs derived from $NbrAtomSpecToMatch
    3370     # simply corresponds to a list of possible matches...
    3371     #
    3372     ($NbrAtomSpec) = @NbrAtomSpecs;
    3373     @NbrAtomSpecAtomicInvariants = split /\./, $NbrAtomSpec;
    3374     $NbrAtomSpecAtomicInvarintCount = scalar @NbrAtomSpecAtomicInvariants;
    3375
    3376     if (!exists $NbrAtomSpecDataMap{$NbrAtomSpecCount}{$NbrAtomSpecAtomicInvarintCount}) {
    3377       @{$NbrAtomSpecDataMap{$NbrAtomSpecCount}{$NbrAtomSpecAtomicInvarintCount}} = ();
    3378     }
    3379     push @{$NbrAtomSpecDataMap{$NbrAtomSpecCount}{$NbrAtomSpecAtomicInvarintCount}}, $Index;
    3380
    3381   }
    3382
    3383   @SortedNbrAtomSpecs = (); @SortedNbrBondSpecs = ();
    3384   @SortedNbrOfNbrAtomSpecs = ();
    3385
    3386   for $NbrAtomSpecCount ( sort { $a <=> $b } keys %NbrAtomSpecDataMap) {
    3387     for $NbrAtomSpecAtomicInvarintCount ( sort { $b <=> $a } keys %{$NbrAtomSpecDataMap{$NbrAtomSpecCount}}) {
    3388       for $Index (@{$NbrAtomSpecDataMap{$NbrAtomSpecCount}{$NbrAtomSpecAtomicInvarintCount}}) {
    3389         push @SortedNbrAtomSpecs, $NbrAtomSpecsRef->[$Index];
    3390         if (defined $NbrBondSpecsRef) {
    3391           push @SortedNbrBondSpecs, $NbrBondSpecsRef->[$Index];
    3392         }
    3393         if (defined $NbrOfNbrAtomSpecsRef) {
    3394           push @SortedNbrOfNbrAtomSpecs, $NbrOfNbrAtomSpecsRef->[$Index];
    3395         }
    3396       }
    3397     }
    3398   }
    3399
    3400   return (\@SortedNbrAtomSpecs, defined $NbrBondSpecsRef ? \@SortedNbrBondSpecs : undef, defined $NbrOfNbrAtomSpecsRef ? \@SortedNbrOfNbrAtomSpecs : undef);
    3401 }
    3402
    3403 # Check whether atom matches supported atom specification...
    3404 #
    3405 sub _DoesAtomSpecificationMatch {
    3406   my($This, $AtomSpecificationToMatch) = @_;
    3407   my($AtomSpecification, $AtomicInvariant, $AtomSpecificationMatched, $AtomicInvariantMatched, $FirstMatch);
    3408
    3409   # Anything to match...
    3410   if (!(defined($AtomSpecificationToMatch) && $AtomSpecificationToMatch)) {
    3411     return 1;
    3412   }
    3413
    3414   # Take out any spaces...
    3415   $AtomSpecificationToMatch =~ s/ //g;
    3416
    3417   # Match specified atom specifications. For multiple atom specifications in a comma delimited string,
    3418   # only one atom specification needs to match for a successful match. It's up to the caller to make
    3419   # sure that the specificaton list is ordered from least to most specific atom specification...
    3420   #
    3421   for $AtomSpecification (split /\,/, $AtomSpecificationToMatch) {
    3422     $AtomSpecificationMatched = 1;
    3423     $FirstMatch = 1;
    3424
    3425     # Match all atom symbol atomic invariants...
    3426     ATOMICINVARIANT: for $AtomicInvariant (split /\./, $AtomSpecification) {
    3427       if ($FirstMatch) {
    3428         # Match atom symbol atomic invariant...
    3429         $FirstMatch = 0;
    3430         $AtomicInvariantMatched = $This->_MatchAtomSymbolAtomicInvariant($AtomicInvariant);
    3431       }
    3432       else {
    3433         # Match non atom symbol atomic invariant...
    3434         $AtomicInvariantMatched = $This->_MatchNonAtomSymbolAtomicInvariant($AtomicInvariant);
    3435       }
    3436
    3437       if (!$AtomicInvariantMatched) {
    3438         # No need to match other atomic invariants...
    3439         $AtomSpecificationMatched = 0;
    3440         last ATOMICINVARIANT;
    3441       }
    3442     }
    3443
    3444     if ($AtomSpecificationMatched) {
    3445       # No need to match other atom specifications...
    3446       return 1;
    3447     }
    3448   }
    3449
    3450   # Nothing matched...
    3451   return 0;
    3452 }
    3453
    3454 # Check whether atom matches atom symbol atomic invariant...
    3455 #
    3456 sub _MatchAtomSymbolAtomicInvariant {
    3457   my($This, $AtomicInvariant) = @_;
    3458   my($NegateMatch, $Status, $AtomicNumber);
    3459
    3460   $Status = 0;
    3461   $NegateMatch = 0;
    3462
    3463   # Does match needs to be negated?
    3464   if ($AtomicInvariant =~ /^!/) {
    3465     $NegateMatch = 1;
    3466     $AtomicInvariant =~ s/^!//;
    3467   }
    3468
    3469   ATOMICINVARIANT: {
    3470     # Any atom match...
    3471     if ($This->_IsWildCardAtomSymbolAtomicInvariant($AtomicInvariant)) {
    3472       $Status = 1;
    3473       last ATOMICINVARIANT;
    3474     }
    3475
    3476     # Atomic number match...
    3477     if ($AtomicInvariant =~ /^#/) {
    3478       $AtomicNumber = $AtomicInvariant; $AtomicNumber =~ s/^#//;
    3479       $Status = ($This->{AtomicNumber} == $AtomicNumber) ? 1 : 0;
    3480       last ATOMICINVARIANT;
    3481     }
    3482
    3483     # Atom symbol match...
    3484     $Status = ($This->{AtomSymbol} =~ /^$AtomicInvariant$/i) ? 1 : 0;
    3485   }
    3486
    3487   if ($NegateMatch) {
    3488     $Status = $Status ? 0 : 1;
    3489   }
    3490
    3491   return $Status;
    3492 }
    3493
    3494 # Is it a wild card atom symbol atomic invariant?
    3495 #
    3496 sub _IsWildCardAtomSymbolAtomicInvariant {
    3497   my($This, $AtomicInvariant) = @_;
    3498
    3499   return ($AtomicInvariant =~ /^(A|\*)$/i) ? 1 : 0;
    3500 }
    3501
    3502 # Check whether atom matches non atom symbol atomic invariants...
    3503 #
    3504 sub _MatchNonAtomSymbolAtomicInvariant {
    3505   my($This, $AtomicInvariant) = @_;
    3506   my($NegateMatch, $Status, $Name, $Value, $UnknownName);
    3507
    3508   ($Status, $NegateMatch, $UnknownName) = ('0') x 3;
    3509
    3510   # Does match needs to be negated?
    3511   if ($AtomicInvariant =~ /^!/) {
    3512     $NegateMatch = 1;
    3513     $AtomicInvariant =~ s/^!//;
    3514   }
    3515
    3516   # Extract atomic invariant name and any value...
    3517   if ($AtomicInvariant =~ /[0-9\*]+/) {
    3518     ($Name, $Value) = $AtomicInvariant =~ /^([a-zA-Z]+)([0-9\-\+\*\>\<\=]+)$/;
    3519   }
    3520   else {
    3521     ($Name, $Value) = ($AtomicInvariant, undef);
    3522   }
    3523
    3524   NAME: {
    3525     # Match number of non-hydrogen atom neighbors
    3526     if ($Name =~ /^X$/i) {
    3527       $Status = (defined($Value) && $This->GetNumOfNonHydrogenAtomNeighbors() == $Value) ? 1 : 0;
    3528       last NAME;
    3529     }
    3530
    3531     # Match total number of atom neighbors including missing hydrogens...
    3532     if ($Name =~ /^T$/i) {
    3533       $Status = (defined($Value) && ($This->GetNumOfNonHydrogenAtomNeighbors() + $This->GetNumOfHydrogens()) == $Value) ? 1 : 0;
    3534       last NAME;
    3535     }
    3536
    3537     # Match formal charge...
    3538     if ($Name =~ /^FC$/i) {
    3539       my $FormalCharge = $This->GetFormalCharge();
    3540       $Status = $This->_MatchNonAtomSymbolAtomicInvariantValue($FormalCharge, $Value);
    3541       last NAME;
    3542     }
    3543
    3544     # Match aromatic annotation indicating whether atom is aromatic...
    3545     if ($Name =~ /^Ar$/i) {
    3546       $Status = $This->IsAromatic() ? 1 : 0;
    3547       last NAME;
    3548     }
    3549
    3550     # Match number of implicit and explicit hydrogens...
    3551     if ($Name =~ /^H$/i) {
    3552       $Status = (defined($Value) && ($This->GetNumOfHydrogens() == $Value)) ? 1 : 0;
    3553       last NAME;
    3554     }
    3555
    3556     # Match ring atom annotation indicating whether atom is in ring...
    3557     if ($Name =~ /^RA$/i) {
    3558       $Status = defined($Value) ? $This->IsInRingOfSize($Value) : ($This->IsInRing() ? 1 : 0);
    3559       last NAME;
    3560     }
    3561
    3562     # Match number of rings for atom..
    3563     if ($Name =~ /^TR$/i) {
    3564       $Status = (defined($Value) && ($Value == $This->GetNumOfRings())) ? 1 : 0;
    3565       last NAME;
    3566     }
    3567
    3568     # Match sum of bond orders to non-hydrogen atom neighbors...
    3569     if ($Name =~ /^BO$/i) {
    3570       $Status = (defined($Value) && $This->GetSumOfBondOrdersToNonHydrogenAtoms() == $Value) ? 1 : 0;
    3571       last NAME;
    3572     }
    3573
    3574     # Match largest bond order of non-hydrogen atom neighbors...
    3575     if ($Name =~ /^LBO$/i) {
    3576       $Status = (defined($Value) && $This->GetLargestBondOrderToNonHydrogenAtoms() == $Value) ? 1 : 0;
    3577       last NAME;
    3578     }
    3579
    3580     # Match number of single bonds to non-hydrogen atom neighbors...
    3581     if ($Name =~ /^SB$/i) {
    3582       $Status = (defined($Value) && $This->GetNumOfSingleBondsToNonHydrogenAtoms() == $Value) ? 1 : 0;
    3583       last NAME;
    3584     }
    3585
    3586     # Match total number of single bonds to atom neighbors including missing and explicit hydrogens...
    3587     if ($Name =~ /^TSB$/i) {
    3588       $Status = (defined($Value) && ($This->GetNumOfSingleBondsToNonHydrogenAtoms() + $This->GetNumOfHydrogens()) == $Value) ? 1 : 0;
    3589       last NAME;
    3590     }
    3591
    3592     # Match number of double bonds to non-hydrogen atom neighbors...
    3593     if ($Name =~ /^DB$/i) {
    3594       $Status = (defined($Value) && $This->GetNumOfDoubleBondsToNonHydrogenAtoms() == $Value) ? 1 : 0;
    3595       last NAME;
    3596     }
    3597
    3598     # Match number of triple bonds to non-hydrogen atom neighbors...
    3599     if ($Name =~ /^TB$/i) {
    3600       $Status = (defined($Value) && $This->GetNumOfTripleBondsToNonHydrogenAtoms() == $Value) ? 1 : 0;
    3601       last NAME;
    3602     }
    3603
    3604     # Match number of aromatic bonds to non-hydrogen atom neighbors...
    3605     if ($Name =~ /^AB$/i) {
    3606       $Status = (defined($Value) && $This->GetNumOfAromaticBondsToNonHydrogenAtoms() == $Value) ? 1 : 0;
    3607       last NAME;
    3608     }
    3609
    3610
    3611     # Match mass number indicating isotope other than most abundant isotope...
    3612     if ($Name =~ /^MN$/i) {
    3613       $Status = (defined($Value) && $This->GetMassNumber() == $Value) ? 1 : 0;
    3614       last NAME;
    3615     }
    3616
    3617     # Match spin multiplicity...
    3618     if ($Name =~ /^SM$/i) {
    3619       my $SpinMultiplicity = $This->GetSpinMultiplicity();
    3620       if (!defined $SpinMultiplicity) { $SpinMultiplicity = 0; }
    3621       $Status = (defined($Value) && defined($SpinMultiplicity) && $Value == $SpinMultiplicity) ? 1 : 0;
    3622       last NAME;
    3623     }
    3624
    3625     $UnknownName = 1;
    3626     carp "Warning: ${ClassName}->_MatchNonAtomSymbolAtomicInvariant: Unknown atomic invariant $AtomicInvariant...";
    3627   }
    3628
    3629   if (!$UnknownName) {
    3630     if ($NegateMatch) {
    3631       $Status = $Status ? 0 : 1;
    3632     }
    3633   }
    3634
    3635   return $Status;
    3636 }
    3637
    3638 # Match atomic invariant value...
    3639 #
    3640 # Specified value format:
    3641 #   . +* : Any positive value
    3642 #   . -* : Any negative value
    3643 #   . >ValidNumber or >=ValidNumber
    3644 #   . <ValidNumber or <=ValidNumber
    3645 #   . Any valid number
    3646 #
    3647 sub _MatchNonAtomSymbolAtomicInvariantValue {
    3648   my($This, $TargetValue, $SpecifiedValue) = @_;
    3649   my($Status);
    3650
    3651   $Status = 0;
    3652
    3653   if (!(defined($TargetValue) && defined($SpecifiedValue))) {
    3654     return $Status;
    3655   }
    3656
    3657   VALUE: {
    3658     if ($SpecifiedValue =~ /^\+\*/) {
    3659       $Status = ($TargetValue > 0) ? 1 : 0;
    3660       last VALUE;
    3661     }
    3662     if ($SpecifiedValue =~ /^\-\*/) {
    3663       $Status = ($TargetValue < 0) ? 1 : 0;
    3664       last VALUE;
    3665     }
    3666     if ($SpecifiedValue =~ /^>/) {
    3667       if ($SpecifiedValue =~ /^>=/) {
    3668         $SpecifiedValue =~ s/^>=//;
    3669         $Status = ($SpecifiedValue >= $TargetValue) ? 1 : 0;
    3670       }
    3671       else {
    3672         $SpecifiedValue =~ s/^>//;
    3673         $Status = ($SpecifiedValue > $TargetValue) ? 1 : 0;
    3674       }
    3675       last VALUE;
    3676     }
    3677     if ($SpecifiedValue =~ /^</) {
    3678       if ($SpecifiedValue =~ /^<=/) {
    3679         $SpecifiedValue =~ s/^<=//;
    3680         $Status = ($SpecifiedValue <= $TargetValue) ? 1 : 0;
    3681       }
    3682       else {
    3683         $SpecifiedValue =~ s/^<//;
    3684         $Status = ($SpecifiedValue < $TargetValue) ? 1 : 0;
    3685       }
    3686       last VALUE;
    3687     }
    3688     # Default is do perform an equality match...
    3689     $Status = ($SpecifiedValue == $TargetValue) ? 1 : 0;
    3690   }
    3691
    3692   return $Status;
    3693 }
    3694
    3695 # Check whether atoms match  bond specifications...
    3696 #
    3697 sub _DoesBondSpecificationMatch {
    3698   my($This, $BondedAtom, $BondSpecificationToMatch) = @_;
    3699   my($BondSpecification, $BondSymbolSpecification, $BondSpecificationMatched);
    3700
    3701   # Anything to match...
    3702   if (!(defined($BondSpecificationToMatch) && $BondSpecificationToMatch)) {
    3703     return 1;
    3704   }
    3705
    3706   # Take out any spaces...
    3707   $BondSpecificationToMatch =~ s/ //g;
    3708
    3709   # Match specified bond specifications. For multiple bond specifications in a comma delimited string,
    3710   # only one bond specification needs to match for a successful match...
    3711   #
    3712   for $BondSpecification (split /\,/, $BondSpecificationToMatch) {
    3713     $BondSpecificationMatched = 1;
    3714
    3715     # Match all specified bond symbol specifications...
    3716     BONDSYMBOL: for $BondSymbolSpecification (split /\./, $BondSpecification) {
    3717       if (!$This->_MatchBondSymbolSpecification($BondedAtom, $BondSymbolSpecification)) {
    3718         # No need to match other bond symbol specifications...
    3719         $BondSpecificationMatched = 0;
    3720         last BONDSYMBOL;
    3721       }
    3722     }
    3723     if ($BondSpecificationMatched) {
    3724       # No need to try matching other bond specifications...
    3725       return 1;
    3726     }
    3727   }
    3728
    3729   # Nothing matched...
    3730   return 0;
    3731 }
    3732
    3733 # Check whether atoms match  bond symbol specification...
    3734 #
    3735 sub _MatchBondSymbolSpecification {
    3736   my($This, $BondedAtom, $BondSymbolSpecification) = @_;
    3737   my($NegateMatch, $Status, $Bond, $BondSymbol, $UnknownBondSymbol);
    3738
    3739   ($Status, $NegateMatch, $UnknownBondSymbol) = ('0') x 3;
    3740
    3741   # Does match needs to be negated?
    3742   if ($BondSymbolSpecification =~ /^!/) {
    3743     $NegateMatch = 1;
    3744     $BondSymbolSpecification =~ s/^!//;
    3745   }
    3746   $BondSymbol = $BondSymbolSpecification;
    3747   $Bond = $This->GetBondToAtom($BondedAtom);
    3748
    3749   BONDSYMBOL: {
    3750     if ($BondSymbol =~ /^(-|1|s|Single)$/i) { $Status = $Bond->IsSingle() ? 1 : 0; last BONDSYMBOL; }
    3751     if ($BondSymbol =~ /^(=|2|d|Double)$/i) { $Status = $Bond->IsDouble() ? 1 : 0; last BONDSYMBOL; }
    3752     if ($BondSymbol =~ /^(#|3|t|Triple)$/i) { $Status = $Bond->IsTriple() ? 1 : 0; last BONDSYMBOL; }
    3753     if ($BondSymbol =~ /^(:|a|Ar|Aromatic)$/i) { $Status = $Bond->IsAromatic() ? 1 : 0; last BONDSYMBOL; }
    3754
    3755     if ($BondSymbol =~ /^(\@|RB|Ring)$/i) { $Status = $Bond->IsInRing() ? 1 : 0; last BONDSYMBOL; }
    3756
    3757     if ($BondSymbol =~ /^(\~|\*|Any)$/i) { $Status = 1; last BONDSYMBOL; }
    3758
    3759     $UnknownBondSymbol = 1;
    3760     carp "Warning: ${ClassName}->_MatchBondSpecification: Unknown bond specification $BondSymbolSpecification...";
    3761   }
    3762
    3763   if (!$UnknownBondSymbol) {
    3764     if ($NegateMatch) {
    3765       $Status = $Status ? 0 : 1;
    3766     }
    3767   }
    3768
    3769   return $Status;
    3770 }
    3771
    3772 # Is it a saturated atom?
    3773 #
    3774 sub IsSaturated {
    3775   my($This) = @_;
    3776
    3777   return !$This->IsUnsaturated();
    3778 }
    3779
    3780 # Is it an unsaturated atom containing at least one non-single bond?
    3781 #
    3782 sub IsUnsaturated {
    3783   my($This) = @_;
    3784   my($NumOfSingleBonds, $NumOfDoubleBonds, $NumOfTripleBonds, $NumOfAromaticBonds);
    3785
    3786   ($NumOfSingleBonds, $NumOfDoubleBonds, $NumOfTripleBonds, $NumOfAromaticBonds) = $This->GetNumOfBondTypesToNonHydrogenAtoms();
    3787
    3788   return ($NumOfDoubleBonds || $NumOfTripleBonds || $NumOfAromaticBonds) ? 1 : 0;
    3789 }
    3790
    3791 # Is atom in a ring?
    3792 #
    3793 sub IsInRing {
    3794   my($This) = @_;
    3795
    3796   # Is this atom in a molecule?
    3797   if (!$This->HasProperty('Molecule')) {
    3798     return undef;
    3799   }
    3800   my($Molecule);
    3801   $Molecule = $This->GetProperty('Molecule');
    3802
    3803   return $Molecule->_IsAtomInRing($This);
    3804 }
    3805
    3806 # Is atom not in a ring?
    3807 #
    3808 sub IsNotInRing {
    3809   my($This) = @_;
    3810
    3811   # Is this atom in a molecule?
    3812   if (!$This->HasProperty('Molecule')) {
    3813     return undef;
    3814   }
    3815   my($Molecule);
    3816   $Molecule = $This->GetProperty('Molecule');
    3817
    3818   return $Molecule->_IsAtomNotInRing($This);
    3819 }
    3820
    3821 # Is atom only in one ring?
    3822 #
    3823 sub IsOnlyInOneRing {
    3824   my($This) = @_;
    3825
    3826   # Is this atom in a molecule?
    3827   if (!$This->HasProperty('Molecule')) {
    3828     return undef;
    3829   }
    3830   my($Molecule);
    3831   $Molecule = $This->GetProperty('Molecule');
    3832
    3833   return $Molecule->_IsAtomInOnlyOneRing($This);
    3834 }
    3835
    3836 # Is atom in a ring of specific size?
    3837 #
    3838 sub IsInRingOfSize {
    3839   my($This, $RingSize) = @_;
    3840
    3841   # Is this atom in a molecule?
    3842   if (!$This->HasProperty('Molecule')) {
    3843     return undef;
    3844   }
    3845   my($Molecule);
    3846   $Molecule = $This->GetProperty('Molecule');
    3847
    3848   return $Molecule->_IsAtomInRingOfSize($This, $RingSize);
    3849 }
    3850
    3851 # Get size of smallest ring containing the atom...
    3852 #
    3853 sub GetSizeOfSmallestRing {
    3854   my($This) = @_;
    3855
    3856   # Is this atom in a molecule?
    3857   if (!$This->HasProperty('Molecule')) {
    3858     return undef;
    3859   }
    3860   my($Molecule);
    3861   $Molecule = $This->GetProperty('Molecule');
    3862
    3863   return $Molecule->_GetSizeOfSmallestAtomRing($This);
    3864 }
    3865
    3866 # Get size of largest ring containing the atom...
    3867 #
    3868 sub GetSizeOfLargestRing {
    3869   my($This) = @_;
    3870
    3871   # Is this atom in a molecule?
    3872   if (!$This->HasProperty('Molecule')) {
    3873     return undef;
    3874   }
    3875   my($Molecule);
    3876   $Molecule = $This->GetProperty('Molecule');
    3877
    3878   return $Molecule->_GetSizeOfLargestAtomRing($This);
    3879 }
    3880
    3881 # Get number of  rings containing the atom...
    3882 #
    3883 sub GetNumOfRings {
    3884   my($This) = @_;
    3885
    3886   # Is this atom in a molecule?
    3887   if (!$This->HasProperty('Molecule')) {
    3888     return undef;
    3889   }
    3890   my($Molecule);
    3891   $Molecule = $This->GetProperty('Molecule');
    3892
    3893   return $Molecule->_GetNumOfAtomRings($This);
    3894 }
    3895
    3896 # Get number of  rings with odd size containing the atom...
    3897 #
    3898 sub GetNumOfRingsWithOddSize {
    3899   my($This) = @_;
    3900
    3901   # Is this atom in a molecule?
    3902   if (!$This->HasProperty('Molecule')) {
    3903     return undef;
    3904   }
    3905   my($Molecule);
    3906   $Molecule = $This->GetProperty('Molecule');
    3907
    3908   return $Molecule->_GetNumOfAtomRingsWithOddSize($This);
    3909 }
    3910
    3911 # Get number of  rings with even size containing the atom...
    3912 #
    3913 sub GetNumOfRingsWithEvenSize {
    3914   my($This) = @_;
    3915
    3916   # Is this atom in a molecule?
    3917   if (!$This->HasProperty('Molecule')) {
    3918     return undef;
    3919   }
    3920   my($Molecule);
    3921   $Molecule = $This->GetProperty('Molecule');
    3922
    3923   return $Molecule->_GetNumOfAtomRingsWithEvenSize($This);
    3924 }
    3925
    3926 # Get number of  rings with specified size containing the atom...
    3927 #
    3928 sub GetNumOfRingsWithSize {
    3929   my($This, $RingSize) = @_;
    3930
    3931   # Is this atom in a molecule?
    3932   if (!$This->HasProperty('Molecule')) {
    3933     return undef;
    3934   }
    3935   my($Molecule);
    3936   $Molecule = $This->GetProperty('Molecule');
    3937
    3938   return $Molecule->_GetNumOfAtomRingsWithSize($This, $RingSize);
    3939
    3940 }
    3941
    3942 # Get number of  rings with size less than specified containing the atom...
    3943 #
    3944 sub GetNumOfRingsWithSizeLessThan {
    3945   my($This, $RingSize) = @_;
    3946
    3947   # Is this atom in a molecule?
    3948   if (!$This->HasProperty('Molecule')) {
    3949     return undef;
    3950   }
    3951   my($Molecule);
    3952   $Molecule = $This->GetProperty('Molecule');
    3953
    3954   return $Molecule->_GetNumOfAtomRingsWithSizeLessThan($This, $RingSize);
    3955 }
    3956
    3957 # Get number of  rings with size greater than specified size containing the atom...
    3958 #
    3959 sub GetNumOfRingsWithSizeGreaterThan {
    3960   my($This, $RingSize) = @_;
    3961
    3962   # Is this atom in a molecule?
    3963   if (!$This->HasProperty('Molecule')) {
    3964     return undef;
    3965   }
    3966   my($Molecule);
    3967   $Molecule = $This->GetProperty('Molecule');
    3968
    3969   return $Molecule->_GetNumOfAtomRingsWithSizeGreaterThan($This, $RingSize);
    3970 }
    3971
    3972 # Get all rings an array of references to arrays containing ring atoms...
    3973 #
    3974 sub GetRings {
    3975   my($This) = @_;
    3976
    3977   # Is this atom in a molecule?
    3978   if (!$This->HasProperty('Molecule')) {
    3979     return undef;
    3980   }
    3981   my($Molecule);
    3982   $Molecule = $This->GetProperty('Molecule');
    3983
    3984   return $Molecule->_GetAtomRings($This);
    3985 }
    3986
    3987 # Get smallest ring as an array containing ring atoms...
    3988 #
    3989 sub GetSmallestRing {
    3990   my($This) = @_;
    3991
    3992   # Is this atom in a molecule?
    3993   if (!$This->HasProperty('Molecule')) {
    3994     return undef;
    3995   }
    3996   my($Molecule);
    3997   $Molecule = $This->GetProperty('Molecule');
    3998
    3999   return $Molecule->_GetSmallestAtomRing($This);
    4000 }
    4001
    4002 # Get largest ring as an array containing ring atoms...
    4003 #
    4004 sub GetLargestRing {
    4005   my($This) = @_;
    4006
    4007   # Is this atom in a molecule?
    4008   if (!$This->HasProperty('Molecule')) {
    4009     return undef;
    4010   }
    4011   my($Molecule);
    4012   $Molecule = $This->GetProperty('Molecule');
    4013
    4014   return $Molecule->_GetLargestAtomRing($This);
    4015 }
    4016
    4017 # Get odd size rings an array of references to arrays containing ring atoms...
    4018 #
    4019 sub GetRingsWithOddSize {
    4020   my($This) = @_;
    4021
    4022   # Is this atom in a molecule?
    4023   if (!$This->HasProperty('Molecule')) {
    4024     return undef;
    4025   }
    4026   my($Molecule);
    4027   $Molecule = $This->GetProperty('Molecule');
    4028
    4029   return $Molecule->_GetAtomRingsWithOddSize($This);
    4030 }
    4031
    4032 # Get even size rings an array of references to arrays containing ring atoms...
    4033 #
    4034 sub GetRingsWithEvenSize {
    4035   my($This) = @_;
    4036
    4037   # Is this atom in a molecule?
    4038   if (!$This->HasProperty('Molecule')) {
    4039     return undef;
    4040   }
    4041   my($Molecule);
    4042   $Molecule = $This->GetProperty('Molecule');
    4043
    4044   return $Molecule->_GetAtomRingsWithEvenSize($This);
    4045 }
    4046
    4047 # Get rings with specified size as an array of references to arrays containing ring atoms...
    4048 #
    4049 sub GetRingsWithSize {
    4050   my($This, $RingSize) = @_;
    4051
    4052   # Is this atom in a molecule?
    4053   if (!$This->HasProperty('Molecule')) {
    4054     return undef;
    4055   }
    4056   my($Molecule);
    4057   $Molecule = $This->GetProperty('Molecule');
    4058
    4059   return $Molecule->_GetAtomRingsWithSize($This, $RingSize);
    4060 }
    4061
    4062 # Get rings with size less than specfied size as an array of references to arrays containing ring atoms...
    4063 #
    4064 sub GetRingsWithSizeLessThan {
    4065   my($This, $RingSize) = @_;
    4066
    4067   # Is this atom in a molecule?
    4068   if (!$This->HasProperty('Molecule')) {
    4069     return undef;
    4070   }
    4071   my($Molecule);
    4072   $Molecule = $This->GetProperty('Molecule');
    4073
    4074   return $Molecule->_GetAtomRingsWithSizeLessThan($This, $RingSize);
    4075 }
    4076
    4077 # Get rings with size greater than specfied size as an array of references to arrays containing ring atoms...
    4078 #
    4079 sub GetRingsWithSizeGreaterThan {
    4080   my($This, $RingSize) = @_;
    4081
    4082   # Is this atom in a molecule?
    4083   if (!$This->HasProperty('Molecule')) {
    4084     return undef;
    4085   }
    4086   my($Molecule);
    4087   $Molecule = $This->GetProperty('Molecule');
    4088
    4089   return $Molecule->_GetAtomRingsWithSizeGreaterThan($This, $RingSize);
    4090 }
    4091
    4092 # Get next object ID...
    4093 sub _GetNewObjectID {
    4094   $ObjectID++;
    4095   return $ObjectID;
    4096 }
    4097
    4098 # Return a string containing vertices, edges and other properties...
    4099 sub StringifyAtom {
    4100   my($This) = @_;
    4101   my($AtomString, $ID, $Name, $AtomSymbol, $AtomicNumber, $XYZVector, $AtomicWeight, $ExactMass, $NumOfNeighbors, $NumOfBonds, $Valence, $MissingHydrogens, $TotalHydrogens, $ImplicitHydrogens, $ExplicitHydrogens, $FormalCharge, $Charge, $SpinMultiplicity, $FreeRadicalElectrons, $StereoCenter, $StereoCenterStatus, $StereoChemistry, $StereochemistryString, $RingAtom, $NumOfRings, $AromaticAtom);
    4102
    4103   $ID = $This->GetID();
    4104   $Name = $This->GetName();
    4105   $AtomSymbol = $This->GetAtomSymbol();
    4106   $AtomicNumber = $This->GetAtomicNumber();
    4107   $XYZVector = $This->GetXYZVector();
    4108
    4109   $AtomicWeight = $This->GetAtomicWeight();
    4110   if (!defined $AtomicWeight) {
    4111     $AtomicWeight = 'undefined';
    4112   }
    4113   $ExactMass = $This->GetExactMass();
    4114   if (!defined $ExactMass) {
    4115     $ExactMass = 'undefined';
    4116   }
    4117   $NumOfNeighbors = $This->GetNumOfNeighbors();
    4118   if (!defined $NumOfNeighbors) {
    4119     $NumOfNeighbors = 'undefined';
    4120   }
    4121   $NumOfBonds = $This->GetNumOfBonds();
    4122   if (!defined $NumOfBonds) {
    4123     $NumOfBonds = 'undefined';
    4124   }
    4125   $Valence = $This->GetValence();
    4126   if (!defined $Valence) {
    4127     $Valence = 'undefined';
    4128   }
    4129
    4130   $MissingHydrogens = $This->GetNumOfMissingHydrogens();
    4131   if (!defined $MissingHydrogens) {
    4132     $MissingHydrogens = 'undefined';
    4133   }
    4134   $TotalHydrogens = $This->GetNumOfHydrogens();
    4135   if (!defined $TotalHydrogens) {
    4136     $TotalHydrogens = 'undefined';
    4137   }
    4138   $ImplicitHydrogens = $This->GetNumOfImplicitHydrogens();
    4139   if (!defined $ImplicitHydrogens) {
    4140     $ImplicitHydrogens = 'undefined';
    4141   }
    4142   $ExplicitHydrogens = $This->GetNumOfExplicitHydrogens();
    4143   if (!defined $ExplicitHydrogens) {
    4144     $ExplicitHydrogens = 'undefined';
    4145   }
    4146
    4147   $FormalCharge = $This->GetFormalCharge();
    4148   if (!defined $FormalCharge) {
    4149     $FormalCharge = 'undefined';
    4150   }
    4151   $Charge = $This->GetCharge();
    4152   if (!defined $Charge) {
    4153     $Charge = 'undefined';
    4154   }
    4155
    4156   $SpinMultiplicity = $This->GetSpinMultiplicity();
    4157   if (!defined $SpinMultiplicity) {
    4158     $SpinMultiplicity = 'undefined';
    4159   }
    4160   $FreeRadicalElectrons = $This->GetFreeRadicalElectrons();
    4161   if (!defined $FreeRadicalElectrons) {
    4162     $FreeRadicalElectrons = 'undefined';
    4163   }
    4164
    4165   $RingAtom = $This->IsInRing();
    4166   if (defined $RingAtom) {
    4167     $RingAtom = $RingAtom  ? 'Yes' : 'No';
    4168     $NumOfRings = $This->GetNumOfRings();
    4169   }
    4170   else {
    4171     $RingAtom = 'undefined';
    4172     $NumOfRings = 'undefined';
    4173   }
    4174
    4175   $AromaticAtom = $This->GetAromatic();
    4176   if (defined $AromaticAtom) {
    4177     $AromaticAtom = $AromaticAtom  ? 'Yes' : 'No';
    4178   }
    4179   else {
    4180     $AromaticAtom = 'undefined';
    4181   }
    4182
    4183   $StereochemistryString = '';
    4184   $StereoCenter = $This->GetStereoCenter();
    4185   if (defined $StereoCenter) {
    4186     $StereoCenterStatus = $This->IsStereoCenter() ? 'Yes' : 'No';
    4187     $StereoChemistry = $This->GetStereochemistry();
    4188     if (!defined $StereoChemistry) {
    4189       $StereoChemistry = 'undefined';
    4190     }
    4191     $StereochemistryString = "StereoCenter: $StereoCenterStatus; Stereochemistry: $StereoChemistry";
    4192   }
    4193
    4194   $AtomString = "Atom: ID: $ID; Name: \"$Name\"; AtomSymbol: \"$AtomSymbol\"; AtomicNumber: $AtomicNumber; XYZ: $XYZVector; AtomicWeight: $AtomicWeight; ExactMass: $ExactMass; NumOfNeighbors: $NumOfNeighbors;  NumOfBonds: $NumOfBonds; Valence: $Valence; MissingHydrogens: $MissingHydrogens; TotalHydrogens: $TotalHydrogens; ImplicitHydrogens: $ImplicitHydrogens; ExplicitHydrogens: $ExplicitHydrogens; FormalCharge: $FormalCharge; Charge: $Charge; SpinMultiplicity: $SpinMultiplicity; FreeRadicalElectrons: $FreeRadicalElectrons; RingAtom: $RingAtom; NumOfAtomRings: $NumOfRings; AromaticAtom: $AromaticAtom";
    4195
    4196   if ($StereochemistryString) {
    4197     $AtomString .= "; $StereochemistryString";
    4198   }
    4199
    4200   return $AtomString;
    4201 }
    4202
    4203 # Load appropriate atom data files from <MayaChemTools>/lib directory used by various
    4204 # object methods in the current class...
    4205 #
    4206 sub _LoadAtomClassData {
    4207   my($MayaChemToolsLibDir);
    4208
    4209   $MayaChemToolsLibDir = GetMayaChemToolsLibDirName();
    4210
    4211   _LoadAtomValenceModelData($MayaChemToolsLibDir);
    4212
    4213 }
    4214
    4215 #
    4216 # Load data for supported valence models...
    4217 #
    4218 sub _LoadAtomValenceModelData {
    4219   my($MayaChemToolsLibDir) = @_;
    4220   my($MDLValenceModelDataFile, $DaylightValenceModelDataFile);
    4221
    4222   %MDLValenceModelDataMap = ();
    4223   %DaylightValenceModelDataMap = ();
    4224
    4225   $MDLValenceModelDataFile = $MayaChemToolsLibDir . "/data/MDLValenceModelData.csv";
    4226   $DaylightValenceModelDataFile = $MayaChemToolsLibDir . "/data/DaylightValenceModelData.csv";
    4227
    4228   if (! -e "$MDLValenceModelDataFile") {
    4229     croak "Error: ${ClassName}::_LoadAtomValenceModelData: MayaChemTools package file, $MDLValenceModelDataFile, is missing: Possible installation problems...";
    4230   }
    4231
    4232   if (! -e "$DaylightValenceModelDataFile") {
    4233     croak "Error: ${ClassName}::_LoadAtomValenceModelData: MayaChemTools package file, $DaylightValenceModelDataFile, is missing: Possible installation problems...";
    4234   }
    4235
    4236   _LoadValenceModelDataFile($MDLValenceModelDataFile, \%MDLValenceModelDataMap);
    4237   _LoadValenceModelDataFile($DaylightValenceModelDataFile, \%DaylightValenceModelDataMap);
    4238
    4239 }
    4240
    4241 #
    4242 # Load valence model data file...
    4243 #
    4244 sub _LoadValenceModelDataFile {
    4245   my($DataFile, $DataMapRef) = @_;
    4246
    4247   # File format:
    4248   #
    4249   # "AtomicNumber","ElementSymbol","FormalCharge","CommomValences"
    4250   #
    4251   my($InDelim, $Line, $NumOfCols, @ColLabels, @LineWords);
    4252
    4253   $InDelim = "\,";
    4254   open DATAFILE, "$DataFile" or croak "Couldn't open $DataFile: $! ...";
    4255
    4256   # Skip lines up to column labels...
    4257   LINE: while ($Line = GetTextLine(\*DATAFILE)) {
    4258     if ($Line !~ /^#/) {
    4259       last LINE;
    4260     }
    4261   }
    4262   @ColLabels= quotewords($InDelim, 0, $Line);
    4263   $NumOfCols = @ColLabels;
    4264
    4265   my($AtomicNumber, $ElementSymbol, $FormalCharge, $CommonValences);
    4266
    4267   # Process element data...
    4268   LINE: while ($Line = GetTextLine(\*DATAFILE)) {
    4269     if ($Line =~ /^#/) {
    4270       next LINE;
    4271     }
    4272     @LineWords = ();
    4273     @LineWords = quotewords($InDelim, 0, $Line);
    4274     if (@LineWords != $NumOfCols) {
    4275       croak "Error: ${ClassName}::_LoadValenceModelDataFile: The number of data fields, @LineWords, in $DataFile must be $NumOfCols.\nLine: $Line...";
    4276     }
    4277
    4278     ($AtomicNumber, $ElementSymbol, $FormalCharge, $CommonValences) = @LineWords;
    4279
    4280     if (exists $DataMapRef->{$AtomicNumber}) {
    4281       # Additional data for an element...
    4282       if (exists $DataMapRef->{$AtomicNumber}{$FormalCharge}) {
    4283         # Duplicate data entries for an element...
    4284         carp "Warning: ${ClassName}::_LoadValenceModelDataFile: Ignoring valence data for element with atomic number $AtomicNumber and formal charge $FormalCharge in data file $DataFile: It has already been loaded.\nLine: $Line...";
    4285         next LINE;
    4286       }
    4287     }
    4288     else {
    4289       # Data for a new element...
    4290       %{$DataMapRef->{$AtomicNumber}} = ();
    4291     }
    4292
    4293     %{$DataMapRef->{$AtomicNumber}{$FormalCharge}} = ();
    4294     $DataMapRef->{$AtomicNumber}{$FormalCharge}{ElementSymbol} = $ElementSymbol;
    4295
    4296     @{$DataMapRef->{$AtomicNumber}{$FormalCharge}{CommonValences}} = ();
    4297     @{$DataMapRef->{$AtomicNumber}{$FormalCharge}{CommonValences}} = sort { $a <=> $b } split /\,/, $CommonValences;
    4298   }
    4299   close DATAFILE;
    4300 }
    4301