# helper functiont that does all the work
overlap_single<-function(sd,a,tar,foil) 
{ 
  nchar(tar)->tarlen
  nchar(foil)->foillen

  if(tarlen!=length(sd)) error("Target Length does not match # of elements in sd");

  ovcalc<-function(sd,len1,len2) # Compute potential overlap of each element from: sd of stimulus, length of stimulus, length of option
  {
    outer(
      1:len1,(1:len2-.5)/len2,  # rescale y centres to 0-1 range
      Vectorize(function(x,y)
        { 
          xloc<-(x-.5)/len1 # rescale x centre
          l<-pnorm(y+.5/len2,xloc,sd[x]/len1)
          r<-pnorm(y-.5/len2,xloc,sd[x]/len1)
          return(abs(l-r))
        }  )
    );
  }

  ovtf<-ovcalc(sd,tarlen,foillen) # for the foil option (target is always stimulus)
  ovtt<-ovcalc(sd,tarlen,tarlen)  # for the correct option (target is always stimulus)

  # identify matching elements
  tfm<-outer(1:tarlen,1:foillen,Vectorize(function(x,y) substr(tar,x,x)==substr(foil,y,y)));  
  ttm<-outer(1:tarlen,1:tarlen,Vectorize(function(x,y) substr(tar,x,x)==substr(tar,y,y))); 

  ffm<-outer(1:foillen,1:foillen,Vectorize(function(x,y) substr(foil,x,x)==substr(foil,y,y))); # foilfoil matches used only to detect repeat

  ov_foil<-sum(ovtf*tfm)^a; 
  ov_target<-sum(ovtt*ttm)^a; 

  ov_target/(ov_foil+ov_target)
}

overlap<-Vectorize(overlap_single,vectorize.args=c("tar","foil"))

# call overlap like overlap(c(.5,.5,.5,.5,.5),2,c("ankle"),c("anxle"))
# parameters are sd(vector),a,target,foil)
# You can add more target-foil pairs of the same target length
# To do more than one length, you'll need a different sd vector anyway
# Return a vector of P(correct) for each trial

overlap_samediff_single<-function(sd,a,thr,tar,reference) 
{ 
  nchar(tar)->tarlen
  nchar(reference)->reflen

  if(tarlen!=length(sd)) error("Target Length does not match # of elements in sd");

  ovcalc<-function(sd,len1,len2) # Compute potential overlap of each element from: sd of stimulus, length of stimulus, length of option
  {
    outer(
      1:len1,(1:len2-.5)/len2,  # rescale y centres to 0-1 range
      Vectorize(function(x,y)
        { 
          xloc<-(x-.5)/len1 # rescale x centre
          l<-pnorm(y+.5/len2,xloc,sd[x]/len1)
          r<-pnorm(y-.5/len2,xloc,sd[x]/len1)
          return(abs(l-r))
        }  )
    );
  }

  ovtr<-ovcalc(sd,reflen,tarlen) # for the reference (target is always stimulus)
#  ovtt<-ovcalc(sd,tarlen,tarlen) # for the target (target is always stimulus)
  # identify matching elements
  trm<-outer(1:tarlen,1:reflen,Vectorize(function(x,y) substr(tar,x,x)==substr(reference,y,y)));  

  ttm<-outer(1:tarlen,1:tarlen,Vectorize(function(x,y) substr(tar,x,x)==substr(tar,y,y))); 

  rrm<-outer(1:reflen,1:reflen,Vectorize(function(x,y) substr(reference,x,x)==substr(reference,y,y))); # refref matches used only to detect repeat

  ov_different<-thr^a; 
  ov_same<-sum(ovtr*trm)^a; 

  ifelse(tar==reference,ov_same,ov_different)/(ov_different+ov_same)
}

overlap_samediff<-Vectorize(overlap_samediff_single,vectorize.args=c("tar","reference"))

# call overlap_samediff like  overlap_samediff(c(.5,.5,.5,.5,.5),2,3,3,"ankle","ankle")
# parameters are sd(vector),a_nonrepeat,a_repeat,soft_threshold,target,reference)
# You can add more target-reference pairs of the same target length
# To do more than one length, you'll need a different sd vector anyway
# Return a vector of P(correct) for each trial

###################################################

#exp 1 data
datafexp1<-read.table("z../datashortExp12AFC.txt")

subset(datafexp1,V7==7)->exp1tar7
exp1tar7$V5->exp1tar7tar
exp1tar7$V6->exp1tar7foil

# call overlap like overlap(c(.5,.5,.5,.5,.5),2,3,c("ankle"),c("anxle"))
overlap(c(.5,.5,.5,.5,.5,.5,.5),2,c(as.character(exp1tar7tar)),c(as.character(exp1tar7foil)))
overlap(c(.5,.5,.5,.5,.5,.5,.5),2,c("dartien"),c("drartien"))

## call overlap like overlap(c(.5,.5,.5,.5,.5),2,3,c("ankle"),c("anxle"))
# parameters are sd(vector),a_nonrepeat,a_repeat,target,foil)

#create a function that runs the overlap for all target lengths,combines all data, calculates sse and 
#includes a summary


exp1fitOV<- function(p,summary=0){
  
  a<-p[1]
  #sd7<-p[2:8]
  #sd8<-p[9:16]
  d<-p[2] 
  r<-p[3]
  #d=1.5;r=1.2
  sd7<-d*(1-exp(-(1:7-.5)/r))
  sd8<-d*(1-exp(-(1:8-.5)/r))
  
  if(any(c(sd7,sd8)<0)) return(1e+16)
  
  subset(datafexp1,V7==7)->tar7
  tar7$V5->tar7tar
  tar7$V6->tar7foil
  subset(datafexp1,V7==8)->tar8
  tar8$V5->tar8tar
  tar8$V6->tar8foil
  overlap(sd7,a,as.character(tar7tar),as.character(tar7foil))->overlap7
  overlap(sd8,a,as.character(tar8tar),as.character(tar8foil))->overlap8
  c(overlap7,overlap8)->pc
  
  sse<-sum((datafexp1$V2-pc)^2)
  
  if(summary>0)
  {
    tempd<-datafexp1
    tempd$V2<-pc
    
    list(fit=sse,
         summary=rbind(data.frame(src="model",tempd),
                       data.frame(src="data",datafexp1))
    )
  }
  else
  {
    sse
  }
  
}


#optim(c(2,.5,.5,.5,.5,.5,.5,.5,.5,.5,.5,.5,.5,.5,.5,.5),
#      fn=exp1fitOV,control=list(maxit=5000,trace=T))->exp1ov

optim(c(2,1.5,1.2),
      fn=exp1fitOV,control=list(maxit=5000,trace=T))->exp1ov

exp1ov$par

#optim(exp1ov$par,
#      fn=exp1fitOV,control=list(maxit=50000,trace=T))->exp1ov2


#round(exp1ov2$par,2)
round(exp1ov$par,2)

#exp1fitOV(exp1ov2$par,summary=1)$summary->ovpredictionsexp1
exp1fitOV(exp1ov$par,summary=1)$summary->ovpredictionsexp1

with(subset(ovpredictionsexp1,src=="model"),xtabs(V2~V3+V7)/xtabs(~V3+V7))
with(subset(ovpredictionsexp1,src=="data"),xtabs(V2~V3+V7)/xtabs(~V3+V7))

#exports

#export OV predictions means Exp1
write.csv(ovpredictionsexp1,"../OVEXP1noBiasRd.txt",quote=FALSE,row.names =FALSE)

