Masking netCDF with np.where

Question

I am not too experienced with programming and completely new to working with netCDF files and xarray.

I am trying to apply a mask to some precipitation data to split it into 7 regions. Here is what the mask looks like: enter image description here

The precipitation data has: precipitation in mm coordinates of time, x and y

The mask has: values 0-7, coordinates of x and y on the exact same grid as the precip data

I have already done this on a different mask which had the same structure using

prp1950_masked = prp1950.where(ice_mask==3,drop=True)

where prp1950 is the precipitation data and ice_mask is the ice mask with values 0-3. But when I try to do it for this regional mask, I get:

ValueError: Cannot apply_along_axis when any iteration dimensions are 0

After a whole day of trying, my best guess is that it has something to do with the fact that in prp1950 and ice_mask, x appears first and then y, whereas in the regional mask regs y appears first in the list of coordinates. I have already tried various ways of transposing regs but nothing has worked. See images below.

enter image description here

enter image description here

enter image description here

I would be so very grateful for any solutions to this/analyses of what could be going wrong! (Bearing in mind that most of the code will look like absolute gibberish to me)