library(stockassessment)

oldwd<-setwd("data")

catch.no<-read.ices('cn.dat')
catch.mean.weight<-read.ices('cw.dat')
dis.mean.weight<-read.ices('dw.dat')
land.mean.weight<-read.ices('lw.dat')
stock.mean.weight<-read.ices('sw.dat')
prop.mature<-read.ices('mo_raw.dat')
natural.mortality<-read.ices('nm.dat')
surveys<-read.ices('survey.dat')
land.no<-read.ices('lf.dat')
dis.no<-catch.no-land.no 
prop.f<-read.ices('pf.dat')
prop.m<-read.ices('pm.dat')


# redirect output to avoid mgcv message generating warning icons
zz=file("dummy",open="wt")
sink(zz)
sink(zz,type="message")
library(mgcv)
sink(zz,type="message")
sink()
file.remove("dummy")

# Smooth maturity

skipYears=c(1:10)
columnsToSmooth=1:5
mo=prop.mature[-c(skipYears),]

for(cc in columnsToSmooth){
    ww = mo[,cc];
    tt = 1:length(ww)
    tmp = gam(ww ~ s(tt))
    mo[,cc] = predict(tmp);
}
prop.mature[-c(skipYears),]=mo

head=readLines("mo_raw.dat",n=5);

file="mo.dat"
cat(paste(head[1],"smoothed"),"\n",file=file)
cat(head[2],"\n",file=file,append=TRUE)
cat(head[3],"\n",file=file,append=TRUE)
cat(head[4],"\n",file=file,append=TRUE)
cat(head[5],"\n",file=file,append=TRUE)
write.table(round(prop.mature,4),file=file,row.names=FALSE,col.names=FALSE,append=TRUE)


setwd(oldwd)

# Modify to 6+ data
cutage<-6
low<-1
GE<-which(as.numeric(colnames(catch.no))>=cutage)
E<-which(as.numeric(colnames(catch.no))==cutage)
w<-catch.no[,GE]/rowSums(catch.no[,GE])
wex<-rbind(w,w[nrow(w),])
wD<-dis.no[,GE]/ifelse(rowSums(dis.no[,GE])>0,rowSums(dis.no[,GE]),1)
wL<-land.no[,GE]/rowSums(land.no[,GE])
catch.no[,E]<-rowSums(catch.no[,GE])
catch.no<-catch.no[,low:E]
prop.mature[,E]<-rowSums(prop.mature[,GE]*wex)
prop.mature<-prop.mature[,low:E]
stock.mean.weight[,E]<-rowSums(stock.mean.weight[,GE]*wex)
stock.mean.weight<-stock.mean.weight[,low:E]
catch.mean.weight[,E]<-rowSums(catch.mean.weight[,GE]*w)
catch.mean.weight<-catch.mean.weight[,low:E]
dis.mean.weight[,E]<-rowSums(dis.mean.weight[,GE]*wD)
dis.mean.weight<-dis.mean.weight[,low:E]
land.mean.weight[,E]<-rowSums(land.mean.weight[,GE]*wL)
land.mean.weight<-land.mean.weight[,low:E]
natural.mortality[,E]<-rowSums(natural.mortality[,GE]*wex)
natural.mortality<-natural.mortality[,low:E]
land.no[,E]<-rowSums(land.no[,GE])
land.no<-land.no[,low:E]
land.frac<-ifelse(catch.no>0,land.no/catch.no,1)
prop.f<-prop.f[,low:E]
prop.m<-prop.m[,low:E]

dat<-setup.sam.data(surveys=surveys,
                    residual.fleet=catch.no, 
                    prop.mature=prop.mature, 
                    stock.mean.weight=stock.mean.weight, 
                    catch.mean.weight=catch.mean.weight, 
                    dis.mean.weight=dis.mean.weight, 
                    land.mean.weight=land.mean.weight,
                    prop.f=prop.f, 
                    prop.m=prop.m, 
                    natural.mortality=natural.mortality, 
                    land.frac=land.frac)
                    
save(dat, file="run/data.RData")