You could overlay a base legend with transparent bty='n' legends and fine-tune using adj=.
plot(1, type="n", col=2)
legend('topright', legend=rep('', 9), lty=c(rep(NA, 3), rep(1, 6)), title='',
ncol=3, col=c(rep(NA, 3), 2, 5, 3, 6, 4, 7))
legend('topright', legend=c(paste0('T=', 1:3), rep('', 6)), title='', bty='n',
adj=.3, ncol=3, col=c(rep(1, 3), rep(NA, 6)))
legend('topright', legend=c('', 'foo', ' bar'), adj=.6, bty='n', ncol=3,
col=c(rep(1, 3), rep(NA, 6)))

Update
To be even more flexible you could use the par()$usr coordinates and define three adjustment parameters p*. For sake of stability I strongly recommend to use the png device or similar with fixed width and height.
png('foo.png', width=480, height=480)
plot(matrix(1:12, 3, 4), type="n", col=2)
pu <- par()$usr
p1 <- 2.1; p2 <- 1.9; p3 <- 1.4
legend(pu[3] - p1, pu[4], legend=rep('', 21), lty=rep(1, 21), title='',
ncol=7, col=rep(NA, 21))
legend(pu[3] - p1, pu[4], legend=c(paste0('T=', 1:3), rep('', 6)), title='', bty='n',
adj=.3, ncol=3, col=c(rep(1, 3), rep(NA, 6)))
legend(pu[3] - p2, pu[4], legend=rep('', 3), lty=1, title='', col=c(2, 5, 3), bty='n')
legend(pu[3] - p3, pu[4], legend=rep('', 3), lty=1, title='', col=c(6, 4, 7), bty='n')
legend(pu[3] - p2, pu[4], legend='fooooooooooo', bty='n', ncol=3, adj=.12)
legend(pu[3] - p3, pu[4], legend='baaaaaaaaaar', bty='n', ncol=3, adj=.12)
box()
dev.off()
