You ask about "the HW model and/or more general ARIMA models". Note that ARIMA and smoothing methods are not in general equivalent or subsets of each other.
For state space models like ETS, which generalize exponential smoothing, you can calculate PIs from predictive densities. In simple cases, you can actually calculate this parametrically, but more often, you will simulate. Take a look at Hyndman et al., Forecasting with Exponential Smoothing, Chapter 6.
The article by Yar & Chatfield (1990, IJF) that Bruna w links to (+1) covers the specific Holt-Winters case (and for completeness, also simpler Single and Double Exponential Smoothing).
For ARIMA models, similar approaches are possible: simple models allow theoretical computations, and sometimes you need to simulate. Brockwell & Davis, Time Series: Theory and Methods (§5.4) gives a few hints, but other theoretical books on ARIMA like Box, Jenkins, Reinsel & Ljung Time Series Analysis: Forecasting and Control (which I don't have at hand right now) may go more deeply.
In any case, prediction intervals and predictive densities should ideally account for multiple sources of uncertainty:
- Model selection uncertainty
- Parameter estimation uncertainty
- Future innovation uncertainty
Few if any approaches cover all three. Model selection uncertainty in particular is often ignored, and future innovations are often forecasted by pretending we know the parameters precisely, while parameter estimation uncertainty is sometimes included (especially in decent simulation-based approaches) and sometimes ignored.