--- a
+++ b/CLI/IntervertebralAnalysis/IntervertebralAnalysis.py
@@ -0,0 +1,92 @@
+#!/usr/bin/env python-real
+
+import sys
+import os
+from datetime import date
+
+import numpy as np
+
+sys.path.append(os.path.dirname(os.path.dirname(os.path.abspath(__file__))))
+from MusculoskeletalAnalysisCLITools import (
+   crop,
+   readImg,
+   readMask,
+   width,
+   writeReport,
+)
+
+
+OUTPUT_FIELDS = [
+    ("Date Analysis Performed", "The current date"),
+    ("Input Volume", "Name of the input volume"),
+    ("Disc Volume (mm^3)", "The volume of the disc"),
+    ("Nucleus Pulposus Volume (mm^3)", "The volume of the NP"),
+    ("Volume Ratio", "The ratio of whole disc volume to NP volume"),
+    ("Annulus Fibrosus Width (mm)", "The width of the AF, calculated by using rotating calipers algorithm on each slice and finding the maximum width"),
+    ("Nucleus Pulposus Width (mm)", "The width of the NP, calculated using the same method as AF width"),
+    ("Disc Height (mm)", "The height of the disc at its center"),
+    ("Disc Height Ratio", "The ratio of disc height to disc width"),
+    ("Voxel Dimension (mm)", "The side length of one voxel, measured in milimeters"),
+]
+
+
+def main(inputImg, inputMask1, inputMask2, voxSize, name, output):
+    imgData = readImg(inputImg)
+    (_, maskData1) = readMask(inputMask1)
+    (_, maskData2) = readMask(inputMask2)
+    # Set the np mask to the smaller one
+    if np.count_nonzero(maskData1) > np.count_nonzero(maskData2):
+        (maskData, npData, imgData) = crop(maskData1, maskData2, imgData)
+    else:
+        (maskData, npData, imgData) = crop(maskData2, maskData1, imgData)
+    (maskData, npData, imgData) = crop(maskData, npData, imgData)
+    if np.count_nonzero(maskData) == 0 or np.count_nonzero(npData) == 0:
+         raise Exception("Segmentation mask is empty.")
+
+    volume = np.count_nonzero(maskData) * voxSize**3
+    npVolume = np.count_nonzero(npData) * voxSize**3
+    vr = npVolume/volume
+
+    print("""<filter-progress>{}</filter-progress>""".format(.33))
+    sys.stdout.flush()
+
+    afWidth= width(maskData) * voxSize
+    npWidth = width(npData) * voxSize
+
+    # Find the center
+    xcent = int(np.rint(np.mean(np.nonzero(maskData)[0])))
+    ycent = int(np.rint(np.mean(np.nonzero(maskData[1]))))
+    center = np.nonzero(maskData[xcent,ycent,:])
+    height = (np.max(center)-np.min(center)) * voxSize
+
+    hw = height/afWidth
+
+    print("""<filter-progress>{}</filter-progress>""".format(.66))
+    sys.stdout.flush()
+
+    fPath = os.path.join(output, "intervertebral.txt")
+    header = [field[0] for field in OUTPUT_FIELDS]
+    data = [
+        date.today(),
+        name,
+        volume,
+        npVolume,
+        vr,
+        afWidth,
+        npWidth,
+        height,
+        hw,
+        voxSize
+    ]
+
+    writeReport(fPath, header, data)
+
+
+
+if __name__ == "__main__":
+    if len(sys.argv) < 7:
+        print(sys.argv)
+        print("Usage: Intervertebral Analysis <input> <mask1> <mask2> <voxelSize> <name> <output>")
+        sys.exit(1)
+    main(sys.argv[1], sys.argv[2], sys.argv[3], float(sys.argv[4]), str(sys.argv[5]), str(sys.argv[6]))
+