Crystal.h

/*  ObjCryst++ Object-Oriented Crystallographic Library
    (c) 2000-2002 Vincent Favre-Nicolin vincefn@users.sourceforge.net
        2000-2001 University of Geneva (Switzerland)

    This program is free software; you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation; either version 2 of the License, or
    (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program; if not, write to the Free Software
    Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
*/
/*   Crystal.h header file for the Crystal object
*
*/
#ifndef _OBJCRYST_CRYSTAL_H_
#define _OBJCRYST_CRYSTAL_H_

#include "ObjCryst/CrystVector/CrystVector.h"

#include "ObjCryst/ObjCryst/General.h"
#include "ObjCryst/RefinableObj/RefinableObj.h"
#include "ObjCryst/ObjCryst/UnitCell.h"
#include "ObjCryst/ObjCryst/ScatteringPower.h"
#include "ObjCryst/ObjCryst/Scatterer.h"

//#include <stdlib.h>
#include <string>
//#include <iomanip>
//#include <cmath>
//#include <typeinfo>
//#include <fstream>
//#include <ctime>

namespace ObjCryst
{
class Scatterer; //forward declaration of another header's class :KLUDGE:
extern const RefParType *gpRefParTypeCrystal;
class NiftyStaticGlobalObjectsInitializer_Crystal
{
   public:
      NiftyStaticGlobalObjectsInitializer_Crystal()
      {
         if (mCount++ == 0)
         {
            gpRefParTypeCrystal=new RefParType (gpRefParTypeUnitCell,"Crystal");
         }
      }
      ~NiftyStaticGlobalObjectsInitializer_Crystal()
      {
         if (--mCount == 0)
         {
            delete gpRefParTypeCrystal;
            gpRefParTypeCrystal=0;
         }
      }
   private:
      static long mCount;
};
static NiftyStaticGlobalObjectsInitializer_Crystal NiftyStaticGlobalObjectsInitializer_Crystal_counter;

//######################################################################
//######################################################################
class Crystal:public UnitCell
{
   public:
      Crystal();
      Crystal(const REAL a, const REAL b, const REAL c,
              const string &SpaceGroupId);
      Crystal(const REAL a, const REAL b, const REAL c, const REAL alpha,
              const REAL beta, const REAL gamma,const string &SpaceGroupId);

      Crystal(const Crystal &oldCryst);
      ~Crystal();
      virtual const string& GetClassName() const;

      void AddScatterer(Scatterer *scatt);
      void RemoveScatterer(Scatterer *scatt, const bool del=true);

      long GetNbScatterer()const;
      Scatterer & GetScatt(const string &scattName);
      const Scatterer & GetScatt(const string &scattName) const;
      Scatterer & GetScatt(const long scattIndex);
      const Scatterer & GetScatt(const long scattIndex) const;

      ObjRegistry<Scatterer>& GetScattererRegistry();
      const ObjRegistry<Scatterer>& GetScattererRegistry()const;
      ObjRegistry<ScatteringPower>& GetScatteringPowerRegistry();
      const ObjRegistry<ScatteringPower>& GetScatteringPowerRegistry()const;

      void AddScatteringPower(ScatteringPower *scattPow);
      void RemoveScatteringPower(ScatteringPower *scattPow, const bool del=true);
      ScatteringPower& GetScatteringPower(const string &name);
      const ScatteringPower& GetScatteringPower(const string &name)const;
      const RefinableObjClock& GetMasterClockScatteringPower()const;

      virtual const ScatteringComponentList& GetScatteringComponentList()const;
      const RefinableObjClock& GetClockScattCompList()const;
      void Print(ostream &os=cout) const;

      CrystMatrix_REAL GetMinDistanceTable(const REAL minDistance=0.1) const;
      void PrintMinDistanceTable(const REAL minDistance=0.1,ostream &os=cout) const;

      ostream& POVRayDescription(ostream &os,const CrystalPOVRayOptions &options)const;

      virtual void GLInitDisplayList(const bool onlyIndependentAtoms=false,
                                     const REAL xMin=-.1,const REAL xMax=1.1,
                                     const REAL yMin=-.1,const REAL yMax=1.1,
                                     const REAL zMin=-.1,const REAL zMax=1.1,
                                     const bool displayNames=false,
                                     const bool hideHydrogens=false)const;

      void CalcDynPopCorr(const REAL overlapDist=1., const REAL mergeDist=.0)const ;
      void ResetDynPopCorr()const ;
      REAL GetDynPopCorr(const Scatterer* pscatt, unsigned int component)const;
      void SetUseDynPopCorr(const int use);
      int GetUseDynPopCorr() const;
      REAL GetBumpMergeCost() const;
      void SetBumpMergeDistance(const ScatteringPower &scatt1,
                                const ScatteringPower &scatt2, const REAL dist=1.5);
      void SetBumpMergeDistance(const ScatteringPower &scatt1,
                                const ScatteringPower &scatt2, const REAL dist,
                                const bool allowMerge);

      void RemoveBumpMergeDistance(const ScatteringPower &scatt1,
                                   const ScatteringPower &scatt2);

      struct BumpMergePar
      {
         BumpMergePar();
         BumpMergePar(const REAL dist, const bool canOverlap=false);
         REAL mDist2;
         bool mCanOverlap;
      };

      typedef std::map<pair<const ScatteringPower*, const ScatteringPower*>,Crystal::BumpMergePar > VBumpMergePar;
      const VBumpMergePar& GetBumpMergeParList()const;
      VBumpMergePar& GetBumpMergeParList();
      const RefinableObjClock& GetClockScattererList()const;

      virtual void XMLOutput(ostream &os,int indent=0)const;
      virtual void XMLInput(istream &is,const XMLCrystTag &tag);
      //virtual void XMLInputOld(istream &is,const IOCrystTag &tag);

      virtual void GlobalOptRandomMove(const REAL mutationAmplitude,
                                       const RefParType *type=gpRefParTypeObjCryst);
      virtual REAL GetLogLikelihood()const;
      virtual void CIFOutput(ostream &os, double mindist = 0.5)const;

      virtual void GetGeneGroup(const RefinableObj &obj,
                                CrystVector_uint & groupIndex,
                                unsigned int &firstGroup) const;
      virtual void BeginOptimization(const bool allowApproximations=false,
                                     const bool enableRestraints=false);
      void AddBondValenceRo(const ScatteringPower&,const ScatteringPower&,const REAL ro);
      void RemoveBondValenceRo(const ScatteringPower&,const ScatteringPower&);
      REAL GetBondValenceCost() const;
      std::map<pair<const ScatteringPower*,const ScatteringPower*>, REAL>& GetBondValenceRoList();
      const std::map<pair<const ScatteringPower*,const ScatteringPower*>, REAL>& GetBondValenceRoList()const;
      void Init(const REAL a, const REAL b, const REAL c, const REAL alpha,
                const REAL beta, const REAL gamma,const string &SpaceGroupId,
                const string& name);

      void SetDeleteSubObjInDestructor(const bool b);
      void ConnectAtoms(const REAL min_relat_dist=0.4, const REAL max_relat_dist=1.3, const bool warnuser_fail=false);
   private:
      void InitOptions();

      int FindScatterer(const string &scattName)const;

      void CalcDistTable(const bool fast)const;

      void CalcBondValenceSum()const;

      ObjRegistry<Scatterer> mScattererRegistry ;

      VBumpMergePar mvBumpMergePar;
      RefinableObjClock mBumpMergeParClock;
      mutable RefinableObjClock mBumpMergeCostClock;
      mutable REAL mBumpMergeCost;
      REAL mBumpMergeScale;


      struct Neighbour
      {
         Neighbour(const unsigned long neighbourIndex,const int sym,
                   const REAL dist2);
         unsigned long mNeighbourIndex;
         unsigned int mNeighbourSymmetryIndex;
         REAL mDist2;
      };
      struct NeighbourHood
      {
         unsigned long mIndex;
         unsigned int mUniquePosSymmetryIndex;
         std::vector<Crystal::Neighbour> mvNeighbour;
      };
      mutable std::vector<NeighbourHood> mvDistTableSq;
      mutable RefinableObjClock mDistTableClock;
      mutable REAL mDistTableMaxDistance;

      mutable ScatteringComponentList mScattCompList;

      RefinableObjClock mLatticeClock;
      RefObjOpt mUseDynPopCorr;

      ObjRegistry<ScatteringPower> mScatteringPowerRegistry;

      //Clocks
      RefinableObjClock mClockScattererList;
      mutable RefinableObjClock mClockScattCompList;
      mutable RefinableObjClock mClockNeighborTable;
      mutable RefinableObjClock mClockDynPopCorr;
      RefinableObjClock mMasterClockScatteringPower;

      RefObjOpt mDisplayEnantiomer;

      map<pair<const ScatteringPower*,const ScatteringPower*>, REAL> mvBondValenceRo;
      RefinableObjClock mBondValenceParClock;
      mutable RefinableObjClock mBondValenceCalcClock;
      mutable RefinableObjClock mBondValenceCostClock;
      mutable REAL mBondValenceCost;
      REAL mBondValenceCostScale;
      mutable std::map<long, REAL> mvBondValenceCalc;
      // Flag indicating whether to delete Scatterers and ScatteringPowers in
      // the destructor (default true). Modified by
      // SetDeleteSubObjInDestructor.
      bool mDeleteSubObjInDestructor;

   #ifdef __WX__CRYST__
   public:
      virtual WXCrystObjBasic* WXCreate(wxWindow*);
      friend class WXCrystal;
   #endif
};

extern ObjRegistry<Crystal> gCrystalRegistry;


}// namespace


#endif //_OBJCRYST_CRYSTAL_H_