00001 #ifndef CCTBX_MAPTBX_REAL_SPACE_GRADIENTS_SIMPLE_H
00002 #define CCTBX_MAPTBX_REAL_SPACE_GRADIENTS_SIMPLE_H
00003
00004 #include <cctbx/maptbx/eight_point_interpolation.h>
00005
00006 namespace cctbx { namespace maptbx {
00007
00008 template <
00009 typename MapFloatType,
00010 typename SiteFloatType>
00011 MapFloatType
00012 real_space_target_simple(
00013 uctbx::unit_cell const& unit_cell,
00014 af::const_ref<MapFloatType, af::c_grid_padded<3> > const& density_map,
00015 af::const_ref<scitbx::vec3<SiteFloatType> > const& sites_cart)
00016 {
00017 MapFloatType result = 0;
00018 for(std::size_t i_site=0;i_site<sites_cart.size();i_site++) {
00019 result += eight_point_interpolation(
00020 density_map,
00021 unit_cell.fractionalize(sites_cart[i_site]));
00022 }
00023 return result;
00024 }
00025
00026 template <
00027 typename MapFloatType,
00028 typename SiteFloatType>
00029 af::shared<scitbx::vec3<SiteFloatType> >
00030 real_space_gradients_simple(
00031 uctbx::unit_cell const& unit_cell,
00032 af::const_ref<MapFloatType, af::c_grid_padded<3> > const& density_map,
00033 af::const_ref<scitbx::vec3<SiteFloatType> > const& sites_cart,
00034 SiteFloatType delta)
00035 {
00036 CCTBX_ASSERT(delta > 0);
00037 typedef scitbx::vec3<SiteFloatType> v3_t;
00038 af::shared<v3_t> result(sites_cart.size(), af::init_functor_null<v3_t>());
00039 v3_t* res = result.begin();
00040 for(std::size_t i_site=0;i_site<sites_cart.size();i_site++,res++) {
00041 v3_t piv = sites_cart[i_site];
00042 v3_t piv_d = piv;
00043 for(unsigned i_axis=0;i_axis<3;i_axis++) {
00044 MapFloatType densities[2];
00045 for(unsigned i_sign=0;i_sign<2;i_sign++) {
00046 piv_d[i_axis] = (i_sign == 0 ? piv[i_axis] + delta
00047 : piv[i_axis] - delta);
00048 fractional<SiteFloatType> site_frac = unit_cell.fractionalize(piv_d);
00049 densities[i_sign] = eight_point_interpolation(
00050 density_map, site_frac);
00051 }
00052 piv_d[i_axis] = piv[i_axis];
00053 (*res)[i_axis] = (densities[0] - densities[1]) / (2 * delta);
00054 }
00055 }
00056 return result;
00057 }
00058
00059 }}
00060
00061 #endif // GUARD
