Highest Density Regions (HDR) misunderstanding

I think the current implementation only works for single mode marginal, as you observed for bimodal where the modes are far away it will gives incorrect answer.
I remember seeing a version that works for multiple mode from @aloctavodia, maybe he can provide a bit more information.