Systems and methods configured to create contour maps of geospatial variables based on hydrometeorological data associated with the variable are described herein. The systems and methods advantageously use simultaneous multi-event universal kriging for spatio-temporal data exploration, analysis and interpolation with the objective of creating a series of related maps, where each map corresponds to a specific sampling event, but wherein some features exhibit spatial relationships persisting over time. In one particular example, water level maps are prepared using the methods, which has the flexibility to allow the conditioning of trend coefficients based on any arbitrary subsets of sample data, and thereby provides a physically based and deterministic rather than wholly-stochastic basis for depicting hydrometeorological data correlations in space and time. The simultaneous MEUK method borrows strength from events possessing a large number of samples to map events possessing fewer data, and allows for implementing “wheel and axle” data analysis.