how to estimate parameters of multimodal distribution

Dear all

Please, is there any function or package for dealing with multimodal distributions? 
I try to fit multimodal distribution or more precisely to find out mixture of normal 
distributions which can lead to my actual data.

I use optim to find (in that case) two parameters but what I want is to let the 
function find out arbitrary number of normal distributions underlaying my actual 
data and estimate all parameters automatically. Actually I want to elaborate some 
function which is similar to procedures used for evaluating DTA, IR or XRD 
spectra. Is the optim way suitable for this task? 

I would be gratefull for **any** hint.

Here is some artificial example.

# making two normals and putting them together

x<-seq(5,100,5)
y1<-abs(dnorm(x,30,5)+rnorm(20,0,.002))
y2<-abs(dnorm(x,60,10)+rnorm(20,0,.002))
y<-y1+y2
y<-y/sum(y)*100

# *** my data actually look like this***
n<-round(y/sum(y)*100)

opt<-optim(c(5,10),fff)
matplot(x,cbind(y,fff2(c(5.3,10))),type="l")

# quite OK but I need a little bit more general solution

fff<-function(p) 
{
p1<-p[1]
p2<-p[2]
v<-dnorm(x,30,p1)+dnorm(x,60,p2)
s<-sum(v)
sum((y-v/s)^2)
}

fff2<-function(p) 
{
p1<-p[1]
p2<-p[2]
v<-dnorm(x,30,p1)+dnorm(x,60,p2)
s<-sum(v)
v/s
}

Thank you very much for any help.
Best regards.

Petr Pikal
petr.pikal at precheza.cz
p.pik at volny.cz