Hydro-infrastructural systems (e.g., flood control dams, stormwater detention basins, and seawalls) are designed to protect the public against the adverse impacts of various hydrologic extremes (e.g., floods, droughts, and storm surges). In their design and safety evaluation, the characteristics of concerned hydrologic extremes affecting the hydrosystem performance often are described by several interrelated random variables—not just one—that need to be considered simultaneously. These multiple random variables, in practical problems, have a mixture of non-normal distributions of which the joint distribution function is difficult to establish. To tackle problems involving multivariate non-normal variables, one frequently adopted approach is to transform non-normal variables from their original domain to multivariate normal space under which a large wealth of established theories can be utilized. This study presents a framework for practical normal transform based on the third-order polynomial in the context of a multivariate setting. Especially, the study focuses on multivariate third-order polynomial normal transform (TPNT) with explicit consideration of sampling errors in sample L-moments and correlation coefficients. For illustration, the modeling framework is applied to establish an at-site rainfall intensity–duration-frequency (IDF) relationship. Annual maximum rainfall data analyzed contain seven durations (1–72 h) with 27 years of useable records. Numerical application shows that the proposed modeling framework can produce reasonable rainfall IDF relationships by simultaneously treating several correlated rainfall data series and is a viable tool in dealing with multivariate data with a mixture of non-normal distributions.