<commit_before>from __future__ import absolute_import, division, print_function
# LIBTBX_PRE_DISPATCHER_INCLUDE_SH export PHENIX_GUI_ENVIRONMENT=1
# LIBTBX_PRE_DISPATCHER_INCLUDE_SH export BOOST_ADAPTBX_FPE_DEFAULT=1
import math
import iotbx.phil
from scitbx import matrix
from cctbx.array_family import flex
from dials.util.options import OptionParser
from dials.util.options import flatten_datablocks, flatten_reflections
from dials.algorithms.indexing.indexer \
     import indexer_base, filter_reflections_by_scan_range

import libtbx.load_env
import logging
logger = logging.getLogger('dials.command_line.search_beam_position')

help_message = '''

Search for a better beam centre using the results of spot finding. Based on
method of Sauter et al., J. Appl. Cryst. 37, 399-409 (2004).

Examples::

  %s datablock.json strong.pickle

''' % libtbx.env.dispatcher_name

phil_scope = iotbx.phil.parse("""
nproc = Auto
  .type = int(value_min=1)
plot_search_scope = False
  .type = bool
max_cell = None
  .type = float
  .help = "Known max cell (otherwise will compute from spot positions)"
scan_range = None
  .help = "The range of images to use in indexing. Number of arguments"
    "must be a factor of two. Specifying \"0 0\" will use all images"
    "by default. The given range follows C conventions"
    "(e.g. j0 <= j < j1)."
  .type = ints(size=2)
  .multiple = True
max_reflections = 10000
  .type = int(value_min=1)
  .help = "Maximum number of reflections to use in the search for better"
          "experimental model. If the number of input reflections is greater"
          "then a random subset of reflections will be used."
mm_search_scope = 4.0
  .help = "Global radius of origin offset search."
  .type = float(value_min=0)
wide_search_binning = 2
  .help = "Modify the coarseness of the wide grid search for the beam centre."
  .type = float(value_min=0)
n_macro_cycles = 1
  .type = int
  .help = "Number of macro cycles for an iterative beam centre search."

seed = 42
  .type = int(value_min=0)

output {
  datablock = optimized_datablock.json
    .type = path
  log = "dials.search_beam_position.log"
    .type = str
  debug_log = "dials.search_beam_position.debug.log"
    .type = str
}
""")

master_params = phil_scope.fetch().extract()


from rstbx.phil.phil_preferences import indexing_api_defs
dps_phil_scope = iotbx.phil.parse('''
include scope rstbx.phil.phil_preferences.indexing_api_defs
d_min = None
  .type = float(value_min=0)
''', process_includes=True)


class better_experimental_model_discovery(object):
  def __init__(self, imagesets, spot_lists, solution_lists,
               amax_lists, horizon_phil, wide_search_binning=1):
    from libtbx import adopt_init_args
    adopt_init_args(self, locals())

  def optimize_origin_offset_local_scope(self):
    """Local scope: find the optimal origin-offset closest to the current overall detector position
       (local minimum, simple minimization)"""

    from libtbx.test_utils import approx_equal
    from rstbx.indexing_api import dps_extended

    detector = self.imagesets[0].get_detector()
    beam = self.imagesets[0].get_beam()
    s0 = matrix.col(beam.get_s0())
    # construct two vectors that are perpendicular to the beam.  Gives a basis for refining beam
    axis = matrix.col((1,0,0))
    beamr0 = s0.cross(axis).normalize()
    beamr1 = beamr0.cross(s0).normalize()
    beamr2 = beamr1.cross(s0).normalize()

    assert approx_equal(s0.dot(beamr1), 0.)
    assert approx_equal(s0.dot(beamr2), 0.)
    assert approx_equal(beamr2.dot(beamr1), 0.)
    # so the orthonormal vectors are self.S0_vector, beamr1 and beamr2

    if self.horizon_phil.indexing.mm_search_scope:
      scope = self.horizon_phil.indexing.mm_search_scope
      plot_px_sz = self.imagesets[0].get_detector()[0].get_pixel_size()[0]
      plot_px_sz *= self.wide_search_binning
      grid = max(1,int(scope/plot_px_sz))
      widegrid = 2 * grid + 1
      scores = flex.double()
      for y in xrange(-grid,grid+1):
        for x in xrange(-grid,grid+1):
          new_origin_offset = x*plot_px_sz*beamr1 + y*plot_px_sz*beamr2
          score = 0
          for i in range(len(self.imagesets)):
            score += self.get_origin_offset_score(new_origin_offset,
                                                  self.solution_lists[i],
                                                  self.amax_lists[i],
                                                  self.spot_lists[i],
                                                  self.imagesets[i])
          scores.append(score)

      def igrid(x): return x - (widegrid//2)

      idxs = [igrid(i)*plot_px_sz for i in xrange(widegrid)]

      # if there are several similarly high scores, then choose the closest
      # one to the current beam centre
      potential_offsets = flex.vec3_double()
      sel = scores > (0.9*flex.max(scores))
      for i in sel.iselection():
        offset = (idxs[i%widegrid])*beamr1 + (idxs[i//widegrid])*beamr2
        potential_offsets.append(offset.elems)
        #print offset.length(), scores[i]
      wide_search_offset = matrix.col(
        potential_offsets[flex.min_index(potential_offsets.norms())])

      #plot_max = flex.max(scores)
      #idx_max = flex.max_index(scores)
      #wide_search_offset = (idxs[idx_max%widegrid])*beamr1 + (idxs[idx_max//widegrid])*beamr2

    else:
      wide_search_offset = None

    # DO A SIMPLEX MINIMIZATION
    from scitbx.simplex import simplex_opt
    class test_simplex_method(object):
      def __init__(selfOO, wide_search_offset=None):
        selfOO.starting_simplex=[]
        selfOO.n = 2
        selfOO.wide_search_offset = wide_search_offset
        for ii in range(selfOO.n+1):
          selfOO.starting_simplex.append(flex.random_double(selfOO.n))
        selfOO.optimizer = simplex_opt(dimension=selfOO.n,
                                       matrix=selfOO.starting_simplex,<commit_msg>raise Sorry if all scores are equal to zero (closes #628)<commit_after>from __future__ import absolute_import, division, print_function
# LIBTBX_PRE_DISPATCHER_INCLUDE_SH export PHENIX_GUI_ENVIRONMENT=1
# LIBTBX_PRE_DISPATCHER_INCLUDE_SH export BOOST_ADAPTBX_FPE_DEFAULT=1
import math
import iotbx.phil
from scitbx import matrix
from cctbx.array_family import flex
from libtbx.utils import Sorry
from dials.util.options import OptionParser
from dials.util.options import flatten_datablocks, flatten_reflections
from dials.algorithms.indexing.indexer \
     import indexer_base, filter_reflections_by_scan_range

import libtbx.load_env
import logging
logger = logging.getLogger('dials.command_line.search_beam_position')

help_message = '''

Search for a better beam centre using the results of spot finding. Based on
method of Sauter et al., J. Appl. Cryst. 37, 399-409 (2004).

Examples::

  %s datablock.json strong.pickle

''' % libtbx.env.dispatcher_name

phil_scope = iotbx.phil.parse("""
nproc = Auto
  .type = int(value_min=1)
plot_search_scope = False
  .type = bool
max_cell = None
  .type = float
  .help = "Known max cell (otherwise will compute from spot positions)"
scan_range = None
  .help = "The range of images to use in indexing. Number of arguments"
    "must be a factor of two. Specifying \"0 0\" will use all images"
    "by default. The given range follows C conventions"
    "(e.g. j0 <= j < j1)."
  .type = ints(size=2)
  .multiple = True
max_reflections = 10000
  .type = int(value_min=1)
  .help = "Maximum number of reflections to use in the search for better"
          "experimental model. If the number of input reflections is greater"
          "then a random subset of reflections will be used."
mm_search_scope = 4.0
  .help = "Global radius of origin offset search."
  .type = float(value_min=0)
wide_search_binning = 2
  .help = "Modify the coarseness of the wide grid search for the beam centre."
  .type = float(value_min=0)
n_macro_cycles = 1
  .type = int
  .help = "Number of macro cycles for an iterative beam centre search."

seed = 42
  .type = int(value_min=0)

output {
  datablock = optimized_datablock.json
    .type = path
  log = "dials.search_beam_position.log"
    .type = str
  debug_log = "dials.search_beam_position.debug.log"
    .type = str
}
""")

master_params = phil_scope.fetch().extract()


from rstbx.phil.phil_preferences import indexing_api_defs
dps_phil_scope = iotbx.phil.parse('''
include scope rstbx.phil.phil_preferences.indexing_api_defs
d_min = None
  .type = float(value_min=0)
''', process_includes=True)


class better_experimental_model_discovery(object):
  def __init__(self, imagesets, spot_lists, solution_lists,
               amax_lists, horizon_phil, wide_search_binning=1):
    from libtbx import adopt_init_args
    adopt_init_args(self, locals())

  def optimize_origin_offset_local_scope(self):
    """Local scope: find the optimal origin-offset closest to the current overall detector position
       (local minimum, simple minimization)"""

    from libtbx.test_utils import approx_equal
    from rstbx.indexing_api import dps_extended

    detector = self.imagesets[0].get_detector()
    beam = self.imagesets[0].get_beam()
    s0 = matrix.col(beam.get_s0())
    # construct two vectors that are perpendicular to the beam.  Gives a basis for refining beam
    axis = matrix.col((1,0,0))
    beamr0 = s0.cross(axis).normalize()
    beamr1 = beamr0.cross(s0).normalize()
    beamr2 = beamr1.cross(s0).normalize()

    assert approx_equal(s0.dot(beamr1), 0.)
    assert approx_equal(s0.dot(beamr2), 0.)
    assert approx_equal(beamr2.dot(beamr1), 0.)
    # so the orthonormal vectors are self.S0_vector, beamr1 and beamr2

    if self.horizon_phil.indexing.mm_search_scope:
      scope = self.horizon_phil.indexing.mm_search_scope
      plot_px_sz = self.imagesets[0].get_detector()[0].get_pixel_size()[0]
      plot_px_sz *= self.wide_search_binning
      grid = max(1,int(scope/plot_px_sz))
      widegrid = 2 * grid + 1
      scores = flex.double()
      for y in xrange(-grid,grid+1):
        for x in xrange(-grid,grid+1):
          new_origin_offset = x*plot_px_sz*beamr1 + y*plot_px_sz*beamr2
          score = 0
          for i in range(len(self.imagesets)):
            score += self.get_origin_offset_score(new_origin_offset,
                                                  self.solution_lists[i],
                                                  self.amax_lists[i],
                                                  self.spot_lists[i],
                                                  self.imagesets[i])
          scores.append(score)

      def igrid(x): return x - (widegrid//2)

      idxs = [igrid(i)*plot_px_sz for i in xrange(widegrid)]

      # if there are several similarly high scores, then choose the closest
      # one to the current beam centre
      potential_offsets = flex.vec3_double()
      if scores.all_eq(0):
        raise Sorry("No valid scores")
      sel = scores > (0.9*flex.max(scores))
      for i in sel.iselection():
        offset = (idxs[i%widegrid])*beamr1 + (idxs[i//widegrid])*beamr2
        potential_offsets.append(offset.elems)
        #print offset.length(), scores[i]
      wide_search_offset = matrix.col(
        potential_offsets[flex.min_index(potential_offsets.norms())])

      #plot_max = flex.max(scores)
      #idx_max = flex.max_index(scores)
      #wide_search_offset = (idxs[idx_max%widegrid])*beamr1 + (idxs[idx_max//widegrid])*beamr2

    else:
      wide_search_offset = None

    # DO A SIMPLEX MINIMIZATION
    from scitbx.simplex import simplex_opt
    class test_simplex_method(object):
      def __init__(selfOO, wide_search_offset=None):
        selfOO.starting_simplex=[]
        selfOO.n = 2
        selfOO.wide_search_offset = wide_search_offset
        for ii in range(selfOO.n+1):
          selfOO.starting_simplex.append(flex.random_double(selfOO.n))
        selfOO.optimizer = simplex_opt(dimension=selfOO.n,
                                       matrix=selfOO.starting_simplex,
{'fraction_non_alphanumeric': 0.09519029694688415, 'fraction_numerical': 0.014554579673776663, 'mean_word_length': 3.101543739279588, 'pattern_counts': {'":': 0, '<': 7, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 18, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '3814114', 'n_tokens_mistral': 4074, 'n_tokens_neox': 3850, 'n_words': 951}