1
/* hkl_info.cpp: class file for hkl_info class + children */
2
//C Copyright (C) 2000-2006 Kevin Cowtan and University of York
4
//L This library is free software and is distributed under the terms
5
//L and conditions of version 2.1 of the GNU Lesser General Public
6
//L Licence (LGPL) with the following additional clause:
8
//L `You may also combine or link a "work that uses the Library" to
9
//L produce a work containing portions of the Library, and distribute
10
//L that work under terms of your choice, provided that you give
11
//L prominent notice with each copy of the work that the specified
12
//L version of the Library is used in it, and that you include or
13
//L provide public access to the complete corresponding
14
//L machine-readable source code for the Library including whatever
15
//L changes were used in the work. (i.e. If you make changes to the
16
//L Library you must distribute those, but you do not need to
17
//L distribute source or object code to those portions of the work
18
//L not covered by this licence.)'
20
//L Note that this clause grants an additional right and does not impose
21
//L any additional restriction, and so does not affect compatibility
22
//L with the GNU General Public Licence (GPL). If you wish to negotiate
23
//L other terms, please contact the maintainer.
25
//L You can redistribute it and/or modify the library under the terms of
26
//L the GNU Lesser General Public License as published by the Free Software
27
//L Foundation; either version 2.1 of the License, or (at your option) any
30
//L This library is distributed in the hope that it will be useful, but
31
//L WITHOUT ANY WARRANTY; without even the implied warranty of
32
//L MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
33
//L Lesser General Public License for more details.
35
//L You should have received a copy of the CCP4 licence and/or GNU
36
//L Lesser General Public License along with this library; if not, write
37
//L to the CCP4 Secretary, Daresbury Laboratory, Warrington WA4 4AD, UK.
38
//L The GNU Lesser General Public can also be obtained by writing to the
39
//L Free Software Foundation, Inc., 59 Temple Place, Suite 330, Boston,
49
Message_fatal message_recip_asu_error( "HKL_info: find_sym reciprocal space ASU error" );
50
Message_ctor message_ctor_hkl_info( " [HKL_info: constructed]" );
53
// methods for 'HKL_info'
56
/*! Update all the lookup tables to be consistent with the modified
58
void HKL_info::update_hkl_list() {
60
hkl_class_lookup.resize( num_reflections() );
61
invresolsq_lookup.resize( num_reflections() );
62
invresolsq_range_ = Range<ftype>();
63
for ( int i = 0; i < num_reflections(); i++ ) {
64
hkl_class_lookup[i] = spacegroup_.hkl_class( hkl_of( i ) );
65
invresolsq_lookup[i] = hkl_of( i ).invresolsq( cell_ );
66
invresolsq_range_.include( invresolsq_lookup[i] );
72
// constructors/destructor
76
Message::message( message_ctor_hkl_info );
80
/*! Construct and initialise HKL_info object. This updates the
81
spacegroup and cell and clears the reflection list. The resolution
82
is used as a rejection criterion for reflections - no HKL will be
83
stored beyond the given limit. Initially there are no reflections in
84
the reflection list: see generate_hkl_list().
86
If any of the parameters have null values, the existing values will
87
be unchanged. The object will only be fully initialised once all
88
parameters are available.
89
\param spacegroup The spacegroup.
90
\param cell The unit cell.
91
\param resolution The resolution limit. */
92
HKL_info::HKL_info( const Spacegroup& spacegroup, const Cell& cell, const Resolution& resolution, const bool& generate )
94
init( spacegroup, cell, resolution, generate );
95
Message::message( message_ctor_hkl_info );
99
/*! Initialise the HKL_info object. This updates the spacegroup and
100
cell and clears the reflection list. The resolution is used as a
101
rejection criterion for reflections - no HKL will be stored beyond
102
the given limit. Initially there are no reflections in the
103
reflection list: see generate_hkl_list().
105
If any of the parameters have null values, the existing values will
106
be unchanged. The object will only be fully initialised once all
107
parameters are available.
108
\param spacegroup The spacegroup.
109
\param cell The unit cell.
110
\param resolution The resolution limit.
111
\param generate If true, a reflection list will be generated for an ASU. */
112
void HKL_info::init( const Spacegroup& spacegroup, const Cell& cell, const Resolution& resolution, const bool& generate )
114
// set spacegroup and cell
115
spacegroup_ = spacegroup;
117
resolution_ = resolution;
120
if ( is_null() ) return;
122
// Create the intergised symops (grid irrelevent)
123
Grid g( 24, 24, 24 );
124
isymop.resize( spacegroup_.num_symops() );
125
for ( int sym = 0; sym < spacegroup_.num_symops(); sym++ )
126
isymop[sym] = Isymop( spacegroup_.symop(sym), g );
130
if ( generate ) generate_hkl_list();
135
/*! Initialise the HKL_info object. This updates the spacegroup and
136
cell and clears the reflection list. The HKL_sampling determines
139
If any of the parameters have null values, the existing values will
140
be unchanged. The object will only be fully initialised once all
141
parameters are available.
142
\param spacegroup The spacegroup.
143
\param cell The unit cell.
144
\param hkl_sampling The resolution limit.
145
\param generate If true, a reflection list will be generated for an ASU. */
146
void HKL_info::init( const Spacegroup& spacegroup, const Cell& cell, const HKL_sampling& hkl_sampling, const bool& generate )
148
// set spacegroup and cell
149
spacegroup_ = spacegroup;
151
hkl_sampling_ = hkl_sampling;
154
if ( spacegroup_.is_null() || cell_.is_null() || hkl_sampling_.is_null() ) return;
156
// estimate resolution
157
resolution_ = hkl_sampling.resolution( cell );
159
// Create the intergised symops (grid irrelevent)
160
Grid g( 24, 24, 24 );
161
isymop.resize( spacegroup_.num_symops() );
162
for ( int sym = 0; sym < spacegroup_.num_symops(); sym++ )
163
isymop[sym] = Isymop( spacegroup_.symop(sym), g );
168
std::vector<HKL> add;
169
// make a reflection list for the given refln sampling
170
HKL lim = hkl_sampling.hkl_limit();
172
// make a list of valid reflections
173
for (rfl.h()=-lim.h(); rfl.h()<=lim.h(); rfl.h()++)
174
for (rfl.k()=-lim.k(); rfl.k()<=lim.k(); rfl.k()++)
175
for (rfl.l()=-lim.l(); rfl.l()<=lim.l(); rfl.l()++)
176
if ( spacegroup_.recip_asu(rfl) &&
177
hkl_sampling.in_resolution( rfl ) &&
178
!( spacegroup_.hkl_class(rfl).sys_abs() ) ) add.push_back(rfl);
179
// update the reflection data lists
187
/*! \return true if the object has not been initalised. */
188
bool HKL_info::is_null() const
189
{ return ( spacegroup_.is_null() || cell_.is_null() || resolution_.is_null() ); }
192
/*! Using current cell, spacegroup, resolution. */
193
void HKL_info::generate_hkl_list()
195
std::vector<HKL> add;
196
// make a reflection list for the given cell, symm and resolution
197
/* TO MAKE A BOX TO HOLD A SPHERE IN REAL OR RECIPROCAL SPACE
198
In reciprocal space, use the real space cell dimension / resolution
199
In real space, use the recip space dimension * radius */
201
int hmax = int(cell_.descr().a()/resolution_.limit());
202
int kmax = int(cell_.descr().b()/resolution_.limit());
203
int lmax = int(cell_.descr().c()/resolution_.limit());
204
ftype s_lim = resolution_.invresolsq_limit();
205
// make a list of valid reflections
206
for (rfl.h()=-hmax; rfl.h()<=hmax; rfl.h()++)
207
for (rfl.k()=-kmax; rfl.k()<=kmax; rfl.k()++)
208
for (rfl.l()=-lmax; rfl.l()<=lmax; rfl.l()++)
209
if ( spacegroup_.recip_asu(rfl) && rfl.invresolsq(cell_) < s_lim &&
210
!( spacegroup_.hkl_class(rfl).sys_abs() ) ) add.push_back(rfl);
211
// update the reflection data lists
217
/*! The new HKLs are transformed to the default reciprocal ASU, and
218
added to the reflection list. Duplicates and reflections outside the
219
resoluution limit are ignored. Then the fast lookup tables for HKL,
220
invresolsq, and reflection class are rebuilt.
221
\param add The list of new reflections to add. */
222
void HKL_info::add_hkl_list( const std::vector<HKL>& add ) {
223
HKL equiv; int sym; bool friedel;
224
for ( int i = 0; i < add.size(); i++ ) {
225
if ( add[i].invresolsq( cell_ ) <= resolution_.invresolsq_limit() ) {
226
equiv = find_sym( add[i], sym, friedel );
227
if ( lookup.index_of( equiv ) < 0 ) hkl.push_back( equiv );
234
/*! Returns the index of the reflection, the sym no. and Friedel flag.
236
HKL HKL_info::find_sym( const HKL& rfl, int& sym, bool& friedel ) const
238
// find the symmetry operator mapping hkl into the reflection
239
// list and determine the friedel flag
242
// now find the symop which gives a reflection from the list
243
for (sym = 0; sym < spacegroup_.num_primops(); sym++) {
244
equiv = rfl.transform(isymop[sym]);
245
if ( spacegroup_.recip_asu(equiv) ) { friedel = false; return equiv; }
247
if ( spacegroup_.recip_asu(equiv) ) { friedel = true ; return equiv; }
249
Message::message( message_recip_asu_error );
254
void HKL_info::debug() const
256
std::cout << "Num reflns " << hkl.size() << "\n";
260
} // namespace clipper