Building on @P-Lapointe solution, but making it extremely easy, you could use the maximum values from your data using max() and then you re-use those maximum values to set the legend xy coordinates. To make sure you don't get beyond the borders, you set up ylim slightly over the maximum values.
a=c(rnorm(1000))
b=c(rnorm(1000))
par(mfrow=c(1,2))
plot(a,ylim=c(0,max(a)+1))
legend(x=max(a)+0.5,legend="a",pch=1)
plot(a,b,ylim=c(0,max(b)+1),pch=2)
legend(x=max(b)-1.5,y=max(b)+1,legend="b",pch=2)
