Wow, that's a very interesting model! I can think of two approaches to your problem. One is simply to use the Solver routine in Excel. The other is more theoretical (at first), but might yield a good solution.
First we can simplify your model expression a tad:
\begin{align*}
S_n(x) =& \frac{a_n}{\sigma_1\sqrt{2\pi}}\left[ \exp\left(\frac{-(x-\mu_1)^2}{2\sigma_1^2}\right) - \exp\left(\frac{-(x-\mu'_1)^2}{2\sigma_1^2}\right)\right]\\
+ &\frac{b_n}{\sigma_2\sqrt{2\pi}}\left[ \exp\left(\frac{-(x-\mu_2)^2}{2\sigma_2^2}\right) - \exp\left(\frac{-(x-\mu'_2)^2}{2\sigma_2^2}\right)\right]\\
+ &\frac{c_n}{\sigma_3\sqrt{2\pi}}\left[ \exp\left(\frac{-(x-\mu_3)^2}{2\sigma_3^2}\right) - \exp\left(\frac{-(x-\mu'_3)^2}{2\sigma_3^2}\right)\right].
\end{align*}
(With exponentials that complicated, it's usually better typesetting practice to use the exp function.)
Now, for each $n,$ we have a series of data points $\{x_j, S_{nj}\}$ to which we wish to fit the model. We'll try least-squares regression (this is not a linear model, because the coefficients we're trying to find don't show up linearly in the model) with a cost function
$$C_n=\sum_{j=1}^N(S_n(x_j)-S_{nj})^2. $$
Here $S_n(x_j)$ refers to the model evaluated at the $x_j$ data point, and $S_{nj}$ refers to the known (measured) $S$ value of the $n$th spectrum at $x_j$.
To solve for $\{\sigma_i,\mu_i, \mu_i'\},\; 1\le i\le 3,$ the usual procedure is to take the partial derivatives of $C_n$ and set them equal to zero. This will be a system of nine equations. Now, we will need to assemble these complicated ingredients one step at a time. Here are some intermediate calculations:
\begin{align*}
\frac{\partial S_n(x_j)}{\partial \mu_1}&=\frac{a_n}{\sigma_1\sqrt{2\pi}}\,\exp\left(\frac{-(x-\mu_1)^2}{2\sigma_1^2}\right)\left(\frac{2(x-\mu_1)}{2\sigma_1^2}\right)=\frac{a_n(x-\mu_1)}{\sigma_1^3\sqrt{2\pi}}\,\exp\left(\frac{-(x-\mu_1)^2}{2\sigma_1^2}\right) \\
\frac{\partial S_n(x_j)}{\partial \mu_1'}&=-\frac{a_n(x-\mu_1')}{\sigma_1^3\sqrt{2\pi}}\,\exp\left(\frac{-(x-\mu_1')^2}{2\sigma_1^2}\right).
\end{align*}
The other partials w.r.t. the means $\mu_2, \mu_2', \mu_3,$ and $\mu_3'$ will be analogous. The standard deviation partials are much more complicated. First, let's absorb non-$\sigma$ constants in one term together, like this:
$$_{1}S_{n}(x)=\frac{_1A_n}{\sigma_1}\left[ \exp\left(-\frac{B_1}{\sigma_1^2}\right) - \exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right],$$
where $_1A_n=a_n/\sqrt{2\pi}, \; B_1=(x-\mu_1)^2/2,$ and $B_1'=(x-\mu_1')^2/2.$ I'm using the notation $_1S_n(x)$ to denote the first term in $S_n(x)$. Now we calculate the partial derivative thus:
\begin{align*}
\frac{\partial (_1S_n(x))}{\partial \sigma_1}&=-\frac{_1A_n}{\sigma_1^2}\left[ \exp\left(-\frac{B_1}{\sigma_1^2}\right) - \exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right]
+\frac{_1A_n}{\sigma_1}\left[ \left(\frac{2B_1}{\sigma_1^3}\right)\exp\left(-\frac{B_1}{\sigma_1^2}\right) - \left(\frac{2B_1'}{\sigma_1^3}\right)\exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right] \\
&=-\frac{\sigma_1^2\,_1A_n}{\sigma_1^4}\left[ \exp\left(-\frac{B_1}{\sigma_1^2}\right) - \exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right]
+\frac{2\,_1A_n}{\sigma_1^4}\left[ B_1\exp\left(-\frac{B_1}{\sigma_1^2}\right) - B_1'\exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right] \\
&=\frac{_1A_n}{\sigma_1^4} \left[(2B_1-\sigma_1^2)\exp\left(-\frac{B_1}{\sigma_1^2}\right)-(2B_1'-\sigma_1^2)\exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right] \\
&=\frac{a_n}{\sigma_1^4\sqrt{2\pi}} \left[((x-\mu_1)^2-\sigma_1^2)\exp\left(-\frac{(x-\mu_1)^2}{2\sigma_1^2}\right)-((x-\mu_1')^2-\sigma_1^2)\exp\left(-\frac{(x-\mu_1')^2}{2\sigma_1^2}\right)\right].
\end{align*}
A more standard (pun intended, of course) way to write this would probably (pun intended) be:
$$\frac{\partial (_1S_n(x))}{\partial \sigma_1}=
\frac{a_n}{\sigma_1^2\sqrt{2\pi}} \left[\left(\frac{(x-\mu_1)^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x-\mu_1)^2}{2\sigma_1^2}\right)-\left(\frac{(x-\mu_1')^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x-\mu_1')^2}{2\sigma_1^2}\right)\right] $$
Whew! And those were only the intermediate calculations! What we really need to calculate is $\partial C_n/\partial\mu_1,$ and so on. Thankfully, we've already done most of the hard work. We have
\begin{align*}
\frac{\partial C_n}{\partial\mu_1}&=2\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj})\frac{a_n(x_j-\mu_1)}{\sigma_1^3\sqrt{2\pi}}\,\exp\left(\frac{-(x_j-\mu_1)^2}{2\sigma_1^2}\right)\right] \\
\frac{\partial C_n}{\partial\mu_1'}&=-2\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj})\frac{a_n(x_j-\mu_1')}{\sigma_1^3\sqrt{2\pi}}\,\exp\left(\frac{-(x_j-\mu_1')^2}{2\sigma_1^2}\right)\right] \\
\frac{\partial C_n}{\partial\sigma_1}&=2\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj})\frac{a_n}{\sigma_1^2\sqrt{2\pi}} \left[\left(\frac{(x_j-\mu_1)^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x_j-\mu_1)^2}{2\sigma_1^2}\right)-\left(\frac{(x_j-\mu_1')^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x_j-\mu_1')^2}{2\sigma_1^2}\right)\right]\right].
\end{align*}
The formatting there is off, because MathJax on M.SE only has so much width it can use. Thankfully, these expressions simplify a bit when we set them equal to zero:
\begin{align*}
0&=\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj})(x_j-\mu_1)\,\exp\left(\frac{-(x_j-\mu_1)^2}{2\sigma_1^2}\right)\right] \\
0&=\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj})(x_j-\mu_1')\,\exp\left(\frac{-(x_j-\mu_1')^2}{2\sigma_1^2}\right)\right] \\
0&=\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj}) \left[\left(\frac{(x_j-\mu_1)^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x_j-\mu_1)^2}{2\sigma_1^2}\right)-\left(\frac{(x_j-\mu_1')^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x_j-\mu_1')^2}{2\sigma_1^2}\right)\right]\right].
\end{align*}
Plus, of course, we have the six other analogous equations.
These are certainly complicated formulae, but they're no insurmountable. It would be possible to code it all up in Python, for example, and use, say, a gradient descent algorithm to solve the system. See Andrew Ng's Machine Learning course for more information.
