How to use the pyuvsim.spherical_coordinates_basis_transformation.r_hat function in pyuvsim

Parameters
        ----------
        telescope_location : astropy EarthLocation object
            Location of the telescope.

        Returns
        -------
        array of floats
            Rotation matrix that defines the mapping (RA,Dec) <--> (Alt,Az),
            shape (3, 3, Ncomponents).
        """
        # Find mathematical points and vectors for RA/Dec
        theta_radec = np.pi / 2. - self.dec.rad
        phi_radec = self.ra.rad
        radec_vec = scbt.r_hat(theta_radec, phi_radec)
        assert radec_vec.shape == (3, self.Ncomponents)

        # Find mathematical points and vectors for Alt/Az
        theta_altaz = np.pi / 2. - self.alt_az[0, :]
        phi_altaz = self.alt_az[1, :]
        altaz_vec = scbt.r_hat(theta_altaz, phi_altaz)
        assert altaz_vec.shape == (3, self.Ncomponents)

        R_avg = self._calc_average_rotation_matrix(telescope_location)

        R_exact = np.zeros((3, 3, self.Ncomponents), dtype=np.float)
        for src_i in range(self.Ncomponents):
            intermediate_vec = np.matmul(R_avg, radec_vec[:, src_i])

            R_perturb = scbt.vecs2rot(r1=intermediate_vec, r2=altaz_vec[:, src_i])

Returns
        -------
        array of floats
            Rotation matrix that defines the mapping (RA,Dec) <--> (Alt,Az),
            shape (3, 3, Ncomponents).
        """
        # Find mathematical points and vectors for RA/Dec
        theta_radec = np.pi / 2. - self.dec.rad
        phi_radec = self.ra.rad
        radec_vec = scbt.r_hat(theta_radec, phi_radec)
        assert radec_vec.shape == (3, self.Ncomponents)

        # Find mathematical points and vectors for Alt/Az
        theta_altaz = np.pi / 2. - self.alt_az[0, :]
        phi_altaz = self.alt_az[1, :]
        altaz_vec = scbt.r_hat(theta_altaz, phi_altaz)
        assert altaz_vec.shape == (3, self.Ncomponents)

        R_avg = self._calc_average_rotation_matrix(telescope_location)

        R_exact = np.zeros((3, 3, self.Ncomponents), dtype=np.float)
        for src_i in range(self.Ncomponents):
            intermediate_vec = np.matmul(R_avg, radec_vec[:, src_i])

            R_perturb = scbt.vecs2rot(r1=intermediate_vec, r2=altaz_vec[:, src_i])

            R_exact[:, :, src_i] = np.matmul(R_perturb, R_avg)

        return R_exact