diff --git a/autotest/gdrivers/pds4.py b/autotest/gdrivers/pds4.py
index 987eb65fa006..9439953e7790 100755
--- a/autotest/gdrivers/pds4.py
+++ b/autotest/gdrivers/pds4.py
@@ -1800,3 +1800,24 @@ def test_pds4_oblique_cylindrical_write():
     check_pds4_oblique_cylindrical(filename)
 
     gdal.GetDriverByName("PDS4").Delete(filename)
+
+
+###############################################################################
+
+
+def test_pds4_read_right_to_left(tmp_path):
+
+    numpy = pytest.importorskip("numpy")
+    pytest.importorskip("osgeo.gdal_array")
+
+    tmp_filename = str(tmp_path / "tmp.xml")
+    ref_ds = gdal.Open("data/byte.tif")
+    gdal.Translate(tmp_filename, ref_ds, format="PDS4")
+    xml_content = open(tmp_filename, "rt").read()
+    # Generate a fake Right to Left oriented image
+    open(tmp_filename, "wt").write(
+        xml_content.replace("Left to Right", "Right to Left")
+    )
+    ds = gdal.Open(tmp_filename)
+    # Test that we flip the image along the horizontal axis
+    assert numpy.all(ds.ReadAsArray()[::, ::-1] == ref_ds.ReadAsArray())
diff --git a/frmts/pds/pds4dataset.cpp b/frmts/pds/pds4dataset.cpp
index 559868b77975..f91481f49fce 100644
--- a/frmts/pds/pds4dataset.cpp
+++ b/frmts/pds/pds4dataset.cpp
@@ -1584,6 +1584,13 @@ PDS4Dataset *PDS4Dataset::OpenInternal(GDALOpenInfo *poOpenInfo)
         "");
     const bool bBottomToTop = EQUAL(pszVertDir, "Bottom to Top");
 
+    const char *pszHorizDir = CPLGetXMLValue(
+        psProduct,
+        "Observation_Area.Discipline_Area.Display_Settings.Display_Direction."
+        "horizontal_display_direction",
+        "");
+    const bool bRightToLeft = EQUAL(pszHorizDir, "Right to Left");
+
     auto poDS = std::make_unique<PDS4Dataset>();
     poDS->m_osXMLFilename = osXMLFilename;
     poDS->eAccess = eAccess;
@@ -1999,14 +2006,21 @@ PDS4Dataset *PDS4Dataset::OpenInternal(GDALOpenInfo *poOpenInfo)
 
             for (int i = 0; i < l_nBands; i++)
             {
+                vsi_l_offset nThisBandOffset = nOffset + nBandOffset * i;
+                if (bBottomToTop)
+                {
+                    nThisBandOffset +=
+                        static_cast<vsi_l_offset>(nLines - 1) * nLineOffset;
+                }
+                if (bRightToLeft)
+                {
+                    nThisBandOffset +=
+                        static_cast<vsi_l_offset>(nSamples - 1) * nPixelOffset;
+                }
                 auto poBand = std::make_unique<PDS4RawRasterBand>(
-                    poDS.get(), i + 1, poDS->m_fpImage,
-                    (bBottomToTop) ? nOffset + nBandOffset * i +
-                                         static_cast<vsi_l_offset>(nLines - 1) *
-                                             nLineOffset
-                                   : nOffset + nBandOffset * i,
-                    nPixelOffset, (bBottomToTop) ? -nLineOffset : nLineOffset,
-                    eDT,
+                    poDS.get(), i + 1, poDS->m_fpImage, nThisBandOffset,
+                    bRightToLeft ? -nPixelOffset : nPixelOffset,
+                    bBottomToTop ? -nLineOffset : nLineOffset, eDT,
                     bLSBOrder ? RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN
                               : RawRasterBand::ByteOrder::ORDER_BIG_ENDIAN);
                 if (!poBand->IsValid())