r"""Calculate local elastic matrix or stiffness matrix for local DNA segment

        .. note:: Here local DNA segment referred to less than 5 base-pair long.

        In case of :ref:`base-step-image`: Shift (:math:`Dx`), Slide (:math:`Dy`), Rise (:math:`Dz`),
        Tilt (:math:`\tau`), Roll (:math:`\rho`) and Twist (:math:`\omega`), following elastic matrix is calculated.

        .. math::

            \mathbf{K}_{base-step} = \begin{bmatrix}
            K_{Dx}        & K_{Dx,Dy}      & K_{Dx,Dz}      & K_{Dx,\tau}      & K_{Dx,\rho}      & K_{Dx,\omega} \\
            K_{Dx,Dy}     & K_{Dy}         & K_{Dy,Dz}      & K_{Dy,\tau}      & K_{Dy,\rho}      & K_{Dy,\omega} \\
            K_{Dx,Dz}     & K_{Dy,Dz}      & K_{Dz}         & K_{Dz,\tau}      & K_{Dz,\rho}      & K_{Dz,\omega} \\
            K_{Dx,\tau}   & K_{Dy,\tau}    & K_{Dz,\tau}    & K_{\tau}         & K_{\tau, \rho}   & K_{\tau,\omega} \\
            K_{Dx,\rho}   & K_{Dy,\rho}    & K_{Dz,\rho}    & K_{\tau, \rho}   & K_{\rho}         & K_{\rho,\omega} \\
            K_{Dx,\omega} & K_{Dy,\omega}  & K_{Dz,\omega}  & K_{\tau, \omega} & K_{\rho, \omega} & K_{\omega} \\
            \end{bmatrix}


        In case of :ref:`helical-base-step-image`: x-displacement (:math:`dx`), y-displacement (:math:`dy`), h-rise (:math:`h`),
        inclination (:math:`\eta`), tip (:math:`\theta`) and twist (:math:`\Omega`), following elastic matrix is calculated.

        .. math::

            \mathbf{K}_{helical-base-step} = \begin{bmatrix}
            K_{dx}        & K_{dx,dy}      & K_{dx,h}      & K_{dx,\eta}      & K_{dx,\theta}      & K_{dx,\Omega} \\
            K_{dx,dy}     & K_{dy}         & K_{dy,h}      & K_{dy,\eta}      & K_{dy,\theta}      & K_{dy,\Omega} \\
            K_{dx,h}      & K_{dy,h}       & K_{h}         & K_{h,\eta}       & K_{h,\theta}       & K_{h,\Omega} \\
            K_{dx,\eta}   & K_{dy,\eta}    & K_{h,\eta}    & K_{\eta}         & K_{\eta, \theta}   & K_{\eta,\Omega} \\
            K_{dx,\theta} & K_{dy,\theta}  & K_{h,\theta}  & K_{\eta, \theta} & K_{\theta}         & K_{\theta,\Omega} \\
            K_{dx,\Omega} & K_{dy,\Omega}  & K_{h,\Omega}  & K_{\eta, \Omega} & K_{\theta, \Omega} & K_{\Omega} \\
            \end{bmatrix}


        Parameters
        ----------
        bp : list
            List of two base-steps forming the DNA segment.
            For example: with ``bp=[5, 50]``, 5-50 base-step segment will be considered.

        frames : list
            List of two trajectory frames between which parameters will be extracted. It can be used to select portions
            of the trajectory. For example, with ``frames=[100, 1000]``, 100th to 1000th frame of the trajectory will be
            considered.

        helical : bool
            If ``helical=True``, elastic matrix for **helical base-step** parameters are calculated. Otherwise,
            by default, elastic matrix for **base-step** parameters are calculated.

        unit : str
            Unit of energy. Allowed units are: ``'kT', 'kJ/mol' and 'kcal/mol'``.

        Return
        ------
        mean : numpy.ndarray
            Value of parameters at which energy is zero. Minimum point on energy landscape.

            if ``helical=False``

            .. math::
                \begin{bmatrix}
                    Dx_0  &  Dy_0 & Dz_0 & \tau_0 & \rho_0 & \omega_0
                \end{bmatrix}

            if ``helical=True``

             .. math::
                \begin{bmatrix}
                    dx_0  &  dy_0 & h_0 & \eta_0 & \theta_0 & \Omega_0
                \end{bmatrix}

        result : numpy.ndarray
            Elastic matrix.

def calculateLocalElasticity(self, bp, frames=None, helical=False, unit='kT'):
        r"""Calculate local elastic matrix or stiffness matrix for local DNA segment

        .. note:: Here local DNA segment referred to less than 5 base-pair long.

        In case of :ref:`base-step-image`: Shift (:math:`Dx`), Slide (:math:`Dy`), Rise (:math:`Dz`),
        Tilt (:math:`\tau`), Roll (:math:`\rho`) and Twist (:math:`\omega`), following elastic matrix is calculated.

        .. math::

            \mathbf{K}_{base-step} = \begin{bmatrix}
            K_{Dx}        & K_{Dx,Dy}      & K_{Dx,Dz}      & K_{Dx,\tau}      & K_{Dx,\rho}      & K_{Dx,\omega} \\
            K_{Dx,Dy}     & K_{Dy}         & K_{Dy,Dz}      & K_{Dy,\tau}      & K_{Dy,\rho}      & K_{Dy,\omega} \\
            K_{Dx,Dz}     & K_{Dy,Dz}      & K_{Dz}         & K_{Dz,\tau}      & K_{Dz,\rho}      & K_{Dz,\omega} \\
            K_{Dx,\tau}   & K_{Dy,\tau}    & K_{Dz,\tau}    & K_{\tau}         & K_{\tau, \rho}   & K_{\tau,\omega} \\
            K_{Dx,\rho}   & K_{Dy,\rho}    & K_{Dz,\rho}    & K_{\tau, \rho}   & K_{\rho}         & K_{\rho,\omega} \\
            K_{Dx,\omega} & K_{Dy,\omega}  & K_{Dz,\omega}  & K_{\tau, \omega} & K_{\rho, \omega} & K_{\omega} \\
            \end{bmatrix}


        In case of :ref:`helical-base-step-image`: x-displacement (:math:`dx`), y-displacement (:math:`dy`), h-rise (:math:`h`),
        inclination (:math:`\eta`), tip (:math:`\theta`) and twist (:math:`\Omega`), following elastic matrix is calculated.

        .. math::

            \mathbf{K}_{helical-base-step} = \begin{bmatrix}
            K_{dx}        & K_{dx,dy}      & K_{dx,h}      & K_{dx,\eta}      & K_{dx,\theta}      & K_{dx,\Omega} \\
            K_{dx,dy}     & K_{dy}         & K_{dy,h}      & K_{dy,\eta}      & K_{dy,\theta}      & K_{dy,\Omega} \\
            K_{dx,h}      & K_{dy,h}       & K_{h}         & K_{h,\eta}       & K_{h,\theta}       & K_{h,\Omega} \\
            K_{dx,\eta}   & K_{dy,\eta}    & K_{h,\eta}    & K_{\eta}         & K_{\eta, \theta}   & K_{\eta,\Omega} \\
            K_{dx,\theta} & K_{dy,\theta}  & K_{h,\theta}  & K_{\eta, \theta} & K_{\theta}         & K_{\theta,\Omega} \\
            K_{dx,\Omega} & K_{dy,\Omega}  & K_{h,\Omega}  & K_{\eta, \Omega} & K_{\theta, \Omega} & K_{\Omega} \\
            \end{bmatrix}


        Parameters
        ----------
        bp : list
            List of two base-steps forming the DNA segment.
            For example: with ``bp=[5, 50]``, 5-50 base-step segment will be considered.

        frames : list
            List of two trajectory frames between which parameters will be extracted. It can be used to select portions
            of the trajectory. For example, with ``frames=[100, 1000]``, 100th to 1000th frame of the trajectory will be
            considered.

        helical : bool
            If ``helical=True``, elastic matrix for **helical base-step** parameters are calculated. Otherwise,
            by default, elastic matrix for **base-step** parameters are calculated.

        unit : str
            Unit of energy. Allowed units are: ``'kT', 'kJ/mol' and 'kcal/mol'``.

        Return
        ------
        mean : numpy.ndarray
            Value of parameters at which energy is zero. Minimum point on energy landscape.

            if ``helical=False``

            .. math::
                \begin{bmatrix}
                    Dx_0  &  Dy_0 & Dz_0 & \tau_0 & \rho_0 & \omega_0
                \end{bmatrix}

            if ``helical=True``

             .. math::
                \begin{bmatrix}
                    dx_0  &  dy_0 & h_0 & \eta_0 & \theta_0 & \Omega_0
                \end{bmatrix}

        result : numpy.ndarray
            Elastic matrix.

        """
        acceptedUnit = ['kT', 'kJ/mol', 'kcal/mol']

        if unit not in acceptedUnit:
            raise ValueError(" {0} not accepted. Use any of the following: {1} ".format(unit, acceptedUnit))

        frames = self._validateFrames(frames)

        name = '{0}-{1}-{2}-{3}-local-{4}'.format(bp[0], bp[1], frames[0], frames[1], int(helical))

        if bp[1]-bp[0]+1 > 4:
            raise ValueError("Selected span {0} is larger than 4, and therefore, not recommended for local elasticity".format(bp[1]-bp[0]+1))

        if name not in self.esMatrix:
            time, array = self.extractLocalParameters(self.dna, bp, helical=helical, frames=frames)
            mean = np.mean(array, axis = 1)
            esMatrix = self.getElasticMatrix(array)
            self.esMatrix[name] = esMatrix
            self.minimumPoint[name] = mean
        else:
            esMatrix = self.esMatrix[name]
            mean = self.minimumPoint[name]

        if unit == 'kJ/mol':
            result = 2.4946938107879997 * esMatrix      # (1.38064852e-23 * 300 * 6.023e23 / 1000 )  kT.NA/1000
        elif unit == 'kcal/mol':
            result = 0.5962461306854684 * esMatrix      # (1.38064852e-23 * 300 * 6.023e23 / 1000 / 4.184)  kT.NA/1000
        else:
            result = esMatrix

        return mean, result