This will easier in Matlab than Origin I think, but in principle you can use any software that has a least-squares or similar fitting algorithms implemented. XPS peaks are Voigt functions, a convolution of a Gaussian (instrumental) and Lorentzian (natural) peak shape with a particular mixing ratio. True Voigt functions are computationally demanding so most software using sum or profuct approximations with a fixed mixing ratio (typically 70% Lorentzian).See http://onlinelibrary.wiley.com/doi/10.1002/sia.2527/full for formulae. Note you will be also need to do a manual subtraction of the background (usually fit with a Shirley background but linear are often used as well).

You'll need to script these formulae and use a least-squares fitting procedure to determine the position, area and width (FWHM) of the fitted peak cf. the data. Obviously XP spectra usually contain many peaks, so you will need to do a fitting of a set of different peaks, using as many constraints as possible (fitting generally degrades rapidly the more variables you ask of it). For instance, FWHM of several different  species can sometimes be approximated as equal (or within a narrow range of equal), and the area ratio of spin-orbit split peaks can be fixed (1:2 for p orbitals, 2:3 for d orbitals, 3:4 for f orbitals). FWHM of spin-orbit  split peaks can also be fixed at a ratio (often close to equal, but e.g. for Ti 2p it's closer to double) as can the binding energy separation if you have several chemical species of the doublets. Because you'll be doing this manually you will need to be very careful how you constrain the fitting procedure.

This is doable, but time consuming; I've done this in Matlab when analysing thousands of similar spectra as part of time-resolved experiments (measuring small shifts in binding energy). However, software like CASA XPS make fitting many different core levels very easy; fitting yourself will require you to modify the fitting procedure every time you look at a new peak. Also, for metals and high-conductivity materials like graphite, assymetric peaks are required which will be more difficult to implement (easiest way is to multiply your Voigt approximation by an exponential to give it a tail). Commercial software have dozens of different lineshapes available including the more proper Doniach Sunjic lineshape. 

Good luck!

I am only familiar with Origin.

You can do two things here:

1. use the tools in -> -> . Here the problem is that you cannot set the peak values of your Gaussian bells

2. press Ctrl+Y, then follow the dialog. You need to select e.g. Gaussian functions, then go to / and select the number of Replicas (i.e. how many peak functions you want to fit). Then you can go to , you enter the xc-values for each of the peaks you expect and set these parameters as fixed, then try to get your best fit. Play a little around with some simple or reference data to get a hang of that too. It is quite good but needs some practicing to use it properly.

thank you Prof.Ben F. Spencer  and Prof.Florian J. Stadler , I will do as you suggested