diff --git a/qiskit/quantum_info/analysis/distance.py b/qiskit/quantum_info/analysis/distance.py
index 761dd0440468..8d282a20e7dc 100644
--- a/qiskit/quantum_info/analysis/distance.py
+++ b/qiskit/quantum_info/analysis/distance.py
@@ -14,13 +14,55 @@
 import numpy as np
 
 
+def _hellinger_distance(dist_p, dist_q):
+    """Computes the Hellinger distance between
+    two counts distributions.
+
+    Parameters:
+        dist_p (dict): First dict of counts.
+        dist_q (dict): Second dict of counts.
+
+    Returns:
+        float: Distance
+
+    References:
+        `Hellinger Distance @ wikipedia <https://en.wikipedia.org/wiki/Hellinger_distance>`_
+    """
+    p_sum = sum(dist_p.values())
+    q_sum = sum(dist_q.values())
+
+    p_normed = {}
+    for key, val in dist_p.items():
+        p_normed[key] = val/p_sum
+
+    q_normed = {}
+    for key, val in dist_q.items():
+        q_normed[key] = val/q_sum
+
+    total = 0
+    for key, val in p_normed.items():
+        if key in q_normed.keys():
+            total += (np.sqrt(val) - np.sqrt(q_normed[key]))**2
+            del q_normed[key]
+        else:
+            total += val
+    total += sum(q_normed.values())
+
+    dist = np.sqrt(total)/np.sqrt(2)
+
+    return dist
+
+
 def hellinger_fidelity(dist_p, dist_q):
     """Computes the Hellinger fidelity between
     two counts distributions.
 
-    The fidelity is defined as 1-H where H is the
-    Hellinger distance.  This value is bounded
-    in the range [0, 1].
+    The fidelity is defined as :math:`\\left(1-H^{2}\\right)^{2}` where H is the
+    Hellinger distance.  This value is bounded in the range [0, 1].
+
+    This is equivalent to the standard classical fidelity
+    :math:`F(Q,P)=\\left(\\sum_{i}\\sqrt{p_{i}q_{i}}\\right)^{2}` that in turn
+    is equal to the quantum state fidelity for diagonal density matrices.
 
     Parameters:
         dist_p (dict): First dict of counts.
@@ -49,27 +91,10 @@ def hellinger_fidelity(dist_p, dist_q):
             res2 = execute(qc, sim).result()
 
             hellinger_fidelity(res1.get_counts(), res2.get_counts())
-    """
-    p_sum = sum(dist_p.values())
-    q_sum = sum(dist_q.values())
-
-    p_normed = {}
-    for key, val in dist_p.items():
-        p_normed[key] = val/p_sum
-
-    q_normed = {}
-    for key, val in dist_q.items():
-        q_normed[key] = val/q_sum
-
-    total = 0
-    for key, val in p_normed.items():
-        if key in q_normed.keys():
-            total += (np.sqrt(val) - np.sqrt(q_normed[key]))**2
-            del q_normed[key]
-        else:
-            total += val
-    total += sum(q_normed.values())
 
-    dist = np.sqrt(total)/np.sqrt(2)
-
-    return 1-dist
+    References:
+        `Quantum Fidelity @ wikipedia <https://en.wikipedia.org/wiki/Fidelity_of_quantum_states>`_
+        `Hellinger Distance @ wikipedia <https://en.wikipedia.org/wiki/Hellinger_distance>`_
+    """
+    dist = _hellinger_distance(dist_p, dist_q)
+    return (1-dist**2)**2
diff --git a/releasenotes/notes/fix_hellinger_fidelity_def-7f60d53e3577fdaf.yaml b/releasenotes/notes/fix_hellinger_fidelity_def-7f60d53e3577fdaf.yaml
new file mode 100644
index 000000000000..02a18587ffbf
--- /dev/null
+++ b/releasenotes/notes/fix_hellinger_fidelity_def-7f60d53e3577fdaf.yaml
@@ -0,0 +1,6 @@
+---
+fixes:
+  - |
+    The definition of the Hellinger fidelity from has been corrected from the
+    previous defition of :math`1-H(P,Q)` to :math:`[1-H(P,Q)**2]**2` so that it
+    is equal to the quantum state fidelity of P, Q as diagonal density matrices.