mleHMM <-
function(x,m,G=as.list(1:m),hmm=NULL,...){
  
#   #number of states
#   m = length(unlist(G))
  
  if(!is.null(hmm)){
    stateParameters = hmm@stateParameters
    gamma = hmm@transitionMatrix
    
    distribution = hmm@internVariables$distribution
  }
  else{
    
    #univariate?
    if(class(x)=="numeric"){
      distribution="normal"  
    }
    #multivariate?
    if(class(x)=="matrix"){
      if(dim(x)[2]==1)
        distribution="normal"    
      else
        distribution="mv_normal"  
    }  
    
    
  }


  if(distribution=="normal"){
    
    #if one of the stating parameters is missing calculate values with RHmm  
    if(is.null(hmm)){
      cat("No Starting Parameter supplied.\n")
      cat("Use EM-algorithm of package RHmm to estimate HMM...")
      emFit = HMMFit(obs=x,nStates=m,...)
      muIn = emFit$HMM$distribution$mean
      sigmaIn = sqrt(emFit$HMM$distribution$var)
      stateParameters = list(mu = muIn, sigma = sigmaIn)
      gamma = emFit$HMM$transMat
    }
          
    #transformation for state dep. parameters
    n2w = function(x){c(x[1],log(x[2]))}
    w2n = function(x){c(x[1],exp(x[2]))}    
    
    #density
    f0 = function(x,par){
      mu    = par[1]
      sigma = par[2]
      return(dnorm(x,mean=mu,sd=sigma))
    }

    #state dep. parameter as list
    startP = extern2intern(parameter=stateParameters,distribution=distribution,m=m)
    
    
    cat(" done. \n")
    cat("Use nlm to improve fit...")
    
    ergMLE = mleHMMintern(x=x,m=m,G=G,f=f0,pn2pw=n2w,pw2pn=w2n,
                        startParameters=startP,startGamma=gamma,...)
    cat(" done. \n")
    
  }
  
  if(distribution=="mv_normal"){
    
    #if one of the stating parameters is missing calculate values with RHmm  
    if(is.null(hmm)){
      cat("No Starting Parameter supplied.\n")
      cat("Use EM-algorithm of package RHmm to estimate HMM...")
      
      emFit = HMMFit(obs=x,nStates=m,...)
      
      sigmaIn = array(dim=c(dim(x)[2],dim(x)[2],m))
      muIn = matrix(nrow=dim(x)[2],ncol=m)
      for(i in 1:m){
        sigmaIn[,,i] = emFit$HMM$distribution$cov[[i]]
        muIn[,i] = emFit$HMM$distribution$mean[[i]]
      }
      stateParameters = list(mu = muIn, sigma = sigmaIn)
      gamma = emFit$HMM$transMat
    
      cat(" done. \n")
      cat("Use nlm to improve fit...\n")
    }

  
    #get dimension of dataset -> needed for transformation w2n()
    D = dim(x)[2]
    
    #par is list with entries mu (vector) an sigma (matrix)
    n2w = function(par){
      #print(par)   
      #get parameters out of list
      mu = par$mu
      sigma = par$sigma
      
      #cholesky decomposition
      sigmaChol = t(chol(sigma))
      
      #write parameter to vector
      para = diag(sigmaChol)
      for(d in 1:(D-1)){
        para = c(para,sigmaChol[(d+1):D,d])
      }
      
      para = c(para,mu)
      return(para)
      
    }
    
    w2n = function(par){
      
      sigmaTMP = diag(par[1:D])  
      par = par[-(1:D)]
      for(d in 1:(D-1)){
        sigmaTMP[(d+1):D,d] = par[1:(D-d)]
        par = par[-c(1:(D-d))]
      }
      
      sigma = sigmaTMP %*% t(sigmaTMP)
      
      return(list(mu=tail(x=par,n=D),sigma=sigma))
      
    }

    #density
    f0 = function(x,par){
      mu    = par$mu
      sigma = par$sigma
      return(dmnorm(x=x,mean=mu,varcov=sigma))
    }
    
    
    #state dep. parameter as a list
    startP = extern2intern(parameter=stateParameters,distribution=distribution,m=m)
    

    
    ergMLE = mleHMMintern(x=x,m=m,G=G,f=f0,pn2pw=n2w,pw2pn=w2n,
                        startParameters=startP,startGamma=gamma,...)
    
  }
  
  #transform state dep. to nicer output
  ergMLE$mle$stateDepParameters = intern2extern(ergMLE$mle$stateDepParameters,distribution,m)

  par_est = setHMM(transitionMatrix=ergMLE$mle$gamma,mu=ergMLE$mle$stateDepParameters$mu,sigma=ergMLE$mle$stateDepParameters$sigma)
  
  ergMLE$mle = par_est
  
  return(ergMLE)
}
