code: import itertools
  import numpy as np
  
  from collections import defaultdict
  
  
  def lagged_diff(x, k):
      """
      Returns the sequence of x[i] - x[i - k], as an array with the same size as x.
      Boundary conditions are handled as follows:
          x[i] = x[0]    if i < 0
          x[i] = x[n-1]  if i >= n, where n = len(x)
      """
      s = np.roll(x, k)
      if k >= 0:
          s[:k] = x[0]
      else:
          s[k:] = x[-1] # NOTE: lag is negative here
      return x - s
  
  
  def outlier_mask(x, threshold=3.0):
      """
      Returns an outlier mask for array x, based on Tukey's rule and assuming that the inlier
      distribution of x (the distribution of 'good' values) is Gaussian. 'threshold' represents a
      number of Gaussian sigmas.
      """
      q1, med, q3 = np.percentile(x, [25, 50, 75])
      std = (q3 - q1) / 1.349
      return (x - med) > threshold * std
  
  
  def genlags(radius, geofactor=1.5):
      lag = 1
      while lag <= radius:
          yield lag
          yield -lag
          lag = max(int(geofactor * lag), lag + 1)
  
  
  def iqrm_mask(x, radius=5, threshold=3.0):
      """
      Compute the IQRM mask for one-dimensional input data x.
      The input 'x' is expected to represent a per-channel statistic that measures RFI contamination
      in a block of time-frequency data. Any statistic can be used, but an important requirement is
      that larger values must indicate higher levels of RFI contamination.
  
      Parameters
      ----------
      x : list or ndarray
          Input data (1-dimensional)
      radius : int, optional
          Radius in number of elements. If a float is passed, it is truncated. A recommended value
          is 10% of the number of frequency channels
      threshold : float, optional
          Flagging threshold in number of Gaussian sigmas
  
      Returns
      -------
      mask : ndarray
          Boolean mask with the same size as the input 'x', where 'True' denotes an outlier
      votes_cast : dict of sets
          Dictionary of sets, where the keys are input array indices i that cast at least one vote,
          and the values are the set of array indices that received a vote from i.
      """
      x = np.asarray(x)
      n = len(x)
      radius = int(radius)
  
      if not radius > 0:
          raise ValueError("radius must be > 0")
  
      threshold = float(threshold)
      if not threshold > 0:
          raise ValueError("threshold must be > 0")
  
      # These data structures both represent a directed graph
      # votes_cast[i] contains the recipients of votes cast by i
      # votes_received[i] contains the casters of votes received by i
      votes_cast = defaultdict(set)
      votes_received = defaultdict(set)
  
      for lag in genlags(radius):
          d = lagged_diff(x, lag)
          m = outlier_mask(d, threshold)
  
          # m[i] = True  <=> point j = i - lag cast a vote on i
          #              <=> point i received a vote from j = i - lag
          I = np.where(m)[0]
          J = np.clip(I - lag, 0, n - 1)
  
          for i, j in zip(I, J):
              votes_cast[j].add(i)
              votes_received[i].add(j)
  
      mask = np.zeros_like(x, dtype=bool)
      
      # i gets masked by j if both the following conditions are True:
      # 1) j has cast a vote on i
      # 2) j has cast strictly less votes in total than i has received in total
      for i, casters in votes_received.items():
          for j in casters:
              if j in votes_cast and len(votes_cast[j]) < len(votes_received[i]):
                  mask[i] = True
                  break
  
      return mask, dict(votes_cast)
   import numpy as np
  from pytest import raises
  
  from iqrm import iqrm_mask
  
  
  def generate_noise(nchan=1024, seed=0):
      # IMPORTANT: set the random seed for reproducible results
      np.random.seed(seed)
      return np.random.normal(size=nchan)
  
  
  def generate_noise_with_outlier_range(start, end, nchan=1024, seed=0):
      s = generate_noise(nchan=nchan, seed=seed)
      s[start:end] += 100
      return s
  
  
  def test_param_checks():
      nchan = 1024
      s = np.zeros(nchan)
  
      with raises(ValueError): # radius must be > 0
          iqrm_mask(s, radius=0)
  
      with raises(ValueError): # threshold must be > 0
          iqrm_mask(s, threshold=0)
  
  
  def test_masking_noise():
      s = generate_noise()
  
      for radius in range(1, 6):
          mask, __ = iqrm_mask(s, radius=radius, threshold=4.0)
          assert np.alltrue(~mask)
  
  
  def test_masking_single_outlier():
      nchan = 1024
      indices = [0, 1, 42, 213, 740, 1022, 1023]
  
      for index in indices:
          # NOTE: when using radius = 1, if the either the first or last element are the sole 
          # outlier, they won't be detected (the single vote they receive is not valid). 
          # We thus start at radius=2.
          for radius in (2, 3, 4, 6, 9):
              s = generate_noise_with_outlier_range(index, index+1, nchan=nchan)
              mask, __ = iqrm_mask(s, radius=radius, threshold=4.0)
              assert mask[index] == True
  
  
  def test_masking_outlier_range():
      # The idea here is to generate data that looks like a top hat, i.e. noise plus a contiguous
      # range of high ouliers with similar values.
  
      # If the edges of the top-hat lie "far away" from the edges of the input array, then we expect
      # all outliers to be masked as long as:
      # max trial lag value > width
  
      # NOTE: for a top-hat that lies at the edge of the input array, the situation is different,
      # and the radius required to mask all outliers is such that:
      # max trial lag value > 2*width
  
      nchan = 1024
      indices = [67, 213, 486, 740, 959]
      trial_lag_sequence = (1, 2, 3, 4, 6, 9, 13)
  
      for index in indices:
          for jj, width in enumerate(trial_lag_sequence[:-1]):
              s = generate_noise_with_outlier_range(index, index+width, nchan=nchan)
              radius = trial_lag_sequence[jj+1]
              mask, __ = iqrm_mask(s, radius=radius, threshold=4.0)
              assert np.alltrue(mask[index:index+width])
apis: numpy.zeros_like numpy.roll numpy.asarray numpy.clip numpy.where numpy.percentile
apis: numpy.random.normal numpy.random.seed numpy.alltrue numpy.zeros