Thank you, @ptrblck. I read this post. In my view, the solution (inclusion-exclusion version) provided in this post works because the sum of the mask_region is wanted. However, in my question, I want to obtain the raw values in the mask region. It seems that the inclusion-exclusion version is not very suitable in this case. I am not very sure… Can you give me some tips in using inclusion-exclusion method?