Pointless: hkl_datatypes.hh Source File

00001 // hkl_datatypes.hh
00002 //
00003 // Data type definitions
00004 // Just reflection datatypes
00005 
00006 
00007 #ifndef HKL_DATATYPES_HEADER
00008 #define HKL_DATATYPES_HEADER
00009 
00010 // Clipper
00011 #include <clipper/clipper.h>
00012 #include "clipper/core/clipper_precision.h"
00013 using clipper::ftype;
00014 using clipper::Vec3;
00015 using clipper::Mat33;
00016 using clipper::String;
00017 using clipper::Metric_tensor;
00018 typedef clipper::Vec3<double> DVect3;
00019 typedef clipper::Mat33<double> DMat33;
00020 typedef clipper::Vec3<float> FVect3;
00021 typedef clipper::Mat33<float> FMat33;
00022 typedef clipper::Vec3<int> IVect3;
00023 
00024 #include "cmtzlib.h"    // CCP4 MTZlib headers (namespace CMtz)
00025 #include "csymlib.h"    // CCP4 symmetry stuff
00026 #include "matvec_utils.hh"  // Matrix & vector utilities
00027 #include "util.hh"
00028 
00029 typedef float  Rtype;
00030 typedef double Dtype;
00031 
00032 typedef std::pair<int,int> IndexPair;
00033 typedef std::pair<float,float> RPair;
00034 typedef std::pair<double,double> DPair;
00035 
00036 
00037 namespace scala
00038 {
00039   enum Chirality {UNKNOWN, CHIRAL, NONCHIRAL, CENTROSYMMETRIC};
00040   std::string Chiral_as_string(const Chirality& chiral);
00041 
00042   class SpaceGroup; // forward definition
00043 
00044   enum AnomalousClass {ALL, BOTH, IPLUS, IMINUS};
00045 
00046   //--------------------------------------------------------------
00047   // Orientation things
00048   //
00049   // Returns angles in degrees between principle reciprocal axes a*, b*, c*
00050   // and the rotation spindle (along z)
00051   // On entry:
00052   //  DUB is setting matrix [D][U][B]
00053   DVect3 SpindleToPrincipleAxes(const DMat33& DUB);
00054   //--------------------------------------------------------------
00055   // Returns index of smallest component of vector
00056   int SmallestComponent(const DVect3& v);
00057   //--------------------------------------------------------------
00058   // Format closet reciprocal axis to spindle
00059   std::string FormatSpindleToPrincipleAxes(const DMat33& DUB);
00060   //--------------------------------------------------------------
00061   class IntRange;  // forward definition
00062   //================================================================
00063   typedef clipper::datatypes::I_sigI<float> ClipperI_sigI;
00064   class IsigI
00066   {
00067   public:
00068     IsigI() {clipper::Util::set_null(I_); clipper::Util::set_null(sigI_); }
00069     IsigI(const float& I, const float& sigI) {I_=I; sigI_=sigI;}
00070     IsigI(const IsigI& Is) {
00071       if (this != &Is) {I_=Is.I(); sigI_=Is.sigI();}
00072     }
00073     // Construct from clipper I_sigI
00074     IsigI(const ClipperI_sigI& cIsig) {I_=cIsig.I(); sigI_=cIsig.sigI();}
00075 
00076     // Accessors
00077     const float&    I() const {return  I_;}
00078     const float& sigI() const {return sigI_;}
00079     // write access
00080     float&    I() {return    I_;}
00081     float& sigI() {return sigI_;}
00082 
00083     void scale(const float& a) {I_ *= a; sigI_ *=a;}
00084     void scale(const double& a) {I_ *= float(a); sigI_ *=float(a);}
00085 
00086     IsigI scaleIs(const float& a) const {return IsigI(I_*a, sigI_*a);}
00087     IsigI scaleIs(const double& a) const {return IsigI(I_*a, sigI_*a);}
00088 
00089     bool missing() const { return (clipper::Util::is_nan(I_) || clipper::Util::is_nan(sigI_)); }
00090     bool missingI() const { return (clipper::Util::is_nan(I_)); }
00091     bool missingsigI() const { return (clipper::Util::is_nan(sigI_)); }
00092 
00093 
00094   private:
00095     float I_, sigI_;
00096   };
00097 
00098   //======================================================================
00099   //--------------------------------------------------------------
00100   class ReindexOp : public RTop<double>
00102 
00104   {
00105   public:
00106     ReindexOp() : RTop<double>(RTop<double>::identity()),
00107       strict(false), deviation(0.0) {}
00109     ReindexOp(const Mat33<double>& Hin,
00110               const Vec3<double>& Vin = Vec3<double>(0.0,0.0,0.0))
00111       : RTop<double>(Hin, Vin),
00112       strict(false), deviation(0.0) {}
00114     ReindexOp(const RTop<double>& RTin)
00115       : RTop<double>(RTin),
00116       strict(false), deviation(0.0) {}
00118     ReindexOp(const ReindexOp& Rdxin)
00119       : RTop<double>(Rdxin)    {
00120       strict = Rdxin.strict;
00121       deviation = Rdxin.deviation;
00122     }
00123 
00125     ReindexOp(const std::string& Operator);
00126 
00128     void SetStrict(const bool& Strict) {strict = Strict;}
00130     bool Strict() const {return strict;}
00132     void SetDeviation(const double& Deviation) {deviation = Deviation;}
00134     double Deviation() const {return deviation;}
00135     
00137     Mat33<double> transpose() const {return rot().transpose();}
00139     ReindexOp inverse() const {return ReindexOp(RTop<double>::inverse());}
00140     
00142     std::vector<double> as_Dvec() const {return MVutil::SetVMat33(rot());}
00143     
00145     bool IsIdentity() const {return MVutil::is_rtop_ident(*this);}
00147     bool Positive() const {
00148       return (rot().det() > 0.0);}
00149     
00151     std::string as_hkl() const;
00153     std::string as_matrix() const;
00155     std::string as_XML() const;
00157     std::string as_hkl_XML() const;
00159     bool equals(const ReindexOp& b,
00160                 const double& tol = 1.0e-6) const;
00161     
00163     void FindSimplest(const SpaceGroup& symm);
00164 
00166     friend bool operator < (const ReindexOp& a,const ReindexOp& b);
00167     
00168     friend bool operator == (const ReindexOp& a,const ReindexOp& b);
00169     friend bool operator != (const ReindexOp& a,const ReindexOp& b);
00170     
00171     friend ReindexOp operator * (const ReindexOp& a,const ReindexOp& b);
00172     friend Vec3<double> operator * (const Vec3<double> v, const ReindexOp& R);
00173     friend ReindexOp operator * (const RTop<double>& a,const ReindexOp& b);
00174     friend ReindexOp operator * (const ReindexOp& a,const RTop<double>& b);
00175     
00176      
00177   private:
00178     bool strict;
00179     double deviation;
00180   };
00181   //--------------------------------------------------------------
00183 
00192   class MetricTensor
00193   {
00194   public:
00196     inline MetricTensor() {} 
00197 
00198     //  angles in degrees
00199     MetricTensor( const ftype& a, const ftype& b, const ftype& c, const ftype& alph, const ftype& beta, const ftype& gamm );
00201     MetricTensor(const std::vector<Dtype> cell);
00203     MetricTensor(const Mat33<Dtype> MetricMatrix);
00205     inline ftype lengthsq( const Vec3<>& v ) const
00206     { return ( v[0]*(v[0]*m00 + v[1]*m01 + v[2]*m02) +
00207                v[1]*(v[1]*m11 + v[2]*m12) + v[2]*(v[2]*m22) ); }
00209     inline ftype lengthsq( const Vec3<int>& v ) const
00210     {   ftype h = ftype(v[0]); ftype k = ftype(v[1]); ftype l = ftype(v[2]);
00211       return h*(h*m00 + k*m01 + l*m02) + k*(k*m11 + l*m12) + l*(l*m22); }
00212 
00214     std::vector<Dtype> cell() const;
00215 
00217     MetricTensor inverse() const;
00218 
00220     Mat33<Dtype> matrix() const;
00221 
00222     String format() const;  
00223   private:
00224     ftype m00, m11, m22, m01, m02, m12;
00225   };
00226 
00227   //======================================================================
00228   class Scell
00230 
00234   {
00235   public:
00236     Scell() {cell_.assign(6,0.0);}
00238     //  derive reciprocal cell and orthogonalisation matrix B
00239     Scell(const std::vector<Dtype>& rcell);
00241     Scell(const float* cell);
00243     Scell(const MetricTensor& metric_tensor, const bool real=true);
00245     Scell(const clipper::Cell& ccell);
00247     Scell(const double& a,const double& b,const double& c,
00248           const double& alpha, const double& beta, const double& gamma);
00250     void init(const std::vector<Dtype>& rcell);
00251 
00252     std::string formatPrint(const bool& newline=true) const;  
00253     void dump() const;   
00254 
00256     inline Mat33<Dtype> Bmat() const {return Bmat_;}
00258     inline MetricTensor rec_metric_tensor() const
00259     {return recip_metric_tensor_;}
00260 
00262     inline std::vector<Dtype>  UnitCell() const {return cell_;}
00264     inline std::vector<Dtype>  ReciprocalCell() const {return rec_cell_;}
00266     inline const Dtype& operator[](const int& i) const {return cell_[i];}
00268 
00273     Scell change_basis(const ReindexOp& reindex_op) const;
00274 
00276     std::string format(const int w=7, const int p=4) const;
00278     std::string xml() const;
00280     clipper::Cell ClipperCell() const
00281     {return clipper::Cell(clipper::Cell_descr
00282           (cell_[0],cell_[1],cell_[2],cell_[3],cell_[4],cell_[5]));}
00284     double Volume() const;
00285 
00287     bool AngleTest(const double& AlphaTest,
00288                    const double& BetaTest,
00289                    const double& GammaTest,
00290                    const double Tol=0.2) const;
00292     bool equals(const Scell& other, const double& tolA=0.02, const double& tolD=0.02) const;
00294     bool equalsTol(const Scell& other, const double& AngularTolerance) const;
00295 
00297     double Difference(const Scell& other) const;
00298 
00300     bool null() const;
00301 
00302   private:
00303     std::vector<Dtype> cell_;
00304     std::vector<Dtype> rec_cell_;
00305     Mat33<Dtype> Bmat_;
00306     MetricTensor  metric_tensor_;
00307     MetricTensor  recip_metric_tensor_;
00308   };
00309   //======================================================================
00310   class UnitCellSet {
00312   public:
00313     UnitCellSet(){}
00315     UnitCellSet(const std::vector<Scell>& Cells);
00317     void init(const std::vector<Scell>& Cells);
00318 
00320     void AddCell(const Scell& Cell);
00321 
00323     int Number() const {return cells.size();}
00324 
00326     std::vector<Scell> Cells() const {return cells;}
00327 
00328     // return average cell
00329     Scell AverageCell() const {return averagecell;}
00330 
00332     std::vector<double> Deviations() const;
00334     double WorstDeviation() const;
00336     std::vector<double> RmsD() const;
00337 
00339     Scell Average(const int& idxexclude=-1) const;
00340 
00341   private:
00342     std::vector<Scell> cells;
00343     Scell averagecell;
00344 
00345   }; // UnitCellSet
00346   //======================================================================
00348 
00349   class Hkl : public Vec3<int>
00350   {
00351   public:
00352     inline Hkl() {}                
00353     inline explicit Hkl( const Vec3<int>& v ) :
00354       Vec3<int>( v ) {}            
00355     inline Hkl( const int& h, const int& k, const int& l ) :
00356       Vec3<int>( h, k, l ) {}      
00357     inline Hkl( const Vec3<double>& d ) :
00358       Vec3<int>(Nint(d[0]), Nint(d[1]), Nint(d[2])) {} 
00359     inline Hkl( const clipper::HKL& chkl) :
00360       Vec3<int>(chkl.h(),chkl.k(),chkl.l()) {} 
00361     inline const int& h() const { return (*this)[0]; }  
00362     inline const int& k() const { return (*this)[1]; }  
00363     inline const int& l() const { return (*this)[2]; }  
00364     inline int& h() { return (*this)[0]; }  
00365     inline int& k() { return (*this)[1]; }  
00366     inline int& l() { return (*this)[2]; }  
00367     bool IsGeneral() const; 
00368 
00369     clipper::HKL HKL() const {return clipper::HKL(*this);}
00370 
00372     inline Vec3<Dtype> real() const
00373     {return Vec3<Dtype>(Dtype((*this)[0]),Dtype((*this)[1]),Dtype((*this)[2]));}
00375     inline Dtype invresolsq( const Scell& cell ) const
00376     {return cell.rec_metric_tensor().lengthsq(real());}
00378     Hkl change_basis(const ReindexOp& reindex_op) const;
00380     bool change_basis(Hkl& Newhkl, const ReindexOp& reindex_op) const;
00382     Vec3<Dtype> orth(const Scell& cell) const
00383     {return (cell.Bmat() * real());}
00384     std::string format() const;  
00385 
00386     //  Note that this uses 10-bit packing, so is not guaranteed to
00387     //  produce a unique code, but it is good enough for some purposes
00388     int code() const;
00389 
00390     friend inline Hkl operator -(const Hkl& h1)
00391     { return Hkl( -h1.h(), -h1.k(), -h1.l() ); }
00392     friend inline Hkl operator +(const Hkl& h1, const Hkl& h2)
00393     { return Hkl( h1.h()+h2.h(), h1.k()+h2.k(), h1.l()+h2.l() ); }
00394     friend inline Hkl operator -(const Hkl& h1, const Hkl& h2)
00395     { return Hkl( h1.h()-h2.h(), h1.k()-h2.k(), h1.l()-h2.l() ); }
00396     friend inline Hkl operator *(const int& s, const Hkl& h1)
00397     { return Hkl( s*h1.h(), s*h1.k(), s*h1.l() ); }
00398   };
00399   //--------------------------------------------------------------
00400   Hkl HklDecode(const int& kode);
00401 
00402   //======================================================================
00403   class PxdName
00405   {
00406   public:
00407     PxdName() : pname_(""), xname_(""), dname_("") {}
00409     PxdName(const std::string& pname_in,
00410             const std::string& xname_in,
00411             const std::string& dname_in)
00412       :  pname_(pname_in), xname_(xname_in), dname_(dname_in) {}
00413   
00414   
00415     const std::string& pname() const {return pname_;} 
00416     const std::string& xname() const {return xname_;} 
00417     const std::string& dname() const {return dname_;} 
00418     
00419     std::string& pname() {return pname_;} 
00420     std::string& xname() {return xname_;} 
00421     std::string& dname() {return dname_;} 
00422 
00424     void update(const PxdName& NewPxd);
00425   
00426     std::string formatPrint() const; 
00427     std::string format() const; 
00428 
00430     friend inline bool operator ==(const PxdName& pxd1, const PxdName& pxd2)
00431     {
00432       return (pxd1.pname()==pxd2.pname() &&
00433               pxd1.xname()==pxd2.xname() &&
00434               pxd1.dname()==pxd2.dname());
00435     }
00437     bool is_blank() const
00438     {
00439       return (pname_=="" &&
00440               xname_=="" &&
00441               dname_=="");
00442     }
00443 
00444   private:
00445     std::string pname_, xname_, dname_;
00446   };
00447   //======================================================================
00448   class Xdataset
00450 
00452   {
00453   public:
00454     Xdataset(){}
00456     Xdataset(const PxdName& pxdname, const Scell& cell,
00457              const float& wavel, const int& setid);
00458   
00460     void add_batch(const int& batch_num);
00461     void ClearBatchList() {batches.clear();} 
00462     int num_batches() const {return batches.size();} 
00463 
00464     void change_basis(const ReindexOp& reindex_op);
00465 
00466     void AddRunIndex(const int& RunIndex); 
00467     void ClearRunList() {run_index_list.clear();} 
00468     std::vector<int> RunIndexList() const {return run_index_list;} 
00469   
00470     PxdName pxdname() const {return pxdname_;} 
00471     int setid() const {return setid_;} 
00472     int& setid() {return setid_;}  
00473 
00474     Scell cell() const {return cell_;} 
00475     Scell& cell() {return cell_;} 
00476     float wavelength() const {return wavel_;} 
00477     float& wavelength() {return wavel_;} 
00478 
00479     float Mosaicity() const {return av_mosaic;} 
00480     float& Mosaicity() {return av_mosaic;}  
00481    
00482     std::string formatPrint() const; 
00483 
00485     void AddCellWavelength(const Scell& newcell, const float& wavel);
00486 
00488     void AverageCellWavelength(); 
00490     int NumberofCells() const {return allcells_.Number();}
00491 
00492     std::string formatAllCells() const; 
00493 
00494     double WorstDeviation() const;
00495 
00497     friend bool operator == (const Xdataset& a,const Xdataset& b);
00498     
00499   private:
00500     PxdName pxdname_;
00501     int setid_;
00502     Scell cell_;
00503     float wavel_;
00504     float av_mosaic;
00505     // List of all batch numbers belonging to this dataset
00506     std::vector<int> batches;
00507     // List of runs (indices)
00508     std::vector<int> run_index_list;
00509     // List of unit cells if multiple runs
00510     UnitCellSet allcells_;
00511     // List of wavelengths if multiple runs
00512     std::vector<float> allwavel_;
00513   };
00514   //======================================================================
00515   class BatchNumber
00517   {
00518   public:
00519     BatchNumber() : number(0), original_number(0) {}
00521     BatchNumber(const int& n) : number(n), original_number(n) {}
00523     BatchNumber(const int& n, const int& norig):  number(n), original_number(norig) {}
00524 
00525     int Number() const {return number;}  
00526     int OriginalNumber() const {return original_number;} 
00527 
00528   private:
00529     int number;
00530     int original_number;
00531   };
00532   //======================================================================
00533   class Batch
00535 
00536   {
00537   public:
00538     Batch();  // a dummy batch
00540     Batch(const CMtz::MTZBAT& batch, const bool& accept, const int& idataset);
00542     void change_basis(const ReindexOp& reindex_op);
00544     void SetCellConstraint(const std::vector<int>& lbcell);
00546     void SetCell(const Scell& newcell);
00548     int& FileNumber() {return file_num;}
00550     int  FileNumber() const {return file_num;}
00551 
00553     void OffsetPhi(const float& Phioffset);
00554     float PhiOffset() const {return phioffset;} 
00555     void SetRunIndex(const int& RunIndex); 
00556     int RunIndex() const {return run_index;} 
00557 
00558     int num() const {return batchinfo.num;} 
00559     int& num() {return batchinfo.num;}  
00560 
00562     void OffsetNum(const int& Offset);
00564     int  BatchNumberOffset() const {return offset;}  // get
00565 
00566     int index() const {return Xdataset_index;} 
00567     int& index() {return Xdataset_index;} 
00568     int DatasetID() const {return batchinfo.nbsetid;} 
00569     int& DatasetID() {return batchinfo.nbsetid;}  
00570 
00572     float Phi1() const {return valid_phi ? batchinfo.phistt : 0.0f;}
00574     float Phi2() const {return valid_phi ? batchinfo.phiend : 0.0f;}
00576     float MidPhi() const {return valid_phi ?  0.5*(batchinfo.phistt+batchinfo.phiend) : 0.0f;}
00577     bool ValidPhi() const {return valid_phi;} 
00578     float PhiRange() const
00579     {return Phi2() - Phi1();} 
00580 
00581     // valid_time: > 0 valid time information (time1, time2)
00582     //             < 0 time inferred from phi
00583     //              -1 set = phi
00584     //              -2 inverted from descending phi
00585     //             = 0 no time information
00587     float Time1() const {return (valid_time!=0) ? batchinfo.time1 : 0.0f;}
00589     float Time2() const {return (valid_time!=0) ? batchinfo.time2 : 0.0f;}
00591     float MidTime() const {return (valid_time!=0) ? 0.5*(batchinfo.time1+batchinfo.time2) : 0.0f;}
00593     bool IsValidTime() const {return (valid_time!=0);}
00595     bool IsTimePhi() const {return (valid_time<0);}
00597     void SetTimeColumnPresent() {valid_time = +1;}
00599     void SetTimeColumnAbsent() {valid_time = 0;}
00601     void SetPhiAsProxyTime() {valid_time = -1;}
00603     void SetTimeReversedFromPhi() {valid_time = -2;}
00605     void OffsetTime(const float& Timeoffset);
00607     void StoreTimeRange(const float& Time1, const float& Time2) {
00608       batchinfo.time1 =  Time1;
00609       batchinfo.time2 =  Time2;
00610     }
00611 
00612 
00613     bool Accepted() const {return accepted;}  
00614     void SetAccept(const bool& accept) {accepted = accept;} 
00615 
00616     Scell cell() const {return bcell;} 
00617     float Mosaicity() const {return batchinfo.crydat[0];} 
00618     float Wavelength() const {return batchinfo.alambd;} 
00619 
00621     bool ValidOrientation() const {return valid_Umat;}
00622     // Orientation information:
00623     //   x = [R] [D] [U] [B] h
00624     //  where h   = (h k l)
00625     //       [B]  orthogonalisation matrix from unit cell
00626     //       [U]  orientation matrix
00627     // Two cases for rotation:
00628     //   1) rotation is about the outer goniostat axis e1
00629     //      (a) single axis Phi around e1
00630     //          [R] = [Phi]      [D] = [Phi0], often = [I]
00631     //      (b) 3-axis goniostat, spindle rotation by omega around e1
00632     //          [R] = [Omega]    [D] = [Chi/Kappa][Phi]
00633     //   2) 3-axis goniostat, spindle rotation by phi around e3
00634     //          [R] = [Omega][Chi/Kappa][Phi]     [D] = [Phi0]
00635     //              = [E1E2] [Phi]  (phi varying)
00636     //              
00637     // Transformation operations:
00638     //   four useful reference frames
00639     //    (a) reciprocal lattice indices h
00640     //    (b) camera frame (Ewald sphere frame)
00641     //        s(r)  = x = [R][D][U][B]h
00642     //    (c) zero rotation angle frame (phi/omega = 0)
00643     //        s(r0) = [D][U][B]h = [R]^-1 s(r)
00644     //    (d) polar orthogonalised crystal frame to put a defined
00645     //      reciprocal lattice axis ~ along z
00646     //        s(pole) = [P][B]h
00647     //                = [P][DU]^-1 s(r0)
00648     //       [P] is a cyclic permutation matrix (so [P]^-1 = [P]transpose)
00649     //
00650     // NB in the following FtoF rountines, all vectors (apart from hkl) are
00651     // in dimensionless reciprocal lattice units (ir s(r), s(r0), s(pole)
00652 
00654     DVect3 HtoSr(const Hkl& hkl, const float& phi) const;
00656     DVect3 HtoSr0(const Hkl& hkl) const;
00658     DVect3 SrtoSr0(const DVect3& sr, const float& phi) const;
00660 
00661     DVect3 Sr0toP(const DVect3& sr0) const;
00662 
00663     DMat33 Umat() const {return U;} 
00664 
00666     std::string SpindleToPrincipleAxis() const;
00667 
00669     DVect3 Source() const {return s0;}
00671     bool PhiScan() const {return phiscan;}
00672 
00674 
00676     void SetPole(const int& Pole);
00678     int Pole() const {return pole;}
00679 
00681 
00682     CMtz::MTZBAT batchdata() const {return batchinfo;}
00683 
00685     std::string format() const;
00687     void print() const;
00688 
00690     friend bool operator < (const Batch& a,const Batch& b);
00691 
00692   private:
00693     // initialise batchinfo to zeroes
00694     void initBatchInfo();
00695 
00696     CMtz::MTZBAT batchinfo;   // All in original frame
00697     int Xdataset_index;      // index into crystal/dataset list (-1 if rejected)
00698     int run_index;           // index into run list
00699     bool accepted;
00700     bool valid_cell;
00701     bool valid_Umat;
00702     // valid_time: > 0 valid time information (time1, time2)
00703     //             < 0 time inferred from phi
00704     //              -1 set = phi
00705     //              -2 inverted from descending phi
00706     //             = 0 no time information
00707     int valid_time;
00708     bool valid_phi;  // valid phi information (phi1, phi2)
00709     float phioffset;   // Offset in phi from original
00710     float timeoffset;   // Offset in time from original
00711     int offset;      // number offset from original
00712     int file_num;    // serial number of file from which this batch came
00713 
00714     DVect3 spindle; // goniostat vector for principle rotation axis (= e1 or e3)
00715 
00716     // Stored information, in Cambridge frame
00717     //  Cambridge frame has
00718     //    source vector s0 (anti-parallel to beam) approximately along -x
00719     //      ie x is along beam (if perpendicular to rotation axis)
00720     //    principal rotation axis e0 exactly along z
00721     DMat33 Q;  // conversion matrix from input frame to Cambridge
00722     DVect3 s0;  // source vector
00723     DMat33 U;  // [U]     [Orientation]
00724     DMat33 DU; // [D][U]  [Datum][Orientation]
00725     DMat33 DUB; // [D][U][B]  [Datum][Orientation][Orthogonalisation]
00726                    
00727 
00728     int pole;      // absorption pole, 0 if unset
00729     DMat33 PDUinv; // [P][[D][U]]^-1 for crystal frame absorption pole
00730     bool phiscan;  // true for Phi scan in 3-axis system, else false
00731     DMat33 E1E2;   // [[Omega][Chi/Kappa]] matrix for Phi scan in
00732                     // 3-axis system, else [I]
00733 
00734     Scell bcell;
00735   };
00736   //======================================================================
00737   // Return true if dataset setid is in datasets list
00738   //  & return dataset index idataset (-1 if not)
00739   bool in_datasets(const int& setid,
00740                    const std::vector<Xdataset>& datasets,
00741                    int& idataset);
00742   //======================================================================
00743   class BatchSelection
00745 
00751   {
00752   public:
00753     BatchSelection();
00754     void clear(); 
00755 
00757     bool Null() const {return (batchlist.size()==0 && flaglist.size()==0);}
00758 
00760 
00762     void AddBatch(const int& batch, const int& fileSeriesList);
00764 
00766     void AddRange(const int& batch1, const int& batch2,
00767                   const int& fileSeriesList, const int& flag=-1);
00768 
00769 
00771     int FindInSelection(const int& batch,
00772                         const int& fileSeriesTest) const;
00773 
00775     bool InSelection(const int& batch, const int& fileSeriesTest) const;
00776 
00781     //<!
00782     //<!  2) Specified selection on original file numbering from 2nd or
00783     //<!  subsequent file, fileSeriesList > 0, only test if
00784     //<!  fileSeriesTest > 0
00785 
00787     int FlagNumber(const int& batch);
00788 
00790     std::vector<int> UniqueFlags() const;
00791 
00793     std::vector<IntRange> BatchRanges(const int& flag);
00794 
00796     void CheckSeries(const int& NumFileSeries) const;
00797 
00798   private:
00799     std::vector<int> batchlist;         // list of batch number
00800     std::vector<int> fileseries_list;   // file series numbers for list 
00801     std::vector<IntRange> batchranges;  // batch ranges
00802     std::vector<int> fileseries_range;  // file series numbers for ranges
00803     std::vector<int> flaglist;          // a numerical flag eg a run number
00804   };
00805 //--------------------------------------------------------------
00806   class FileRead
00808   {
00809   public:
00810     FileRead() :read(false), nexcluded(0) {}
00812 
00817     FileRead(const bool& Opened, const bool& Read, 
00818              const bool& SymmetrySet, const int& Nexcluded)
00819       : opened(Opened), read(Read),
00820         symmetrySet(SymmetrySet), nexcluded(Nexcluded)   {}
00821     
00823     bool Opened() const {return opened;}
00825     bool Read() const {return read;}
00827     bool SymmetrySet() const {return symmetrySet;}
00829     int Nexcluded() const {return nexcluded;}
00830 
00831   private:
00832     bool opened;        // true if successfully opened
00833     bool read;          // true if file was read
00834     int  nexcluded;     // number of observation(parts) excluded
00835     bool symmetrySet;   // true if the file contained symmetry information
00836                         // (eg lattice type)
00837   };
00838 }  // namespace scala
00839 #endif
00840 
00841 // end  hkl_datatypes.hh