We deal with multiple image warping, which computes deformation fields between an image and a collection of images, as an extension of variational image registration. Using multiple image warping, we develop a variational method for the computation of average images of biological organs in three-dimensional Euclidean space. The average shape of three-dimensional biological organs is an essential feature to discriminate abnormal organs from normal organs. There are two kinds of volumetric image sets in medical image analysis. The first one is a collection of static volumetric data of an organ and/or organs. The other is a collection of temporal volumetric data of an organ and/or organs. A collection of temporal volumetric beating hearts is an example of temporal volumetric data. For spatiotemporal volumetric data, we can compute (1) the temporal average, which is the average of a heart during a cycle, (2) the frame average, which is the average of hearts at a frame, and (3) the temporal average of frame averages.