### Function plot only the legend of a markov random field, like in drawmaps
plot_legend <- function(ylim_spat, xlimits = c(1,40), ylimits = c(1,2), 
                               cex.legend = 2,  digits = 2, par_new = FALSE, 
                               shift_text = c(8.5,-7.5), col_map = NULL) {
    
    op <- par(no.readonly = TRUE)
    
    lowerlimit <- ylim_spat[1]
    upperlimit <- ylim_spat[2]
    
    
    if(par_new) par(mar=c(0,0,0,0),oma=c(0,0,0,0))
    
    plot(xlimits, ylimits, type = "n", axes = FALSE, col = "white", 
         xlab = "", ylab = "", main="")
    
    legend.colors <- col_map
    nrcolors <- length(legend.colors)
    
    ylo <- yro <- ylimits[1] + (0.70 * (ylimits[2] - ylimits[1]))/1
    ylu <- yru <- ylimits[1] + (0.30 * (ylimits[2] - ylimits[1]))/1
    tylu <- tyru <- ylimits[1] + 0.35
    xlu <- xlo <- xlimits[1] + 0.01 * (xlimits[2] - xlimits[1])
    xru <- xro <- xlimits[1] + 0.99 * (xlimits[2] - xlimits[1])
    step <- (xru - xlu)/nrcolors
    for (i in 0:(nrcolors - 1)) {
        polygon(c(xlo + step * i, xlo + step * (i + 1), 
                  xlu + step * (i + 1), xlu + step * i), 
                c(ylo, yro, yru, ylu), col = legend.colors[i + 1], 
                border = legend.colors[i + 1])
    }
    lines(c(xlo, xro, xru, xlu, xlo), c(ylo, yro, yru, 
                                        ylu, ylo), col = "black")
    text(xlu + shift_text[1] * step, tylu, round(lowerlimit, digits=digits), 
         col = "black", cex = cex.legend)
    text(xru + shift_text[2] * step, tyru, round(upperlimit, digits=digits), 
         col = "black", cex = cex.legend)
    par(op)
}

