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