diff --git a/lib/ants/cli/ancil_create_shapefile.py b/lib/ants/cli/ancil_create_shapefile.py
index 7081626..3cd817a 100755
--- a/lib/ants/cli/ancil_create_shapefile.py
+++ b/lib/ants/cli/ancil_create_shapefile.py
@@ -9,53 +9,377 @@
 
 Creates and saves a shapefile from a list of pairs of longitude, latitude
 points defining a single polygon in a specified polygon file.
+
+Rotated pole domains can be specified using the land sea mask argument,
+where the longitude, latitude pairs are transformed to the new pole.
+Unless a source cube is provided with the source coordinate reference system,
+it assumed that the points defined in the json file are on a standard spherical
+unrotated geodetic coordinate reference system.
 """
 import argparse
 import json
+import logging
+import os
+import warnings
 
+import ants
+import iris.coord_systems
 import numpy as np
+from ants.io.load import load_cube
+from ants.utils.cube import CubeBuilder
 from osgeo import ogr
 from shapely.geometry import Polygon
+from shapely.validation import explain_validity
+
+_LOGGER = logging.getLogger(__name__)
+
+
+def _check_coord_system_type(target_lsm):
+    """
+    Check that target_lsm has a rotated pole coordinate system.
+
+    Parameters
+    ----------
+    target_lsm : :class:`iris.cube.Cube`
+        The lsm cube specifiying the rotated pole coordinate system.
+
+    Raises
+    ------
+    ValueError :
+        If target_lsm does not have a rotated pole coordinate system.
+    """
+
+    is_pole_coords = isinstance(
+        target_lsm.coord_system(), iris.coord_systems.RotatedGeogCS
+    )
+    if not is_pole_coords:
+        raise ValueError(
+            f"target_lsm.coord_system() {target_lsm.coord_system()} is not"
+            f" an instance of {iris.coord_systems.RotatedGeogCS}."
+            f" The landsea mask should specify a valid rotated pole coordinate"
+            f" system."
+        )
+
+
+def _check_polygon_validity(polygon, ccw_expected):
+    """
+    Check if polygon is valid and has the expected orientation.
+
+    A rotation should not change the orientation of the points. If the orientation
+    does change, it is a sign that the transformation has not behaved as expected.
+
+    Parameters
+    ----------
+    polygon : :class:`~shapely.geometry.Polygon`
+        The polygon made from the input json file.
+    ccw_expected : bool
+        Expected polygon orientation before transformation.
+
+    Raises
+    ------
+    ValueError :
+        If the polygon is invalid (for example intersecting edges).
+    ValueError :
+        If the polygon has a different orientation after transformation.
+    """
 
+    if not polygon.is_valid:
+        raise ValueError(f"Polygon is invalid: {explain_validity(polygon)}")
+    if polygon.exterior.is_ccw != ccw_expected:
+        raise ValueError(
+            f"Polygon orientation has changed. Expected is_ccw={ccw_expected}, "
+            f"current polygon has is_ccw={polygon.exterior.is_ccw}."
+        )
 
-def _load_polygon_from_json(json_file):
+
+def _validate_orientation(target_lsm):
     """
-    Load a json file containing a list of pairs of longitude, latitude points
-    to create a polygon from.
+    Validate that the rotated pole coordinate system is not coincident with the
+    geographic north pole.
+
+    Parameters
+    ----------
+    target_lsm : :class:`iris.cube.Cube`
+        The target lsm cube specifying the rotated pole coordinate system.
+
+    Returns
+    -------
+    bool
+        True if target_lsm is a valid rotated pole, False otherwise.
+
+    Warns
+    -----
+    UserWarning
+        If the rotated pole is coincident with the geographic north pole.
+    """
+
+    target_crs = target_lsm.coord_system()
+
+    grid_lon = target_crs.grid_north_pole_longitude
+    grid_lat = target_crs.grid_north_pole_latitude
+    pole_lon_rotation = target_crs.north_pole_grid_longitude
+
+    invalid_crs = ants.utils.ndarray.allclose(
+        [grid_lon, grid_lat, pole_lon_rotation], [0.0, 90.0, 0.0]
+    )
+
+    if invalid_crs:
+        warnings.warn(
+            "target_lsm has a geodetic coordinate system with pole located"
+            f" at grid_longitude={grid_lon}, grid_latitude={grid_lat}."
+            " No transformation will be carried out."
+        )
+    return not invalid_crs
+
+
+def _validate_args(target_lsm_path, source_path):
+    if source_path is not None and target_lsm_path is None:
+        raise ValueError("If --source is passed then --target-lsm must also be given.")
+
+
+def _transform_coordinates(target_lsm, source, points):
+    """
+    Transform the longitude-latitude points to the rotated pole.
+
+    If the transformation results in the points spanning the antimeridian, this
+    can cause incorrect polygons to be created because shapely does not wrap
+    between 180.0 and -180.0 degrees.
+
+    If the distance between sequential transformed longitude points exceeds
+    180.0 degrees, it is assumed these points cross the antimeridian. The points
+    are instead defined in the interval [0.0, 360.0] by adding 360.0 degrees
+    to the negative longitudes.
+
+    Parameters
+    ----------
+    target_lsm : :class:`iris.cube.Cube`
+        The lsm cube specifying the rotated pole coordinate system.
+    source : :class:`iris.cube.Cube`
+        An iris cube specifying the coordinate system of the input json
+        file.
+    points : :class:`numpy.ndarray`
+        An ``(m, 2)`` numpy array of m longitude-latitude pairs.
+
+    Returns
+    -------
+    : :class:`numpy.ndarray`
+        An ``(m, 2)`` numpy array of transformed longitude-latitude pairs.
+    Warns
+    -----
+    UserWarning
+        If the transformed points are assumed to cross the antimeridian.
+    UserWarning
+        If the transformed polygon lies at least partially outside of the
+        domain specified in the target lsm.
+    """
+
+    source_crs = source.coord_system().as_cartopy_crs()
+    target_coord = target_lsm.coord_system()
+    target_crs = target_coord.as_cartopy_crs()
+
+    # We only return longitude and latitude and discard the z coordinate.
+    rotated_points = target_crs.transform_points(
+        source_crs, points[:, 0], points[:, 1]
+    )[:, :2]
+
+    # Create a wrapped array of longitudes and find the point-wise difference.
+    closed_lon = np.vstack([rotated_points, rotated_points[0, :]])
+    lon_diff = np.abs(closed_lon[:-1, 0] - closed_lon[1:, 0])
+
+    # Add 360.0 degrees to negative longitudes if the points cross the antimeridian.
+    if np.any(lon_diff > 180.0):
+        warnings.warn(
+            "The transformed points are assumed to cross the antimeridian. "
+            "The longitudinal points will instead be defined in the interval"
+            " [0, 360.0] degrees."
+        )
+        neg_indices = np.where(rotated_points[:, 0] < 0)
+        rotated_points[neg_indices, 0] += 360.0
+
+    # Find the bounds of the target lsm.
+    bounds_lon = target_lsm.coord(axis="X").bounds
+    bounds_lat = target_lsm.coord(axis="Y").bounds
+    min_lon, max_lon = bounds_lon.min(), bounds_lon.max()
+    min_lat, max_lat = bounds_lat.min(), bounds_lat.max()
+    lons, lats = rotated_points[:, 0], rotated_points[:, 1]
+
+    # Check if the polygon lies in the domain specified by the target lsm.
+    if (
+        lons.min() < min_lon
+        or lons.max() > max_lon
+        or lats.min() < min_lat
+        or lats.max() > max_lat
+    ):
+        warnings.warn("The transformed points lie outside the target lsm domain.")
+
+    _LOGGER.info(
+        "Input json file transformed to new pole rotated coordinate system at "
+        "pole longitude=%s, pole latitude=%s, central rotated longitude=%s.",
+        target_coord.grid_north_pole_longitude,
+        target_coord.grid_north_pole_latitude,
+        target_coord.north_pole_grid_longitude,
+    )
+
+    return rotated_points
+
+
+def _load_cubes(target_lsm_path, source_path):
+    """
+    Load the target lsm and create a source cube if not provided.
+
+    Parameters
+    ----------
+    target_lsm_path : str
+        File path to a land sea mask that provides the new rotated pole.
+    source_path : str
+        File path to a source file specifying the coordinate system of the
+        input json file.
+
+    Returns
+    -------
+    : tuple(:class:`iris.cube.Cube`, :class:`iris.cube.Cube`)
+        A tuple containing the target lsm and the source cube respectively.
+    """
+
+    target_lsm = load_cube(target_lsm_path)
+
+    if source_path is None:
+        crs = iris.coord_systems.GeogCS(6371229.0)
+        source = CubeBuilder(crs, (2, 2))._cube
+    else:
+        source = load_cube(source_path)
+
+    return target_lsm, source
+
+
+def _transform_if_required(target_lsm, source, points):
+    """
+    Perform the transformation to a rotated pole if the coordinate system
+    is valid.
+
+    The points are transformed to the specified rotated pole coordinate system
+    if the provided coordinate system is a rotated pole and is not coincident
+    with the geographic north pole.
+
+    Parameters
+    ----------
+    target_lsm : :class:`iris.cube.Cube`
+        The lsm cube specifying the rotated pole coordinate system.
+    source : :class:`iris.cube.Cube`
+        An iris cube specifying the coordinate system of the input json
+        file.
+    points : :class:`numpy.ndarray`
+        An ``(m, 2)`` numpy array of m longitude-latitude pairs.
+
+    Returns
+    -------
+    : :class:`numpy.ndarray`
+        An ``(m, 2)`` numpy array of transformed longitude-latitude pairs.
+    """
+
+    _check_coord_system_type(target_lsm)
+
+    if _validate_orientation(target_lsm):
+        rotated_points = _transform_coordinates(target_lsm, source, points)
+    else:
+        rotated_points = np.copy(points)
+
+    return rotated_points
+
+
+def _save_json(output, target_lsm_path, source_path, target_lsm, source):
+    """
+    Save a json file containing metadata about the target lsm and source.
+    """
+
+    parent_path = os.path.dirname(output)
+    filename = os.path.splitext(os.path.basename(output))[0]
+
+    target_crs = target_lsm.coord_system()
+    target_proj_params = target_crs.as_cartopy_crs().proj4_params
+    source_crs = source.coord_system()
+    source_proj_params = source_crs.as_cartopy_crs().proj4_params
+
+    prj_metadata = {
+        "target_lsm": {
+            "target_lsm_path": target_lsm_path,
+            "coord_system": type(target_crs).__name__,
+            "grid_north_pole_longitude": target_crs.grid_north_pole_longitude,
+            "grid_north_pole_latitude": target_crs.grid_north_pole_latitude,
+            "north_pole_grid_longitude": target_crs.north_pole_grid_longitude,
+            "proj4_params": target_proj_params,
+        },
+        "source": {
+            "source_path": source_path,
+            "coord_system": type(source_crs).__name__,
+            "proj4_params": source_proj_params,
+        },
+    }
+    json_name = filename + "_prj"
+
+    with open(os.path.join(parent_path, json_name + ".json"), "w") as json_file:
+        json.dump(prj_metadata, json_file, indent=4)
+
+
+def _load_points_from_json(json_file):
+    """
+    Load a json file containing a list of pairs of longitude-latitude points.
 
     Parameters
     ----------
     json_file : str
-        Path to json file
+        Path to json file.
 
     Returns
     -------
-    : :class:`~shapely.geometry.Polygon`
+    : :class:`numpy.ndarray`
     """
+
     with open(json_file, "r") as polygon_json:
-        polygon = json.load(polygon_json)
-    polygon = np.array(polygon)
-    polygon = Polygon(polygon)
-    return polygon
+        points = json.load(polygon_json)
+    points = np.array(points)
 
+    return points
 
-def main(json_file, output):
+
+def main(json_file, output, target_lsm_path, source_path):
     """
+    Create a shape file from pairs of longitude, latitude points.
+
     Loads in a provided json file that defines pairs of longitude, latitude
     points to create a polygon from. That polygon is then used to create a
     shape file that is saved to the specified output location.
 
+    If target_lsm_path is provided, the points are first transformed from a
+    source geodetic coordinate system to a rotated pole coordinate system
+    specified by the lsm. It is assumed that the points in the json file
+    are on an unrotated spherical geodetic grid, unless otherwise specified.
+
     Parameters
     ----------
     json_file : str
         Path to json file
     output : str
         Location to store generated shape file
-
+    target_lsm_path : str
+        File path to a land sea mask that provides the new rotated pole.
+    source_path : str
+        File path to a source file specifying the coordinate system of the
+        input json file.
     """
 
     # Load a json and make a polygon
-    polygon = _load_polygon_from_json(json_file)
+    points = _load_points_from_json(json_file)
+    ccw_expected = Polygon(points).exterior.is_ccw
+
+    # Transform points to a rotated pole if required
+    if target_lsm_path is not None:
+        target_lsm, source = _load_cubes(target_lsm_path, source_path)
+        points = _transform_if_required(target_lsm, source, points)
+        _save_json(output, target_lsm_path, source_path, target_lsm, source)
+
+    polygon = Polygon(points)
+    _check_polygon_validity(polygon, ccw_expected)
 
     # Now convert it to a shapefile with OGR
     driver = ogr.GetDriverByName("Esri Shapefile")
@@ -87,13 +411,28 @@ def _get_parser():
         "json_file", help="Path to json file defining polygon to generate."
     )
     parser.add_argument("output", help="File to save shape file to.")
+    parser.add_argument(
+        "--target-lsm",
+        type=ants.config.filepath_readable,
+        required=False,
+        help="Path to the land sea mask containing the rotated pole"
+        " coordinate system.",
+    )
+    parser.add_argument(
+        "--source",
+        type=ants.config.filepath_readable,
+        required=False,
+        help="Path to a source specifying the coordinate system of the json file.",
+    )
     return parser
 
 
 def cli_interface():
     parser = _get_parser()
     args = parser.parse_args()
-    main(args.json_file, args.output)
+
+    _validate_args(args.target_lsm, args.source)
+    main(args.json_file, args.output, args.target_lsm, args.source)
 
 
 if __name__ == "__main__":
diff --git a/lib/ants/tests/cli/__init__.py b/lib/ants/tests/cli/__init__.py
new file mode 100644
index 0000000..89193a5
--- /dev/null
+++ b/lib/ants/tests/cli/__init__.py
@@ -0,0 +1,4 @@
+# (C) Crown Copyright, Met Office. All rights reserved.
+#
+# This file is part of ANTS and is released under the BSD 3-Clause license.
+# See LICENSE.txt in the root of the repository for full licensing details.
diff --git a/lib/ants/tests/cli/test_ancil_create_shapefile.py b/lib/ants/tests/cli/test_ancil_create_shapefile.py
new file mode 100644
index 0000000..e339307
--- /dev/null
+++ b/lib/ants/tests/cli/test_ancil_create_shapefile.py
@@ -0,0 +1,386 @@
+# (C) Crown Copyright, Met Office. All rights reserved.
+#
+# This file is part of ANTS and is released under the BSD 3-Clause license.
+# See LICENSE.txt in the root of the repository for full licensing details.
+
+import argparse
+import re
+from unittest import mock
+
+import ants.tests
+import iris
+import numpy as np
+from ants.cli.ancil_create_shapefile import (
+    _check_coord_system_type,
+    _check_polygon_validity,
+    _load_cubes,
+    _transform_coordinates,
+    _transform_if_required,
+    _validate_args,
+    _validate_orientation,
+)
+from ants.tests.stock import geodetic
+from ants.utils.cube import CubeBuilder
+from shapely.geometry import Polygon
+
+
+class Test__check_coord_system_type(ants.tests.TestCase):
+
+    def test_unrotated_target_lsm(self):
+        """Test that an error message is raised if the coordinate system in
+        target_lsm is not a rotated pole."""
+
+        target_lsm = geodetic((2, 2))
+        error_msg = re.escape(
+            f"target_lsm.coord_system() {target_lsm.coord_system()} is not"
+            f" an instance of {iris.coord_systems.RotatedGeogCS}."
+            f" The landsea mask should specify a valid rotated pole coordinate"
+            f" system."
+        )
+
+        with self.assertRaisesRegex(ValueError, error_msg):
+            _check_coord_system_type(target_lsm)
+
+
+class Test__check_polygon_validity(ants.tests.TestCase):
+    def test_invalid_polygon(self):
+        """Test that an invalid polygon with intersecting boundaries
+        raises an error."""
+
+        invalid_points = np.array([[1, 1], [-1, -1], [1, -1], [-1, 1]])
+        invalid_poly = Polygon(invalid_points)
+        is_ccw = False
+
+        error_msg = "Polygon is invalid: "
+
+        with self.assertRaisesRegex(ValueError, error_msg):
+            _check_polygon_validity(invalid_poly, is_ccw)
+
+    def test_changed_orientation(self):
+        """Test that a change in polygon orientation raises an error."""
+
+        ccw_points = np.array([[1, 1], [-1, 1], [-1, -1], [1, -1]])
+        ccw_poly = Polygon(ccw_points)
+        is_ccw = False
+
+        error_msg = (
+            "Polygon orientation has changed. Expected is_ccw=False, "
+            "current polygon has is_ccw=True."
+        )
+
+        with self.assertRaisesRegex(ValueError, error_msg):
+            _check_polygon_validity(ccw_poly, is_ccw)
+
+
+class Test__validate_orientation(ants.tests.TestCase):
+    def setUp(self):
+        self.sphere_crs = iris.coord_systems.GeogCS(6371229.0)
+        self.sphere_identity_crs = iris.coord_systems.RotatedGeogCS(
+            90.0, 0.0, ellipsoid=self.sphere_crs
+        )
+
+    def test_warning_raised(self):
+        """Test that a warning is raised with a coordinate system at latitude=90.0,
+        longitude=0.0 is passed. Check the function returns False."""
+
+        target_lsm = CubeBuilder(self.sphere_identity_crs, (2, 2))._cube
+
+        warning_msg = (
+            "target_lsm has a geodetic coordinate system with pole located"
+            " at grid_longitude=0.0, grid_latitude=90.0."
+            " No transformation will be carried out."
+        )
+
+        with self.assertWarnsRegex(UserWarning, warning_msg):
+            valid_crs = _validate_orientation(target_lsm)
+        self.assertFalse(valid_crs)
+
+    def test_true_returned(self):
+        """Test that True is returned when target_lsm has a valid rotated pole."""
+
+        target_lsm = geodetic(
+            (2, 2), north_pole_lat=90.0, north_pole_lon=90.0, crs=self.sphere_crs
+        )
+
+        self.assertTrue(_validate_orientation(target_lsm))
+
+    def test_no_warning_raised_rotated(self):
+        """Test that passing a non-zero central rotated longitude is still
+        accounted for when checking if the coordinate system is rotated."""
+
+        crs = iris.coord_systems.RotatedGeogCS(
+            90.0, 0.0, 180.0, ellipsoid=self.sphere_crs
+        )
+        target_lsm = CubeBuilder(crs, (2, 2))._cube
+
+        self.assertTrue(_validate_orientation(target_lsm))
+
+
+class Test__transform_coordinates(ants.tests.TestCase):
+    def setUp(self):
+        self.sphere_crs = iris.coord_systems.GeogCS(6371229.0)
+        self.sphere_source = CubeBuilder(self.sphere_crs, (2, 2))._cube
+        self.points = np.array([[10, 10], [10, -10], [-10, -10], [-10, 10]])
+
+    def test_identity_rotation_sphere(self):
+        """Test that rotation to a pole at latitude=90.0, longitude=0.0
+        returns the same points.
+
+        By convention, rotated pole coordinate systems will set the prime
+        meridian rotated 180.0 from the specified longitude. To place the
+        prime meridian at 0.0, we apply a further rotation of 180.0, following
+        rotation to the new pole.
+        """
+
+        crs = iris.coord_systems.RotatedGeogCS(
+            90.0, 0.0, 180.0, ellipsoid=self.sphere_crs
+        )
+        target_lsm = CubeBuilder(crs, (2, 2))._cube
+        expected_points = np.copy(self.points)
+
+        rotated_coords = _transform_coordinates(
+            target_lsm, self.sphere_source, self.points
+        )
+
+        self.assertArrayAlmostEqual(expected_points, rotated_coords)
+
+    def test_longitudinal_rotation(self):
+        """Test rotation to a new pole at latitude=90.0, longitude=90.0."""
+
+        crs = iris.coord_systems.RotatedGeogCS(
+            90.0, 90.0, 180.0, ellipsoid=self.sphere_crs
+        )
+        target_lsm = CubeBuilder(crs, (2, 2))._cube
+        expected_points = np.array(
+            [[-80.0, 10.0], [-80.0, -10.0], [-100.0, -10.0], [-100.0, 10]]
+        )
+
+        rotated_coords = _transform_coordinates(
+            target_lsm, self.sphere_source, self.points
+        )
+
+        self.assertArrayAlmostEqual(expected_points, rotated_coords)
+
+    def test_antimeridian_rotation(self):
+        """Test rotation to a new pole at latitude=90.0, longitude=180.0.
+        In this case, the points are assumed to cross the anti meridian
+        and so are instead mapped to the interval [0, 360.0]."""
+
+        crs = iris.coord_systems.RotatedGeogCS(
+            90.0, 180.0, 180.0, ellipsoid=self.sphere_crs
+        )
+        target_lsm = CubeBuilder(crs, (2, 2))._cube
+        expected_points = np.array(
+            [[190.0, 10.0], [190.0, -10.0], [170.0, -10.0], [170.0, 10]]
+        )
+
+        warning_msg = re.escape(
+            "The transformed points are assumed to cross the antimeridian. "
+            "The longitudinal points will instead be defined in the interval"
+            " [0, 360.0] degrees."
+        )
+        with self.assertWarnsRegex(UserWarning, warning_msg):
+            rotated_coords = _transform_coordinates(
+                target_lsm, self.sphere_source, self.points
+            )
+
+        self.assertArrayAlmostEqual(expected_points, rotated_coords)
+
+    def test_local_target_no_domain_warning(self):
+        """Test that no warning is emitted if the transformed polygon
+        lies inside the target domain."""
+
+        crs = iris.coord_systems.RotatedGeogCS(
+            90.0, 10.0, 180.0, ellipsoid=self.sphere_crs
+        )
+        target_lsm = CubeBuilder(crs, (2, 2), xlim=(-25.0, 5.0), ylim=(-12, 12))._cube
+        expected_points = np.array([[0, 10], [0, -10], [-20, -10.0], [-20, 10]])
+
+        rotated_coords = _transform_coordinates(
+            target_lsm, self.sphere_source, self.points
+        )
+
+        self.assertArrayAlmostEqual(expected_points, rotated_coords)
+
+    def test_local_target_domain_warning(self):
+        """Test that a warning is emitted if the transformed polygon
+        lies at least partially outside of the valid target domain."""
+
+        crs = iris.coord_systems.RotatedGeogCS(
+            90.0, 10.0, 180.0, ellipsoid=self.sphere_crs
+        )
+        target_lsm = CubeBuilder(crs, (2, 2), xlim=(-10.0, 5.0), ylim=(-12, 12))._cube
+        expected_points = np.array([[0, 10], [0, -10], [-20, -10.0], [-20, 10]])
+
+        warning_msg = re.escape(
+            "The transformed points lie outside the target lsm domain."
+        )
+        with self.assertWarnsRegex(UserWarning, warning_msg):
+            rotated_coords = _transform_coordinates(
+                target_lsm, self.sphere_source, self.points
+            )
+
+        self.assertArrayAlmostEqual(expected_points, rotated_coords)
+
+
+class Test__validate_args(ants.tests.TestCase):
+
+    def test_error_raised(self):
+        """Test that an error is raised if only --source is passed."""
+
+        args = argparse.Namespace(
+            json_file="json/path",
+            output="output/path",
+            target_lsm=None,
+            source="source/path",
+        )
+        error_msg = "If --source is passed then --target-lsm must also be given."
+
+        with self.assertRaisesRegex(ValueError, error_msg):
+            _validate_args(args.target_lsm, args.source)
+
+
+class Test__load_cubes(ants.tests.TestCase):
+    def test_load_lsm_called(self):
+        """Test that loading the landsea mask is successfully called."""
+
+        with mock.patch("ants.cli.ancil_create_shapefile.load_cube") as mock_load:
+            _ = _load_cubes("target/path", None)
+
+        mock_load.assert_called_once_with("target/path")
+
+    @mock.patch("ants.cli.ancil_create_shapefile.load_cube")
+    def test_cubebuilder_called(self, *args):
+        """Test source cube created if no path given."""
+
+        with mock.patch(
+            "ants.cli.ancil_create_shapefile.CubeBuilder"
+        ) as mock_build_cube:
+            _ = _load_cubes("target/path", None)
+
+        mock_build_cube.assert_called_once()
+
+    def test_load_cube_called(self, *args):
+        """Test load cube called if source path given."""
+
+        with mock.patch("ants.cli.ancil_create_shapefile.load_cube") as mock_load:
+            _ = _load_cubes("target/path", "source/cube")
+
+        self.assertEqual(2, mock_load.call_count)
+
+
+class Test__transform_if_required(ants.tests.TestCase):
+    def setUp(self):
+        self.points = np.array([[10, 10], [10, -10], [-10, -10], [-10, 10]])
+        self.sphere_crs = iris.coord_systems.GeogCS(6371229.0)
+        self.sphere_source = CubeBuilder(self.sphere_crs, (2, 2))._cube
+
+    def test_non_rotated_returned(self):
+        """Test that the input points are returned if a non-rotated pole is
+        provided as the target_lsm."""
+
+        crs = iris.coord_systems.RotatedGeogCS(90.0, 0.0, ellipsoid=self.sphere_crs)
+        target_lsm = CubeBuilder(crs, (2, 2))._cube
+
+        with self.assertWarns(UserWarning):
+            rotated_points = _transform_if_required(
+                target_lsm, self.sphere_source, self.points
+            )
+
+        self.assertArrayEqual(rotated_points, self.points)
+
+    def test_transform_called(self):
+        """Test that transform coordinates is called with the correct arguments
+        when a valid target_lsm is given."""
+
+        crs = iris.coord_systems.RotatedGeogCS(0.0, 0.0, ellipsoid=self.sphere_crs)
+        target_lsm = CubeBuilder(crs, (2, 2))._cube
+
+        with mock.patch(
+            "ants.cli.ancil_create_shapefile._transform_coordinates"
+        ) as mock_transform:
+            _ = _transform_if_required(target_lsm, self.sphere_source, self.points)
+
+        received_lsm, received_source, received_points = mock_transform.call_args.args
+
+        mock_transform.assert_called_once()
+        self.assertEqual(target_lsm, received_lsm)
+        self.assertEqual(received_source, self.sphere_source)
+        self.assertArrayEqual(received_points, self.points)
+
+
+class Test_ite_transform(ants.tests.TestCase):
+
+    def test_uk_pole_rotation(self):
+        """
+        The following test uses points that define a validity polygon over the UK.
+
+        The input points are specified in an unrotated geodetic coordinate system
+        as (lon, lat) pairs. They are then rotated to a pole at lon=177.5,
+        lat=37.5. Note that PROJ applies an additional 180.0 degree rotation
+        to the specified pole location.
+
+        The points have been obtained by unrotating the coordinates under
+        $UMDIR/ancil/data/shapefiles/ite_ukv_polygon/runme.py, which can
+        also be found as part of the ANTS rose-stem test suite under
+        /data/users/ants/sources/ANTS/developer/core/ancil_create_shapefile/
+
+        Notes
+        -----
+        The transformed points differ from the expected points by 360 degrees in
+        longitude. Cartopy treats these longitudes as equivalent, but operations
+        that perform Cartesian comparisons may not consider the resulting polygons
+        equivalent.
+        """
+
+        points = np.array(
+            [
+                [1.63160953, 51.09703216],
+                [-0.29556461, 50.36332573],
+                [-5.37076475, 49.91322904],
+                [-6.01981579, 50.15822986],
+                [-5.16157337, 53.5014088],
+                [-3.80608343, 53.99550534],
+                [-4.04624456, 54.56480733],
+                [-5.10119435, 54.47567927],
+                [-5.95760936, 55.29939785],
+                [-6.998891, 55.85406809],
+                [-7.98467448, 56.73347598],
+                [-7.65168129, 58.36692001],
+                [-3.71924355, 61.09496898],
+                [0.11346329, 61.0768747],
+                [2.0703864, 52.72809453],
+                [1.7594481, 51.23560201],
+            ]
+        )
+        expected_points = np.array(
+            [
+                [362.594, -1.32876],
+                [361.407, -2.1152],
+                [358.15, -2.55],
+                [357.744, -2.28679],
+                [358.417, 1.03058],
+                [359.232, 1.50245],
+                [359.103, 2.07441],
+                [358.488, 2.00291],
+                [358.03, 2.84656],
+                [357.472, 3.43282],
+                [356.986, 4.34795],
+                [357.286, 5.96374],
+                [359.404, 8.6],
+                [361.278, 8.6],
+                [362.766, 0.315629],
+                [362.666, -1.18577],
+            ]
+        )
+
+        sphere_crs = iris.coord_systems.GeogCS(6371229.0)
+        source = CubeBuilder(sphere_crs, (2, 2))._cube
+        crs = iris.coord_systems.RotatedGeogCS(37.5, 177.5, ellipsoid=sphere_crs)
+        target_lsm = CubeBuilder(crs, (2, 2))._cube
+
+        rotated_points = _transform_if_required(target_lsm, source, points)
+
+        # Adding 360.0 to longitude to compare to the expected points.
+        rotated_points[:, 0] += 360.0
+        self.assertArrayAlmostEqual(rotated_points, expected_points)