Calculate harmonics using DFT from real points
The DFT will express your time sequence as a weighted some of complex exponential (basically a set of orthogonal functions, if that's a helpful concept for you). DC is simply a special case of a complex exponential with the frequency being 0. 

Let's assume a sample rate of 1440 Hz (to make the math simple). Then the DFT coefficients mean the following:

X[0]: amplitude at 0 Hz
X[1]: amplitude &amp; phase at 10 Hz
X[2]: amplitude &amp; phase at 20 Hz
...
X[71]: amplitude and phase at 710 Hz
X[72]: amplitude at 720 Hz (Nyquist frequency)


Nyquist is a little bit of an odd-ball. It's a real number so the phase is 0. The value should also be small or close to 0, otherwise you have potentially aliasing.

Since your input is real, you have complex conjugate symmetry as follows

X[73] = X[71]';
X[74] = X[70]';
....
X[143] = X[1]';


Just a quick comment: (X[0] + X[72)/2 is NOT the DC value. X[0] is.