Visualizing Engineering Data with matplotlib

Matplotlib is the cornerstone plotting library for Python-based engineering work. Whether you are analyzing sensor time-series, comparing simulated results against bench measurements, or communicating statistical distributions of manufacturing tolerances, matplotlib provides the primitives needed to produce publication-quality figures. Understanding how to use these primitives correctly — and how to combine them with NumPy and pandas — transforms raw numerical output into clear, defensible visual arguments. This topic covers the full pipeline: drawing basic plots, choosing the right chart type for the data, composing multi-panel figures, annotating and formatting professionally, and overlaying simulation results on measured data for model validation.

Line Plots for Time-Series Data

A line plot is the natural choice whenever the independent axis represents an ordered sequence — most commonly time. The underlying assumption is that adjacent points are meaningfully connected, so drawing a line between them communicates continuity. In matplotlib the workhorse function is ax.plot(x, y) (or equivalently plt.plot(x, y) in the stateful interface). Given arrays of timestamps and corresponding signal values, a single call produces a connected curve:

import numpy as np
import matplotlib.pyplot as plt

t = np.linspace(0, 1.0, 500)          # 500 samples over 1 second
voltage = 3.3 * np.sin(2 * np.pi * 50 * t) + 0.05 * np.random.randn(500)

fig, ax = plt.subplots()
ax.plot(t, voltage, color='steelblue', linewidth=1.2)
ax.set_xlabel('Time (s)')
ax.set_ylabel('Voltage (V)')
ax.set_title('50 Hz Voltage Signal with Noise')
plt.tight_layout()
plt.show()

When multiple signals must appear on the same axes — for instance, three channels from a data acquisition board — each signal gets its own ax.plot() call. Matplotlib automatically cycles through a default color sequence, but for professional figures it is better to specify colors explicitly so the distinction survives black-and-white printing:

ax.plot(t, ch1, color='steelblue',  linewidth=1.2, label='Channel 1')
ax.plot(t, ch2, color='darkorange', linewidth=1.2, label='Channel 2')
ax.plot(t, ch3, color='seagreen',   linewidth=1.2, label='Channel 3')
ax.legend(loc='upper right')

A legend becomes mandatory once more than one series appears on the same axes. Without it a reader cannot determine which curve corresponds to which physical quantity. Pass label='...' to each plot() call and then call ax.legend(); matplotlib collects the labels automatically.

Axis limits deserve explicit attention. By default matplotlib pads the data range by a small margin, which is usually fine, but in engineering contexts you often want to zoom in on a specific transient or steady-state window. Use ax.set_xlim(t_start, t_end) and ax.set_ylim(v_min, v_max). For example, to focus on the first 20 milliseconds of a startup transient:

ax.set_xlim(0, 0.02)
ax.set_ylim(-0.5, 4.0)

Setting limits explicitly also guards against misleading figures where auto-scaling magnifies noise so dramatically that it looks like a large deviation.

Scatter Plots for Correlation Analysis

A scatter plot is the right tool when you want to reveal the relationship between two continuous variables rather than their evolution over time. Each observation becomes a marker positioned at coordinates (x_i, y_i); the absence of connecting lines emphasizes that the ordering of the points is irrelevant. In matplotlib, plt.scatter(x, y) or ax.scatter(x, y) produces this chart:

load    = np.array([10, 20, 35, 50, 65, 80, 100, 120, 140, 160])  # N
deflect = np.array([0.12, 0.25, 0.43, 0.61, 0.80, 0.99, 1.24, 1.49, 1.73, 1.98])  # mm

fig, ax = plt.subplots()
ax.scatter(load, deflect, color='firebrick', s=50, zorder=3, label='Measurements')
ax.set_xlabel('Applied Load (N)')
ax.set_ylabel('Beam Deflection (mm)')

To communicate the general trend, add a best-fit line using NumPy's polynomial fitting routine. np.polyfit(x, y, deg) returns the coefficients of a degree-deg polynomial in least-squares sense; np.polyval(coeffs, x_range) evaluates it:

coeffs = np.polyfit(load, deflect, 1)          # linear fit
x_fit  = np.linspace(load.min(), load.max(), 200)
y_fit  = np.polyval(coeffs, x_fit)

ax.plot(x_fit, y_fit, color='navy', linewidth=1.5,
        linestyle='--', label=f'Linear fit  y={coeffs[0]:.4f}x+{coeffs[1]:.4f}')
ax.legend()

Quantifying the relationship objectively requires a correlation coefficient. NumPy's np.corrcoef(x, y) returns a 2×2 matrix; the off-diagonal element is the Pearson r:

r = np.corrcoef(load, deflect)[0, 1]
print(f'Pearson r = {r:.4f}')

An r value close to ±1 indicates a strong linear relationship; values near 0 suggest no linear association (though non-linear relationships may still exist). With pandas, df['load'].corr(df['deflect']) yields the same result and is convenient when working with DataFrame columns.

When data spans several orders of magnitude — for example, plotting resistivity values from 10⁻⁴ Ω·m to 10⁸ Ω·m — linear axes compress most of the data into an unreadable cluster at one end. Logarithmic scaling fixes this:

ax.set_xscale('log')
ax.set_yscale('log')

On a log-log plot a power-law relationship y = a·x^b appears as a straight line, which is extremely useful for identifying and fitting power-law models common in materials science and fluid mechanics.

Histograms and Distribution Plots

When the goal is to understand how a measured quantity is distributed — not how it varies over time or against another variable — histograms are the standard tool. plt.hist(data, bins=n) partitions the data range into n equal-width bins and draws a bar whose height equals the count of values falling in that bin:

diameters = np.random.normal(loc=25.00, scale=0.05, size=500)  # mm

fig, ax = plt.subplots()
ax.hist(diameters, bins=30, color='slategray', edgecolor='white', linewidth=0.5)
ax.set_xlabel('Shaft Diameter (mm)')
ax.set_ylabel('Count')

For comparison across datasets of different sizes — or for overlaying a theoretical density function — normalize the histogram to a probability density by passing density=True. The area of all bars then sums to 1.0, matching the convention for probability density functions (PDFs):

ax.hist(diameters, bins=30, density=True,
        color='slategray', edgecolor='white', linewidth=0.5, alpha=0.6)

Kernel density estimation (KDE) offers a smooth, continuous alternative to the histogram. It places a small kernel function (typically Gaussian) over each data point and sums them into a single curve. The bandwidth controls the smoothness. Using scipy.stats.gaussian_kde:

from scipy.stats import gaussian_kde

kde = gaussian_kde(diameters, bw_method='scott')
x_range = np.linspace(diameters.min() - 0.1, diameters.max() + 0.1, 400)
ax.plot(x_range, kde(x_range), color='navy', linewidth=2, label='KDE')

KDE is particularly valuable when the sample size is small and histogram bar heights are sensitive to the number and placement of bins. If you are using pandas, series.plot.kde() provides a convenient wrapper.

A common engineering task is checking whether measurement noise is Gaussian — a fundamental assumption in many statistical tests and Kalman filter designs. Overlay a normal distribution curve parameterized by the sample mean and standard deviation:

from scipy.stats import norm

mu, sigma = diameters.mean(), diameters.std()
ax.plot(x_range, norm.pdf(x_range, mu, sigma),
        color='firebrick', linewidth=2, linestyle='--', label='Normal fit')
ax.legend()

If the histogram bars follow the red dashed curve closely, Gaussian noise is a reasonable assumption. Systematic departures — heavy tails, skewness, multimodal humps — indicate that the noise model may need revision, which has direct consequences for any downstream filtering or uncertainty analysis.

Subplots and Multi-Panel Figures

Engineering reports routinely need to show several related plots side-by-side or stacked vertically: for instance, a voltage signal above and its frequency spectrum below, or three different test conditions across a row. Matplotlib's plt.subplots(nrows, ncols) creates a grid of Axes objects sharing a single Figure:

fig, axes = plt.subplots(nrows=3, ncols=1, figsize=(10, 8))

axes[0].plot(t, position, color='steelblue')
axes[0].set_ylabel('Position (m)')

axes[1].plot(t, velocity, color='darkorange')
axes[1].set_ylabel('Velocity (m/s)')

axes[2].plot(t, acceleration, color='seagreen')
axes[2].set_ylabel('Acceleration (m/s²)')
axes[2].set_xlabel('Time (s)')

When the Axes grid is 2-D (multiple rows and columns), the returned axes array is indexed as axes[row, col]. For a 2×2 grid you would write axes[0, 0], axes[0, 1], and so on.

The sharex and sharey parameters are especially valuable when stacking time-series panels. Setting sharex=True links the x-axis limits and ticks across all subplots so that zooming or panning one panel updates all of them simultaneously, and tick labels are suppressed on interior axes to avoid clutter:

fig, axes = plt.subplots(3, 1, sharex=True, figsize=(10, 8))

sharey=True is useful when comparing distributions or signals that should be evaluated on the same scale. Use it when the y-axis meaning is identical across panels; avoid it when each panel has a different physical unit.

Figure size must be scaled deliberately with the number of panels. The default figsize=(6.4, 4.8) inches is fine for a single plot but becomes cramped for a 3×1 or 2×3 grid. A reasonable heuristic is approximately 4 inches of height per row and 5 inches of width per column:

fig, axes = plt.subplots(2, 3, figsize=(15, 8))

Always call plt.tight_layout() (or the newer fig.tight_layout()) after placing all content. It automatically adjusts subplot spacing to prevent axis labels from overlapping adjacent panels — a common problem that makes figures unpresentable without this step.

Formatting and Annotating Plots

A figure without descriptive labels is incomplete regardless of how accurate the data is. At minimum every axis needs a label with its physical quantity and unit, and every multi-series plot needs a legend. The standard pattern using the object-oriented interface is:

ax.set_xlabel('Time (s)', fontsize=12)
ax.set_ylabel('Force (N)', fontsize=12)
ax.set_title('Step Response — Hydraulic Actuator', fontsize=13)

For annotations — drawing a reader's eye to a specific feature such as a peak value, a failure threshold crossing, or a resonance frequency — ax.annotate() provides both a text label and an optional arrow:

ax.annotate(
    'Peak: 4.82 N',
    xy=(0.045, 4.82),           # point to annotate
    xytext=(0.065, 4.60),       # text position
    arrowprops=dict(arrowstyle='->', color='black'),
    fontsize=10
)

The xy argument is the data coordinate of the feature; xytext is where the label text should appear. An arrow connecting the two makes the association unambiguous. This is far more informative than just noting a peak value in a figure caption.

Exporting figures for reports or journal submissions uses plt.savefig(). The file format is inferred from the extension. For raster images such as PNG, set dpi appropriately; 150 dpi is adequate for screen, 300 dpi for print. For vector formats such as PDF or SVG, resolution is irrelevant because the geometry is stored analytically:

plt.savefig('actuator_response.pdf', bbox_inches='tight')
plt.savefig('actuator_response.png', dpi=300, bbox_inches='tight')

bbox_inches='tight' trims excess whitespace around the figure, which is nearly always desirable.

Many engineering organizations and journals enforce style guidelines specifying font sizes (commonly 10–12 pt for axis labels, 8–10 pt for tick labels), minimum line widths (often 1.5 pt for data lines), and restricted color palettes chosen for color-blind accessibility. Rather than setting these properties figure by figure, configure them once via matplotlib's rc system:

import matplotlib as mpl

mpl.rcParams['font.size']       = 11
mpl.rcParams['axes.linewidth']  = 1.2
mpl.rcParams['lines.linewidth'] = 1.5
mpl.rcParams['xtick.major.width'] = 1.2
mpl.rcParams['ytick.major.width'] = 1.2

Place these settings at the top of a script or in a shared configuration module and every subsequent figure in the project will inherit them automatically.

Visualizing Simulation Results Alongside Measured Data

Model validation is one of the most critical activities in engineering, and its visual representation must be unambiguous. The universally accepted convention is to plot measured data as a solid line (or markers) and simulated results as a dashed or dotted line. This signals immediately to any reader which curve was recorded from hardware and which was produced by a computational model:

ax.plot(t_meas, y_meas, color='steelblue',  linewidth=1.5,
        linestyle='solid',  label='Measured')
ax.plot(t_sim,  y_sim,  color='firebrick', linewidth=1.5,
        linestyle='dashed', label='Simulated')
ax.legend()

Before overlaying two signals it is essential to verify that their time bases are aligned. A simulated output typically starts at t=0, but a measured record may carry an absolute timestamp offset or may start recording mid-transient. Failure to align time axes produces apparent phase shifts that can be mistakenly attributed to modeling error. Correct the offset explicitly:

t_meas_aligned = t_meas - t_meas[0]   # shift so measurement starts at t=0

Additionally, if the simulation and measurement have different sample rates, one must be interpolated onto the time grid of the other before computing a difference. NumPy's np.interp() is the simplest tool for this:

y_sim_resampled = np.interp(t_meas_aligned, t_sim, y_sim)

With time alignment confirmed and both signals on the same grid, compute and plot the residual — the point-by-point difference between the simulation output and the measurement. A residual plot is arguably the most diagnostic view in model validation:

residual = y_meas - y_sim_resampled

fig, axes = plt.subplots(2, 1, sharex=True, figsize=(10, 6))

axes[0].plot(t_meas_aligned, y_meas,          color='steelblue',
             linestyle='solid',  label='Measured')
axes[0].plot(t_meas_aligned, y_sim_resampled, color='firebrick',
             linestyle='dashed', label='Simulated')
axes[0].set_ylabel('Output Signal')
axes[0].legend()

axes[1].plot(t_meas_aligned, residual, color='dimgray', linewidth=1.0)
axes[1].axhline(0, color='black', linewidth=0.8, linestyle='--')
axes[1].set_ylabel('Residual (Measured − Simulated)')
axes[1].set_xlabel('Time (s)')

plt.tight_layout()
plt.savefig('validation_plot.pdf', bbox_inches='tight')

Interpreting the residual panel is straightforward: if the residuals are small in magnitude and randomly distributed around zero, the model captures the dominant dynamics well. A systematic offset suggests a steady-state error or gain mismatch. Periodic structure in the residuals indicates unmodeled dynamics at a specific frequency. A drift that grows with time points to an integrating error such as an incorrect initial condition or an unmodeled slow process.

Quantitative summary statistics of the residual — root mean square error (RMSE), maximum absolute error, or mean absolute percentage error — should be computed and displayed, either in the figure title or as an annotation:

rmse = np.sqrt(np.mean(residual**2))
axes[1].set_title(f'RMSE = {rmse:.4f}  |  Max |error| = {np.max(np.abs(residual)):.4f}')

A legend identifying each series is not optional in a validation figure. Anyone reading the report must be able to immediately identify measured versus simulated without referring to surrounding text. Where multiple simulation runs are shown (for example, with different model parameter values), each run should have a distinct line style or color, and all must appear in the legend with descriptive labels.

Combining all of these elements — correct chart type for the data, explicit axis limits and labels, normalized histograms for distribution analysis, shared axes in multi-panel layouts, consistent style settings, and properly formatted validation plots — enables engineers to produce figures that communicate findings clearly, withstand scrutiny, and satisfy the standards of professional technical reporting.