We present global linear and nonlinear simulations of ion temperature gradient instabilities based on a fluid formulation, with an adapted version of the JOREK code. These simulations are performed in realistic global tokamak equilibria based on the solution of the Grad-Shafranov equation. Benchmarking of linear growth rates was successfully completed with respect to previously published data. We find two distinct types of eigenstructures, depending on the magnetic shear. For high shear, when the coupling of poloidal harmonics is strong, ballooning-type eigenmodes are formed, which are up-down asymmetric with a finite ballooning angle, θ0. The poloidal harmonics which form the global eigenmode are found to demonstrate a radial shift, being centered well outside of their corresponding rational surface. Stronger diamagnetic effects increase both θ0 and proportionately shift the m harmonics to larger radii (by as much as two rational surfaces). In the low shear regime, the unstable eigenmodes become narrowly localized between neighboring pairs of rational surfaces, and exhibit no up-down asymmetry. Our simulations also show the generation of finite Reynolds stress due to nonlocal/global profile effects. This stress possesses both poloidally symmetric (n = m = 0) and asymmetric (finite-m) components. Turbulent saturation in nonlinear simulations is demonstrated for both shear regimes.