function (x, y, PDF=F) 
{
    x=sort(x)
    y=sort(y)
    n=length(x)
    m=length(y)
    N=m+n
    if(PDF==T) pdf(file=paste("tiedrankn",n,"m",m,".pdf",sep=""),width=11)
    z = c(x, y)
    rz = rank(z)
    rz=sort(rz)
    rz.u=unique(rz)
    d=NULL
    for(i in 1:length(rz.u)){
        d[i]=sum(rz==rz.u[i]) #for calulating the variance of Ws.star.
    }
    out = combn(rz, length(x), FUN = sum)
    #calculate the breaks br for variable bin width histogram for out 
    out=sort(out)
    out.u = sort(unique(out))
    d.out=diff(out.u)
    kd.out=length(d.out)
    br=c(out.u[1]-d.out[1]/2,out.u[1:kd.out]+d.out/2,
              out.u[kd.out+1]+d.out[kd.out]/2)
    freq=NULL
    for(i in 1:length(out.u)){
     freq[i]=sum(out==out.u[i])
    }   
    mean.Ws.star=n*(N+1)/2
    var.Ws.star=m*n*(N+1)/12-m*n*sum(d^3-d)/(12*N*(N-1))
    freq=freq/length(out)
    mx=min(out)
    Mx=max(out)
    mx3=mx-3
    Mx3=Mx+3
    z=seq(mx3,Mx3,length.out=100)
    dz=dnorm(z,mean.Ws.star,sqrt(var.Ws.star))
    Mfreq=2.2*max(freq)
    hist(out,breaks=br,axes=F,xlab=expression(W[s]^"*"),ylab="relative frequency",
        xlim=c(mx3,Mx3),main="",col=c("grey","grey90"),ylim=c(0,Mfreq))
    rug(out.u,col="black",lwd=2,ticksize=.015)
    text(mx3,.99*Mfreq,
    "The histogram has variable width bins, spanning intervals from midpoint",adj=0)
    text(mx3,.95*Mfreq,"to midpoint of intervals adjacent to respective represented values.",adj=0)
    text(mx3,.91*Mfreq,"These values are marked below each bin and",adj=0)
    text(mx3,.87*Mfreq,"are not necessarily centered on each bin.",adj=0)
    axis(1)
    axis(2)
    abline(h=0)
    lines(z,dz,col="blue")
    if(PDF==T) dev.off()
    readline("hit return\n")
    if(PDF==T) pdf(file=paste("CDFtiedrankn",n,"m",m,".pdf",sep=""),width=11)
    out.step=stepfun(out.u,c(0,cumsum(freq)))
    plot.stepfun(out.step,pch=16,cex=.5,axes=F,ylab=expression(P(W[s]^"*"<=a)),
       xlab="a",main="",xlim=c(mx3,Mx3))
    abline(h=c(0,1))
    pz=pnorm(z,mean.Ws.star,sqrt(var.Ws.star))
    lines(z,pz,col="blue")
    axis(1)
    axis(2)
    text(mx3,.97,"Observations",adj=0)
    text(mx3,.92,paste("x = (",paste(x,collapse=","),"),  ",
      "y = (",paste(y,collapse=","),")"),adj=0)
    text(mx3,.85,"Midranks",adj=0)
    ranksxy=rank(c(x,y))
    Ws.star=sum(ranksxy[1:n])
    prob.Ws.star=sum(freq[out.u<=Ws.star])
    prob.Ws.star.c=sum(freq[out.u>=Ws.star])
    text(mx3,.8,substitute(S^"*"== svec~" ,  "~W[s]^"*"==Ws,list(svec=
         paste("(",paste(ranksxy[1:n],collapse=","),")"),Ws=Ws.star)),adj=0)
    text(mx3,.75,substitute(R^"*"== svec,list(svec=
         paste("(",paste(ranksxy[(n+1):N],collapse=","),")"))),adj=0)
    if(prob.Ws.star<=.5){
       segments(Ws.star,0,Ws.star,prob.Ws.star,col="green")
       segments(Ws.star,prob.Ws.star,mx-5,prob.Ws.star,col="black")
       text(Ws.star+.5,prob.Ws.star/2,substitute(P(W[s]^"*"<= ws)==pws,list(ws=Ws.star,
          pws=signif(prob.Ws.star,4))),adj=0,col="green")}else{
       segments(Ws.star,0,Ws.star,1-prob.Ws.star.c,col="black",lty=3)
       segments(Ws.star,1-prob.Ws.star.c,Ws.star,1,col="green",lty=1)
       segments(Ws.star,1-prob.Ws.star.c,Mx3,1-prob.Ws.star.c,col="black")
       text(Ws.star+2,1-prob.Ws.star.c-.03,substitute(P(W[s]^"*">= ws)==pws,list(ws=Ws.star,
       pws=signif(prob.Ws.star.c,4))),adj=0,col="green")
    }
    if(PDF==T) dev.off()
}