How to create soil texture ternary plot?

I have a data frame with soil texture information and I want to create a soil texture plot similar to this.

How can I do this with:

1. different colors per area and
2. labels for classes within triangles like this ?

example data:

area <- c('S1','S2','S3','S4','S5','S6','S7','S8')
sand <- c(76.4,56.9,61.7,64.5,71,70.1,60.5,53.7)
silt<-c(9.3,23.1,23,17.4,13.5,13.4,21.1,30.6)
clay<-c(14.3, 20,15.4,18,15.5,16.6,18.4,15.7)
my_data<-data.frame(area,sand,silt,clay)

Basic plot with the ggtern package:

theme_set(theme_bw())
my_data %>%
    ggtern(aes(
        x = sand,
        y = clay,
        z = silt )) +
    geom_point()

with your data:

library(plotrix)

soil.texture(my_data[,2:4],col.symbols=1:8,bg.symbols=1:8,point.labels=my_data$area,pch=21)
legend( x=1, 
        legend=my_data$area,
        col=1:8, 
        fill=1:8 )

which give you:

EDIT

To avoid the background of the label you have to look inside the source code of soil.texture. If you do, you will notice that the function uses boxed.labels to plot the text in the triangle plot, so you can't avoid this behavior by playing with the function's arguments. Instead, you can write your own version of the function, changing boxed.labels with the simpler text function, as follow:

my_soil_texture <- function (soiltexture = NULL, main = "", at = seq(0.1, 0.9, 
                                                  by = 0.1), axis.labels = c("percent sand", "percent silt", 
                                                                             "percent clay"), tick.labels = list(l = seq(10, 90, by = 10), 
                                                                                                                 r = seq(10, 90, by = 10), b = seq(10, 90, by = 10)), show.names = TRUE, 
          show.lines = TRUE, col.names = "gray", bg.names = par("bg"), 
          show.grid = FALSE, col.axis = "black", col.lines = "gray", 
          col.grid = "gray", lty.grid = 3, show.legend = FALSE, label.points = FALSE, 
          point.labels = NULL, col.symbols = "black", pch = par("pch"), 
          ...) 
{
  par(xpd = TRUE)
  triax.plot(x = NULL, main = main, at = at, axis.labels = axis.labels, 
             tick.labels = tick.labels, col.axis = col.axis, show.grid = show.grid, 
             col.grid = col.grid, lty.grid = lty.grid)
  arrows(0.12, 0.41, 0.22, 0.57, length = 0.15)
  arrows(0.78, 0.57, 0.88, 0.41, length = 0.15)
  arrows(0.6, -0.1, 0.38, -0.1, length = 0.15)
  if (show.lines) {
    triax.segments <- function(h1, h3, t1, t3, col) {
      segments(1 - h1 - h3/2, h3 * sin(pi/3), 1 - t1 - 
                 t3/2, t3 * sin(pi/3), col = col)
    }
    h1 <- c(85, 70, 80, 52, 52, 50, 20, 8, 52, 45, 45, 65, 
            45, 20, 20)/100
    h3 <- c(0, 0, 20, 20, 7, 0, 0, 12, 20, 27, 27, 35, 40, 
            27, 40)/100
    t1 <- c(90, 85, 52, 52, 43, 23, 8, 0, 45, 0, 45, 45, 
            0, 20, 0)/100
    t3 <- c(10, 15, 20, 7, 7, 27, 12, 12, 27, 27, 55, 35, 
            40, 40, 60)/100
    triax.segments(h1, h3, t1, t3, col.lines)
  }
  if (show.names) {
    xpos <- c(0.5, 0.7, 0.7, 0.73, 0.73, 0.5, 0.275, 0.275, 
              0.27, 0.27, 0.25, 0.135, 0.18, 0.055, 0.49, 0.72, 
              0.9)
    ypos <- c(0.66, 0.49, 0.44, 0.36, 0.32, 0.35, 0.43, 
              0.39, 0.3, 0.26, 0.13, 0.072, 0.032, 0.024, 0.18, 
              0.15, 0.06) * sin(pi/3)
    snames <- c("clay", "silty", "clay", "silty clay", "loam", 
                "clay loam", "sandy", "clay", "sandy clay", "loam", 
                "sandy loam", "loamy", "sand", "sand", "loam", "silt loam", 
                "silt")
    text(xpos, ypos, labels=snames)# here I switched from boxed.labels(xpos, ypos, snames, border = FALSE, xpad = 0.5) to text(xpos, ypos, labels=snames)
  }
  par(xpd = FALSE)
  if (is.null(soiltexture)) 
    return(NULL)
  soilpoints <- triax.points(soiltexture, show.legend = show.legend, 
                             label.points = label.points, point.labels = point.labels, 
                             col.symbols = col.symbols, pch = pch, ...)
  invisible(soilpoints)
}

So now you can just use your own function for the plot:

my_soil_texture(my_data[,2:4],col.symbols=1:8,bg.symbols=1:8,point.labels=my_data$area,pch=21,col.grid=3)
legend( x=1, 
        legend=my_data$area,
        col=1:8, 
        fill=1:8 )

which gives you:

I want to highlight a minor issue with the change just made. Since you can change the background color of the main triangle, boxed labels() tries to work out whether the white or black text will be more easily read based on the background color and displays the text accordingly. So switching from boxed.labels() to text() can cause problems if you want to change the background of your plot