We define an effective divisor on $\bR^N$ to be a function with finite support $\mu:\bR^N\to\bZ_{\geq 0}$. Its mass, denoted by $\bm(\mu)$, is the nonnegative integer
$$\bm(\mu)=\sum_{\bp\in\bR^N} \mu(\bp). $$
We denote by $\Div_+(\bR^N)$ the set of effective divisors. Note that $\Div_+(\bR^N)$ has a natural structure of Abelian semigroup.
For any $\bp\in\bR^N$ we denote by $\delta_\bp$ the Dirac divisor of mass $1$ and supported at $\bp$. The Dirac divisors generate the semigroup $\Div_+(\bR^N)$. We have a natural topology on $\Div_+(\bR^N)$ where $\mu_n\to \mu$ if and only if
$$ \bm(\mu_n)\to \bm(\mu),\;\; {\rm dist}\,\bigr(\;\supp(\mu_n),\; \supp(\mu)\;\bigr)\to 0, $$
where ${\rm dist}$ denotes the Haudorff distance.
A center of mass is a map
$$\eC:\Div_+(\bR^N)\to\Div_+(\bR^N) $$
satisfying the following conditions.
1. (Localization) For any divisor $\mu$ the support of $\eC(\mu)$ consists of a single point $\bc(\mu)$.
2. (Conservation of mass)
$$\bm(\mu)=\bm\bigl(\;\eC(\mu)\;\bigr),\;\;\forall\mu \in\Div_+(\bR^N),$$
so that
$$\eC(\mu)=\bm(\mu)\delta_{\bc(\mu)},\;\;\forall\mu \in\Div_+(\bR^N).$$
3. (Normalization)
$$\bc(m\delta_\bp)=\bp,\;\;\bc(\delta_\bp+\delta_\bq)=\frac{1}{2}(\bp+\bq),\;\;\forall \bp,\bq\in\bR^N,\;\;m\in\bZ_{>0}. $$
4. (Additivity)
$$\eC(\mu_1+\mu_2)= \eC\bigl(\,\eC(\mu_1)+\eC(\mu_2)\,\bigr),\;\;\forall \mu_1,\mu_2\in \Div_+(\bR^N). $$
For example, the correspondence
$$\Div_+ \ni \mu\mapsto \eC_0(\mu)=\bm(\mu)\delta_{\bc_0(\mu)}\in\Div_+,\;\;\bc_0(\mu):=\frac{1}{\bm(\mu)}\sum_\bp \mu(\bp)\bp $$
is a center-of-mass map. I want to show that this is the only center of mass map.
Proposition If $\eC:\Div_+(\bR^N)\to \Div_+(\bR^N)$ is a center-of-mass map, then $\eC=\eC_0$.
Proof. We carry the proof in several steps.
Step 1 (Rescaling). We can write the additivity property as
$$ \bc(\mu_1+\mu_2) =\bc\bigl(\; \bm(\mu_1)\delta_{\bc(\mu_1)} +\bm(\mu_2)\delta_{\bc(\mu_2)}\;\bigr). $$
In particular, this implies that the rescaling property
$$\bc( k\mu)=\bc(\mu),\;\;\forall\mu \in\Div_+,\;\; k\in\bZ_{>0}. \tag{R}\label{R}$$
This follows by induction $k$. For $k=1$ it is obviously true. In general
$$ \bc( k\mu)=\bc\bigl(\;(k-1)\bm(\mu)\delta_{\bc(\;(k-1)\mu)}+\bm(\mu)\delta_{\bc(\mu)}\;\bigr) =\bc\bigl(\; k\bm(\mu)\delta_{\bc(\mu)}\;\bigr)={\bc(\mu)} $$
Step 2. (Equidistribution) For any $n>0$ and any collinear points $\bp_1,\dotsc,\bp_n$ such that
$$|\bp_1-\bp_2|=\cdots=|\bp_{n-1}-\bp_n| $$
we have
$$\eC\Bigl(\sum_{k=1}^n\delta_{\bp_k}\;\Bigr)=\eC_0\Bigl(\sum_{k=1}^n\delta_{\bp_k}\;\Bigr) \tag{E}\label{E}. $$
Equivalently, this means that
$$\bc\Bigl(\sum_{k=1}^n\delta_{\bp_k}\;\Bigr)=\bc_0\bigl(\sum_{k=1}^n\delta_{\bp_k}\;\bigr)={\frac{1}{n}(\bp_1+\cdots+\bp_n)}. $$
We will prove (\ref{E}) arguing by induction on $n$. For $n=1,2$ this follows from the normalization property. Assume that (\ref{E}) is true for any $n< m$. We want to prove it is true for $n=m$.
We distinguish two cases.
(a) $m$ is even, $m= 2m_0$. We set
$$\mu_1=\sum_{j=1}^{m_0} \delta_{\bp_j},\;\;\mu_2=\sum_{j=m_0+1}^{2m_0}\delta_{\bp_j}. $$
Then
$$\bc(\mu_1+\mu_2)= \bc\bigl(\; m_0\delta_{\bc(\mu1)}+m_0\delta{\bc(\mu_2)}\;\bigr) =\bc\bigl( \delta_{\bc(\mu_1)}+\delta_{\bc(\mu_2)}\;\bigr). \tag{1}\label{2}$$
By induction
$$ \bc(\mu_1)=\bc_0(\mu_1),\;\;\bc(\mu_2)=\bc_0(\mu_2). $$
The normalization condition now implies that
$$ \bc\bigl( \delta_{\bc(\mu_1)}+\delta_{\bc(\mu_2)}\;\bigr)=\bc_0\bigl( \delta_{\bc_0(\mu_1)}+\delta_{\bc_0(\mu_2)}\;\bigr). $$
Now run the arguments in (\ref{2}) in reverse, with $\bc$ replaced by $\bc_0$.
(b) $m$ is odd, $m=2m_0+1$. Define
$$\mu_1=\delta_{\bp_{m_0+1}},\;\;\mu_2'=\sum_{j<m_0+1}\delta_{\bp_j},\;\;\mu_2''=\sum_{j>m_0+1}\delta_{\bp_j},\;\;\mu_2=\mu_2'+\mu_2''. $$
(Observe that $\bp_{m_0+1}$ is the mid-point in the string of equidistant collinear points $\bp_1,\dotsc,\bp_{2m_0+1}$. ) We have
$$\eC(\mu_2'+\mu_2'')=\eC\bigl( \; \eC(\mu_2')+\eC(\mu_2'')\;\bigr). $$
By induction
$$ \eC(\mu_2')+\eC(\mu_2'')= \eC_0(\mu_2')+\eC_0(\mu_2'') =m_0\delta_{\bc_0(\mu_2')}+m_0\delta_{\bc_0(\mu_2'')}=m_0\bigl(\;\delta_{\bc_0(\mu_2')}+\delta_{\bc_0(\mu_2'')}\;\bigr). $$
Observing that
$$\frac{1}{2}\bigl(\bc_0(\mu_2')+\bc_0(\mu_2'')\;\bigr)=\bp_{m_0+1} $$
we deduce
$$\eC(\mu_2)= \eC(\mu_2'+\mu_2'')=m_0\eC\bigl( \delta_{\bc_0(\mu_2')}+\delta_{\bc_0(\mu_2'')}\;\bigr)=m_0\eC_0\bigl( \delta_{\bc_0(\mu_2')}+\delta_{\bc_0(\mu_2'')}\;\bigr)=2m_0\delta_{\bp_{m_0+1}}=2m_0\mu_1. $$
Finally we deduce
$$\eC(\mu)=\eC\bigl(\;\eC(\mu_1)+\eC(\mu_2)\;\bigr)=\eC\bigl(\;\eC(\mu_1)+2m_0\eC(\mu_1)\;\bigr)= (2m_0+1)\delta_{\bp_{m_0+1}}=\eC_0(\mu). $$
Step 3. (Replacement) We will show that for any distinct points $\bq_1,\bq_2$ and any positive integers $m_1,m_2$ we can find $(m_1+m_2)$ equidistant points $\bp_1,\dotsc,\bp_{m_1+m_2}$ on the line determined by $\bq_1$ and $\bq_2$ such that
$$m_1\delta_{\bq_1}=\eC_0\Bigl(\sum_{j=1}^{m_1} \delta_{\bp_j}\Bigr)=\eC\Bigl(\sum_{j=1}^{m_1} \delta_{\bp_j}\Bigr),\;\;\;m_2\delta_{\bq_2}=\eC_0\Bigl(\sum_{j=m_1+1}^{m_1+m_2} \delta_{\bp_j}\Bigr)=\eC\Bigl(\sum_{j=m_1+1}^{m_1+m_2} \delta_{\bp_j}\Bigr). $$
This is elementary. Without restricting the generality we can assume that $\bq_1$ and $\bq_2$ lie on an axis (or geodesic) $\bR$ of $\bR^N$, $\bq_0=0$ and $\bq_2=q>0$. Clearly we can find real numbers $x_0, r$, $r>0$, such that
$$
\frac{1}{m_1}\sum_{j=1}^{m_1}(x_0+j r)=0,\;\;\frac{1}{m_2}\sum_{j=m_1+1}^{m_1+m_2}(x_0+jr)=q. $$
Indeed, the above two equalities can be rewritten as
$$ x_0+\frac{m_1+1}{2} r=0, $$
$$ q=x_0 +m_1 r+\frac{m_2+1}{2}=x_0+\frac{m_1+1}{2} r+\frac{m_1+m_2}{2} r. $$
Now place the points $\bp_j$ at the locations $x_0+jr$.
Step 4. (Conclusion) We argue on by induction on mass that
$$\eC(\mu)=\eC_0(\mu),\;\;\forall \mu\in \Div_+\tag{2}\label{3}$$
Clearly, the normalization condition shows that (\ref{3}) is true if $\supp\mu$ consists of a single point, or if $\bm(\mu)\leq 2$.
In general if $\bm(\mu)>2$ we write $\mu=\mu_1+\mu_2$ where $m_1=\bm(\mu_1),m_2=\bm(\mu_2)<\bm(\mu)$.
By induction we have
$$ \eC(\mu)= \eC\bigl( \eC(\mu_1)+\eC(\mu_2)\bigr)=\eC(\;\eC_0(\mu_1)+\eC_0(\mu_2)\;\bigr). $$
If $\bc_0(\mu_1)=\bc_0(\mu_2)$ the divisors $\eC_0(\mu_1)$ $\eC_0(\mu_2)$ are supported at the same point and we are done. Suppose that $\bq_1=\bc_0(\mu_1)\neq\bc_0(\mu_2)=\bq_2$. By Step 3, we can find equidistant points $\bp_1,\dotsc,\bp_{m_1+m_2}$ such that
$$ m_1\delta_{\bq_1}=\eC\Bigl(\sum_{j=1}^{m_1} \delta_{\bp_j}\Bigr)= \eC_0\Bigl(\sum_{j=1}^{m_1} \delta_{\bp_j}\Bigr) $$
$$ m_2\delta_{\bq_2}=\eC\Bigl(\sum_{j=m_1+1}^{m_1+m_2} \delta_{\bp_j}\Bigr)=\eC_0\Bigl(\sum_{j=m_1+1}^{m_1+m_2} \delta_{\bp_j}\Bigr). $$
We deduce that
$$\eC(\mu)=\eC\Bigl(\sum_{k=1}^{m+1+m_2}\delta_{\bp_k}\Bigr),\;\; \eC_0(\mu)=\eC_0\Bigl(\sum_{k=1}^{m+1+m_2}\delta_{\bp_k}\Bigr). $$
The conclusion now follows from (\ref{E}). q.e.d
Remark. The above proof does not really use the linear structure. If we uses only the fact that any two points in $\bR^N$ determine a unique geodesic. The Normalization condition can be replaced by the equivalent one
$$ \bc(\delta_\bp+\delta_\bq)= \mbox{the midpoint of the geodesic segment $[\bp,\bq]$}. $$
If we replace $\bR^N$ with a hyperbolic space the same arguments show that there exists at most one center of mass map.