class_docstr: Differentiable global mutual information loss via Parzen windowing method.
  
  Reference:
      https://dspace.mit.edu/handle/1721.1/123142, Section 3.1, equation 3.1-3.5, Algorithm 1
class_docstr_tok: ["Differentiable","global","mutual","information","loss","via","Parzen","windowing","method",".","Reference",":","https",":","\/\/dspace.mit.edu\/handle\/1721.1\/123142",",","Section","3.1",",","equation","3.1-3.5",",","Algorithm","1"]
docstr: Parzen windowing with b-spline kernel (adapted from ITK)
  
  Args:
      img: the shape should be B[NDHW].
      order: int.
docstr_tok: ["Parzen","windowing","with","b-spline","kernel","(","adapted","from","ITK",")","Args",":","img",":","the","shape","should","be","B","[","NDHW","]",".","order",":","int","."]
code: def parzen_windowing_b_spline(
      self, img: torch.Tensor, order: int
  ) -> tuple[torch.Tensor, torch.Tensor]:
      """
      Parzen windowing with b-spline kernel (adapted from ITK)
  
      Args:
          img: the shape should be B[NDHW].
          order: int.
      """
  
      # Compute binsize for the histograms.
      #
      # The binsize for the image intensities needs to be adjusted so that
      # we can avoid dealing with boundary conditions using the cubic
      # spline as the Parzen window.  We do this by increasing the size
      # of the bins so that the joint histogram becomes "padded" at the
      # borders. Because we are changing the binsize,
      # we also need to shift the minimum by the padded amount in order to
      # avoid minimum values filling in our padded region.
      #
      # Note that there can still be non-zero bin values in the padded region,
      # it's just that these bins will never be a central bin for the Parzen
      # window.
      _max, _min = torch.max(img), torch.min(img)
      padding = 2
      bin_size = (_max - _min) / (self.num_bins - 2 * padding)
      norm_min = torch.div(_min, bin_size) - padding
  
      # assign bin/window index to each voxel
      window_term = torch.div(img, bin_size) - norm_min  # B[NDHW]
      # make sure the extreme values are in valid (non-padded) bins
      window_term = torch.clamp(
          window_term, padding, self.num_bins - padding - 1
      )  # B[NDHW]
      window_term = window_term.reshape(
          window_term.shape[0], -1, 1
      )  # (batch, num_sample, 1)
      bins = torch.arange(
          self.num_bins, device=window_term.device
      ).reshape(
          1, 1, -1
      )  # (1, 1, num_bins)
      sample_bin_matrix = torch.abs(
          bins - window_term
      )  # (batch, num_sample, num_bins)
  
      # b-spleen kernel
      # (4 - 6 * abs ** 2 + 3 * abs ** 3) / 6 when 0 <= abs < 1
      # (2 - abs) ** 3 / 6 when 1 <= abs < 2
      weight = torch.zeros_like(
          sample_bin_matrix, dtype=torch.float
      )  # (batch, num_sample, num_bins)
      if order == 0:
          weight = (
              weight
              + (sample_bin_matrix < 0.5)
              + (sample_bin_matrix == 0.5) * 0.5
          )
      elif order == 3:
          weight = (
              weight
              + (
                  4
                  - 6 * sample_bin_matrix**2
                  + 3 * sample_bin_matrix**3
              )
              * (sample_bin_matrix < 1)
              / 6
          )
          weight = (
              weight
              + (2 - sample_bin_matrix) ** 3
              * (sample_bin_matrix >= 1)
              * (sample_bin_matrix < 2)
              / 6
          )
      else:
          raise ValueError(
              f"Do not support b-spline {order}-order parzen windowing"
          )
  
      weight = weight / torch.sum(
          weight, dim=-1, keepdim=True
      )  # (batch, num_sample, num_bins)
      probability = torch.mean(
          weight, dim=-2, keepdim=True
      )  # (batch, 1, num_bins)
      return weight, probability
  
imports: from __future__ import annotations
  import torch
  from torch.nn import functional
  from torch.nn.modules.loss import _Loss
  from monai.networks.layers import gaussian_1d, separable_filtering
  from monai.utils import LossReduction
  from monai.utils.module import look_up_option
label_desc: Use image node_id 3 for calling the GlobalMutualInformationLoss obj's underlying member method code with example usage: obj.parzen_windowing_b_spline(img, order) and returns: weight,  probability