Repository: pydap/pydap

Issue:
slices defined in a URL's Constraint expression are not being propagated correctly when retrieving remote data
<!-- Thank you very much for your input! -->

### Description
Slices in the form of `[100:200]` and so on, defined in constraint expressions, affect arrays in both `DAP2` and `DAP4`. Broadly, internally when there is a slice, pydap always assumes it starts from the 0-index. This causes an error in `DAP2`, but there is no error in DAP4.

## Minimal example that reproduces the bug, or a new desired feature
```python
import matplotlib.pyplot as plt
from pydap.client import open_url
import numpy as np

url = "http://test.opendap.org:8080/opendap/netcdf/examples/200803061600_HFRadar_USEGC_6km_rtv_SIO.nc"
```

### DAP4

```python
## get all longitude values
ds = open_url(url+"?dap4.ce=/lon")
lon = ds['lon'][:]
np.shape(lon)
>>> (701,)

## subset in DAP4 - get 100 points away from edge data values
ds_ce = open_url(url+'?dap4.ce=/lon[100:1:199]')
# URL of array is correct though!
ds_ce['lon'] 
>>> <BaseType with data Dap4BaseProxy('http://test.opendap.org:8080/opendap/netcdf/examples/200803061600_HFRadar_USEGC_6km_rtv_SIO.nc?dap4.ce=/lon[100:1:199]', 'lon', dtype('>f4'), (100,), (slice(None, None, None),))>
lon_ce1 = ds_ce['lon'][:].data
```
Even though the URL and the shape are correct,`slice(None, None, None)` contains no information about the slice in the CE.

The following comparisons FAILS but should be True:

```python
lon[100] == lon_ce1[0]
>>> np.False
assert (np.equal(lon[100:200], lon_ce1)).all()
>>> ---------------------------------------------------------------------------
AssertionError                            Traceback (most recent call last)
Cell In[41], line 1
----> 1 assert (np.equal(lon[100:200], lon_ce1)).all()

AssertionError:
```

### DAP2
I get something slightly different, but the issue remains the same. Starting with the constraint expression:
```python
ds_ce = open_url(url+'?lon[100:1:199]', protocol='dap2')
# The following shows that there is a problem with the pydap-defined slice
ds_ce['lon']
>>><BaseType with data BaseProxy('http://test.opendap.org:8080/opendap/netcdf/examples/200803061600_HFRadar_USEGC_6km_rtv_SIO.nc', 'lon', dtype('>f4'), (100,), (slice(100, 100, 1),))>
```

The slice is defined as `slice(100, 100, 1)`. Even attempting to download data results an error:
```python
ds_ce['lon'][:].data
>>> HTTPError: 400 Bad Request
Error { 
    code = 400;
    message = "libdap error transmitting DataDDS: The start value of an array index is past the stop value.";
}
```


### proposed solution?

I think I know what is going on. THis is only an issue when the slicing defined in the Constraint Expression does NOT start from `0`. This is , somehow `pydap` fails to identify the start of the Slice. `pydap` correctly identifies the shape of the slice.

Also, repeating the exact steps described above, but with the following URL both DAP2 and DAP4 WORK!

```python
ds_ce = open_url(url+'?dap4.ce=/lon[0:1:99]')
```



Discussion:
Mikejmnez: Added testing for catching scenarios with inconsistent specification of protocol, and syntax of CEs, things pass now. I will merge :)

Patch:
diff --git a/src/pydap/handlers/dap.py b/src/pydap/handlers/dap.py
index dde0500..224e4e3 100644
--- a/src/pydap/handlers/dap.py
+++ b/src/pydap/handlers/dap.py
@@ -85,7 +85,9 @@ class DAPHandler(pydap.handlers.lib.BaseHandler):
 
         self.protocol = self.determine_protocol(protocol)
 
-        self.projection, self.selection = pydap.parsers.parse_ce(self.query)
+        self.projection, self.selection = pydap.parsers.parse_ce(
+            self.query, self.protocol
+        )
         arg = (
             self.scheme,
             self.netloc,
@@ -225,6 +227,17 @@ class DAPHandler(pydap.handlers.lib.BaseHandler):
                 timeout=self.timeout,
             )
 
+        # apply projections to BaseType only
+        # CE for sequences and structs
+        # are not ready (see https://github.com/pydap/pydap/issues/314)
+        for var in self.projection:
+            target = self.dataset
+            while var:
+                token, index = var.pop(0)
+                target = target[token]
+                if isinstance(target, pydap.model.BaseType):
+                    target.data.slice = fix_slice(index, target.shape)
+
     def add_dap2_proxies(self):
         # now add data proxies
         for var in walk(self.dataset, pydap.model.BaseType):
diff --git a/src/pydap/lib.py b/src/pydap/lib.py
index 589dbfb..2fbd34c 100644
--- a/src/pydap/lib.py
+++ b/src/pydap/lib.py
@@ -192,7 +192,7 @@ def fix_slice(slice_, shape):
 
     This is based on this document:
 
-        http://docs.scipy.org/doc/numpy/reference/arrays.indexing.html
+        https://numpy.org/doc/stable/user/basics.indexing.html#slicing-and-striding
 
     """
     # convert `slice_` to a tuple
@@ -211,10 +211,10 @@ def fix_slice(slice_, shape):
     slice_ = tuple(out) + (slice(None),) * expand
 
     out = []
-    for s, n in zip(slice_, shape):
+    for s, N in zip(slice_, shape):
         if isinstance(s, int):
             if s < 0:
-                s += n
+                s += N
             out.append(s)
         else:
             k = s.step or 1
@@ -223,13 +223,13 @@ def fix_slice(slice_, shape):
             if i is None:
                 i = 0
             elif i < 0:
-                i += n
+                i += N
 
             j = s.stop
-            if j is None or j > n:
-                j = n
+            if j is None or j >= (N + i):
+                j = N + i
             elif j < 0:
-                j += n
+                j += N
 
             out.append(slice(i, j, k))
 
diff --git a/src/pydap/parsers/__init__.py b/src/pydap/parsers/__init__.py
index 3efd71d..dd501ee 100644
--- a/src/pydap/parsers/__init__.py
+++ b/src/pydap/parsers/__init__.py
@@ -13,7 +13,7 @@ from ..exceptions import ConstraintExpressionError
 from ..lib import get_var, unquote
 
 
-def parse_projection(input):
+def parse_projection(input, protocol="dap2"):
     """Split a projection into items.
 
     The function takes into account server-side functions, and parse slices
@@ -22,6 +22,10 @@ def parse_projection(input):
     Returns a list of names and slices.
 
     """
+    if protocol == "dap2":
+        key = ","
+    elif protocol == "dap4":
+        key = ";"
 
     def tokenize(input):
         start = pos = count = 0
@@ -30,7 +34,7 @@ def parse_projection(input):
                 count += 1
             elif char == ")":
                 count -= 1
-            elif char == "," and count == 0:
+            elif char == key and count == 0:
                 yield input[start:pos]
                 start = pos + 1
             pos += 1
@@ -110,19 +114,30 @@ def parse_ce(query_string, protocol="dap2"):
     """
     if protocol == "dap2":
         key = "&"
-        tokens = [token for token in unquote(query_string).split(key) if token]
-        if not tokens:
-            projection = []
-            selection = []
-        elif re.search("<=|>=|!=|=~|>|<|=", tokens[0]):
-            projection = []
-            selection = tokens
-        else:
-            projection = parse_projection(tokens[0])
-            selection = tokens[1:]
-    else:
+        if len(query_string) > 0 and query_string[:8] == "dap4.ce=":
+            raise ConstraintExpressionError(
+                "The Constraint Expression %s does not follow the DAP2 "
+                "model specification" % query_string,
+            )
+    elif protocol == "dap4":
+        key = "|"
+        # remove `dap4.ce=` from query string if there is a query
+        if len(query_string) > 0 and query_string[:8] != "dap4.ce=":
+            raise ConstraintExpressionError(
+                "The Constraint Expression %s does not follow the DAP4"
+                "model specification " % query_string,
+            )
+        query_string = query_string[8:]
+    tokens = [token for token in unquote(query_string).split(key) if token]
+    if not tokens:
         projection = []
         selection = []
+    elif re.search("<=|>=|!=|=~|>|<|=", tokens[0]):
+        projection = []
+        selection = tokens
+    else:
+        projection = parse_projection(tokens[0], protocol)
+        selection = tokens[1:]
 
     return projection, selection
 


Test patch:
diff --git a/src/pydap/tests/test_client.py b/src/pydap/tests/test_client.py
index f3a9abf..b423026 100644
--- a/src/pydap/tests/test_client.py
+++ b/src/pydap/tests/test_client.py
@@ -45,7 +45,7 @@ def test_open_url_dap4(remote_url):
     assert data_dap4["s33"][:].data == data_original["s33"][0, 0].data
 
     # subset of vars by indexes
-    var1 = "/s33[0:][0:1:2];"
+    var1 = "/s33[0:1:2][0:1:2];"
     var2 = "/br34[0:1:1][0:1:2][0:1:3];"
     var3 = "/s113[0:1:0][0:1:0][0:1:2]"
     Vars = [var1, var2, var3]
@@ -56,6 +56,15 @@ def test_open_url_dap4(remote_url):
     assert len(dataset.keys()) == len(Vars)
 
 
+def test_open_url_dap4_shape():
+    url = "http://test.opendap.org:8080/opendap/"
+    filename = "netcdf/examples/200803061600_HFRadar_USEGC_6km_rtv_SIO.nc"
+    CE = "?dap4.ce=/lon[100:1:199]"
+    ds_ce = open_url(url + filename + CE)
+    data = ds_ce["lon"][:]
+    assert data.shape == (100,)
+
+
 def test_open_url_seqCE(remote_url):
     seq_url = remote_url + "ff/gsodock.dat"
     data_original = open_url(seq_url)
diff --git a/src/pydap/tests/test_handlers_dap.py b/src/pydap/tests/test_handlers_dap.py
index 3b56fcd..ed2372a 100644
--- a/src/pydap/tests/test_handlers_dap.py
+++ b/src/pydap/tests/test_handlers_dap.py
@@ -174,7 +174,7 @@ class TestDapHandler(unittest.TestCase):
         dataset = DAPHandler("http://localhost:8001/?x[1:1:2]", self.app1).dataset
 
         self.assertEqual(dataset.x.data.shape, (2,))
-        self.assertEqual(dataset.x.data.slice, (slice(1, 2, 1),))
+        self.assertEqual(dataset.x.data.slice, (slice(1, 3, 1),))
 
     def test_grid_array_with_projection(self):
         """Test that a grid array can be properly pre sliced."""
diff --git a/src/pydap/tests/test_lib.py b/src/pydap/tests/test_lib.py
index c205993..75acef8 100644
--- a/src/pydap/tests/test_lib.py
+++ b/src/pydap/tests/test_lib.py
@@ -135,6 +135,12 @@ class TestFixSlice(unittest.TestCase):
         self.assertEqual(slice2, (slice(2, 8, 1),))
         np.testing.assert_array_equal(x[slice1], x[slice2])
 
+    def test_non_zero_slicestart(self):
+        x = np.arange(0, 100)
+        slice1 = slice(100, 200)
+        slice2 = fix_slice(slice1, x.shape)
+        self.assertEqual(slice2, (slice(100, 200, 1),))
+
 
 class TestCombineSlices(unittest.TestCase):
     """Test the ``combine_slices`` function."""
diff --git a/src/pydap/tests/test_parsers.py b/src/pydap/tests/test_parsers.py
index 79119e2..df25abe 100644
--- a/src/pydap/tests/test_parsers.py
+++ b/src/pydap/tests/test_parsers.py
@@ -123,6 +123,19 @@ class TestParseCe(unittest.TestCase):
         )
         self.assertEqual(selection, [])
 
+    def test_dap4_ce(self):
+        protocol = "dap4"
+        projection, sel = parse_ce("dap4.ce=/lon[100:1:199]", protocol)
+        self.assertEqual(projection, [[("/lon", (slice(100, 200, 1),))]])
+        self.assertEqual(sel, [])
+
+    def test_inconsistent_CE_protocol(self):
+        query1 = "dap4.ce=/lon[100:1:199]"
+        query2 = "lon[100:1:199]"
+        with self.assertRaises(ConstraintExpressionError):
+            parse_ce(query1, "dap2")
+            parse_ce(query2, "dap4")
+
 
 class TestParseHyperslab(unittest.TestCase):
     """Test hyperslab parser."""