.packageName <- "simpleaffy"
"call.exprs" <-
function(x,algorithm="rma",do.log = TRUE,sc=100,method=NA) {
  if(algorithm == "rma-R") {       # use RMA
    if(is.na(method)) {
      method<-"quantiles";
    }
    tmp <- expresso(x, bgcorrect.method="rma",
                      normalize.method=method,pmcorrect.method="pmonly",
                      summary.method="avgdiff");
    if(!do.log) {
      exprs(tmp) <- log2(exprs(tmp));
    }
    return(tmp);
  }
   if(algorithm == "rma") {       # use RMA
    tmp <- rma(x);
    if(!do.log) {
      exprs(tmp) <- 2^exprs(tmp);
    }
    return(tmp);
  }
  else if(algorithm == "gcrma") {       # use GCRMA
    if(is.na(method)) {
      method<-"quantiles";
    }
    tmp <- gcrma(x);
    if(!do.log) {
      exprs(tmp) <- 2^exprs(tmp);
    }
    return(tmp);
  }

  else if(algorithm == "mas5-R") { # use expresso MAS5.0
     if(is.na(method)) {
     tmp1 <- expresso(x, normalize=FALSE,bgcorrect.method="mas",
                        pmcorrect.method="mas",summary.method="mas");
        tmp  <- affy.scalevalue.exprSet(tmp1,sc=sc);
     }
     else {
     tmp1 <- expresso(x, normalize.method=method,bgcorrect.method="mas",
                        pmcorrect.method="mas",summary.method="mas");
        tmp  <- affy.scalevalue.exprSet(tmp1,sc=sc);
     }
     if(do.log) {
       exprs(tmp) <- log2(exprs(tmp));
       description(tmp)@preprocessing$sfs=apply(2^(exprs(tmp) - log2(exprs(tmp1))),2,mean)
       description(tmp)@preprocessing$tgt=sc
     }
     else {
       description(tmp)@preprocessing$sfs=apply(2^(exprs(tmp) - log2(exprs(tmp1))),2,mean)
       description(tmp)@preprocessing$tgt=sc

     }
     return(tmp);
  }
  else if(algorithm == "mas5") { # use Simpleaffy MAS5.0
    tmp <- justMAS(x,tgt=sc);
    if(!do.log) {
      exprs(tmp) <- 2^exprs(tmp);
    }
    return(tmp);
  }
  else {
    stop(paste("Don't know how to compute algorithm '",algorithm,"'."));
  }
}

detection.p.val <- function(x, tau = 0.015,calls=TRUE,alpha1=getAlpha1(cleancdfname(cdfName(x))),alpha2=getAlpha2(cleancdfname(cdfName(x))),ignore.saturated=TRUE) {
  if(alpha1 < 0)      {stop("alpha1 must be  > 0 "); }
  if(alpha1 > alpha2) {stop("alpha2 must be  > alpha1 "); }
  if(alpha2 > 1)      {stop("alpha2 must be  <1 "); }

  cat("Getting probe level data...\n");
  pms <-as.matrix(pm(x));
  mms <-as.matrix(mm(x));

  # Saturation:
  # shouldn't be a problem with new scanners or those that have had an engineer visit
  if(ignore.saturated) { sat <- 46000; }
  else { sat <- -1; }
  
  pns <- probeNames(x);
  unique.pns <- unique(pns);
  cat("Computing p-values\n");
  p<-sapply(1:length(pms[1,]),function(x) { 
    .C("DetectionPValue",as.double(pms[,x]),as.double(mms[,x]),as.character(pns),as.integer(length(mms[,x])),
	as.double(tau),as.double(sat),dpval=double(length(unique.pns)),length(unique.pns),PACKAGE="simpleaffy")$dpval;
  });
  rownames(p) <- unique.pns;
  colnames(p) <- paste(sampleNames(x),".detection.p.val",sep="");
  if(!calls) { 
    l <- list(detection.p.values=p); 
  }
  else       {
    cat("Doing PMA Calls\n");
    calls <- sapply(p,function(y) { if(y < alpha1) { return("P") } else { if(y < alpha2) { return("M") } else { return("A") }}});
    calls <- matrix(calls,nrow=nrow(p),ncol=ncol(p));
    colnames(calls) <- paste(sampleNames(x),".present",sep="");
    rownames(calls) <- rownames(p)
    l<- list(pval=p,call=calls);
    return(l); 
  }

}     

get.annotation <- function(x,cdfname) {
  library(cdfname,character.only=T);
  symb  <- unlist(multiget(x,env=get(paste(cdfname,"SYMBOL",sep=""))),use.names=F);
  desc  <- unlist(multiget(x,env=get(paste(cdfname,"GENENAME",sep=""))),use.names=F);
  accno <- unlist(multiget(x,env=get(paste(cdfname,"ACCNUM",sep=""))),use.names=F);
  uni   <- unlist(multiget(x,env=get(paste(cdfname,"UNIGENE",sep=""))),use.names=F);
  acc.lnk <- paste("=HYPERLINK(\"http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?cmd=search&db=nucleotide&term=",accno,"\",\"",accno,"\")",sep="");
  uni.lnk <- paste("=HYPERLINK(\"http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?cmd=search&db=unigene&term=",uni,"&dopt=unigene\",\"",uni,"\")",sep="")	
  res <- cbind(symb,acc.lnk,uni.lnk,desc);
  colnames(res) <- c("gene name","accession","unigene","description");
  return(res);
}

write.annotation <- function(summary,file="results/annotation.table.xls") {
  write.table(summary,file=file,sep="\t",quote=F,col.names=NA)
}

results.summary <- function(results,cdfname) {
  res <- cbind(means(results),fc(results),sapply(fc(results),function(x) { if(x <0) { -1 * 2 ^ (-1 * x) } else { 2^x } }),tt(results),get.annotation(names(fc(results)),cdfname));
  cns <- colnames(res);
  cns[3]<-"log2(fc)";
  cns[4]<-"fc";
  cns[5]<-"tt";
  colnames(res) <- cns;
  return(res);
}


journalpng <- function(file="figure.png",width=4, height=4,res=300) {
  bitmap(file=file,type="png16m",width=width,height=height,res=res)
}

screenpng <- function(file="figure.png",width=4, height=4,res=72) {
  bitmap(file=file,type="png16m",width=width,height=height,res=res)
}
bg.correct.sa <- function(unnormalised, grid=c(4,4)) {
res         <- unnormalised;
pms         <- unique(unlist(pmindex(res))) - 1 # C counts from 0
mms         <- unique(unlist(mmindex(res))) - 1 # C counts from 0
all         <- c(pms,mms)
intensities <- intensity(res)
rws <- nrow(res)
cls <- ncol(res)
for(no in 1:length(res)){
  this.array <- intensities[,no];
  result <- .C("bgmas",as.integer(as.vector(all)),as.integer(length(all)),
       as.double(as.vector(this.array)),as.integer(length(this.array)),
       as.integer(rws),
       as.integer(cls),
       as.integer(grid[1]),as.integer(grid[2]),
       zonebg=double(grid[1] * grid[2]),
       zonesd=double(grid[1] * grid[2]),corrected=double(length(this.array)),PACKAGE="simpleaffy");
       intensities[,no] <- result$corrected;
  }
  intensity(res) <- intensities;
  return(res);
}

justMAS     <- function(unnormalised,tgt=100,scale=TRUE) {
  ct <- 0.03;
  st <- 10.0;
########################### BACKGROUND
  cat("Background correcting\n");
  bgc <- bg.correct.sa(unnormalised);

  cat("Retrieving data from AffyBatch...");
  pms <-as.matrix(pm(bgc))
  mms <-as.matrix(mm(bgc))
  pns <- probeNames(bgc);
  unique.pns <- unique(pns);
########################### SUMMARIES
  cat("done.\nComputing expression calls... \n");
  expression.calls<-sapply(1:length(pms[1,]),function(x) { 
    cat(".");
    .C("GetExpressionLevels",as.double(pms[,x]),as.double(mms[,x]),as.character(pns),as.integer(length(mms[,x])),
	as.double(ct),as.double(st),exprs=double(length(unique.pns)),length(unique.pns),PACKAGE="simpleaffy")$exprs;
  });
  cat("done.\n");
  rownames(expression.calls) <- unique.pns;
  colnames(expression.calls) <- paste(sampleNames(bgc))
########################### SCALING
  if(scale) {
   cat(paste("scaling to a TGT of",tgt,"..."));
   sfs <- double(length(expression.calls[1,]));
   for(x in 1:length(expression.calls[1,])) {
     vals <- sort(2^expression.calls[,x]);
     l <- length(vals);
     frm <- 0.02 *l;
     to  <- 0.98 *l;
     sf  <- tgt/mean(vals[frm:to]);
     cat(paste("Scale factor for:",sampleNames(unnormalised)[x],sf,"\n"))
     expression.calls[,x] <- log2((2^expression.calls[,x]) * sf)
     sfs[x] <- sf; 
   }
  }
  else {
   res@description@preprocessing$sfs = stop("Arrays were not scaled") 
   res@description@preprocessing$tgt = stop("Arrays were not scaled") 
  }
  res <- new("exprSet", 
             exprs       = expression.calls,
             phenoData   = bgc@phenoData,
             annotation  = bgc@annotation, 
             description = bgc@description, 
             notes       = bgc@notes);

  res@description@preprocessing$sfs = unlist(sfs);
  res@description@preprocessing$tgt = tgt;
  return(res);
}


.mas5     <- function(unnormalised,normalize=TRUE,sc=500,analysis="absolute") {
  ct <- 0.03;
  st <- 10.0;
########################### BACKGROUND
  if(normalize) {
    bgc <- bg.correct.sa(unnormalised);
  }
  pms <-pm(bgc)
  mms <-mm(bgc)
  pns <- probeNames(bgc);
  unique.pns <- unique(pns);
########################### SUMMARIES
  expression.calls<-sapply(1:length(pms[1,]),function(x) { 
    .C("GetExpressionLevels",as.double(pms[,x]),as.double(mms[,x]),as.character(pns),as.integer(length(mms[,x])),
	as.double(ct),as.double(st),exprs=double(length(unique.pns)),length(unique.pns),PACKAGE="simpleaffy")$exprs;
  });
  rownames(expression.calls) <- unique.pns;
  colnames(expression.calls) <- paste(sampleNames(bgc))
########################### SCALING
  sfs <- double(length(expression.calls[1,]));
  for(x in 1:length(expression.calls[1,])) {
    vals <- sort(2^expression.calls[,x]);
    l <- length(vals);
    frm <- 0.02 *l;
    to  <- 0.98 *l;
    sf  <- sc/mean(vals[frm:to]);
    expression.calls[,x] <- log2((2^expression.calls[,x]) * sf)
    sfs[x] <- sf; 
  }
  res <- new("exprSet", 
             exprs       = expression.calls,
             phenoData   = bgc@phenoData,
             annotation  = bgc@annotation, 
             description = bgc@description, 
             notes       = bgc@notes);
  res@description@preprocessing$sfs = unlist(sfs);
  res@description@preprocessing$tgt = sc;
  return(res);
}

library("methods")
library("affy")
# holds the results of a pairwise comparison
setClass("PairComp",representation(means="matrix",fc="numeric",tt="numeric",calls="matrix",group="character",members="character"))

#accessor methods
setGeneric("means", function(object) standardGeneric("means"))
setMethod("means","PairComp",function(object) object@means)

setGeneric("fc", function(object) standardGeneric("fc"))
setMethod("fc","PairComp",function(object) object@fc)

setGeneric("tt", function(object) standardGeneric("tt"))
setMethod("tt","PairComp",function(object) object@tt)

setGeneric("calls", function(object) standardGeneric("calls"))
setMethod("calls","PairComp",function(object) object@calls)

setGeneric("group", function(object) standardGeneric("group"))
setMethod("group","PairComp",function(object) object@group)

setGeneric("members", function(object) standardGeneric("members"))
setMethod("members","PairComp",function(object) object@members)


##subseting. can only happen by gene

setMethod("[", "PairComp", function(x,i,j,...,drop=FALSE) {
  if(nrow(calls(x))>0) {calls <- calls(x)[i,,...,drop=FALSE];}
  else { calls <- matrix(); }
  y <- new ("PairComp",means=means(x)[i,,...,drop=FALSE],fc=fc(x)[i,...,drop=FALSE],tt=tt(x)[i,...,drop=FALSE],calls=calls,group=group(x),members=members(x))
  return(y)
})

setReplaceMethod("[", "PairComp", function(x, i,j,...,value) {
  stop("operation not supported");
})



"get.array.subset.exprset" <-
function(x,group,members) {
  pd <- pData(x);
  grp <- pd[,colnames(pd) == group];
  return(x[,is.element(grp,members)]);
}

"get.array.subset.affybatch" <-
function(x,group,members) {
  pd <- pData(x);
  grp <- pd[,colnames(pd) == group];
  return(x[,is.element(grp,members)]);
}

setGeneric("get.array.subset", function(x,group,members) standardGeneric("get.array.subset"))
setMethod("get.array.subset","AffyBatch",get.array.subset.affybatch);
setMethod("get.array.subset","exprSet",get.array.subset.exprset);

"get.fold.change.and.t.test" <- function(x,group,members,logged = TRUE, a.order=NULL,b.order=NULL,method=c("unlogged","logged","median")) {
  
  a.samples <- exprs(get.array.subset(x,group,members[1]));
  b.samples <- exprs(get.array.subset(x,group,members[2]));
 
  pw <- FALSE;

  if(!is.null(a.order)) { 
    a.samples <- a.samples[,a.order];
    if(!is.null(b.order)) { 
      b.samples <- b.samples[,b.order]; 
      pw <- TRUE;
 
    }
    else {
     stop("Both a.order and b.order must be specified for a paired t-test");
    }
  }
  method <- match.arg(method)
  m <- switch(method,
              logged   = 2,
              unlogged = 1,
              median   = 3);

  a.samples.array <- as.double(t(a.samples));
  b.samples.array <- as.double(t(b.samples));

  nacol <- as.integer(length(a.samples[1,]));
  ngene <- as.integer(length(a.samples[,1]));
  nbcol <- as.integer(length(b.samples[1,]));

  if(class(logged) != "logical") stop("Parameter 'logged' should be TRUE or FALSE")
  if((nacol == 1) | (nbcol == 1))  warning("There was only one sample in one (or both) of your sample groups. Not computing t-tests - instead, returning 0.0 for p-scores...");

  c.res <- .C("FCMandTT",a.samples.array,b.samples.array,nacol,nbcol,ngene,as.logical(logged),pw,as.integer(m),ma = double(ngene),mb = double(ngene),fc = double(ngene),tt = double(ngene),PACKAGE="simpleaffy")
  means <- cbind(c.res$ma,c.res$mb);
  colnames(means) <- members;
  fc <- c.res$fc;
  tt <- c.res$tt;
  names(fc)       <- rownames(a.samples);
  rownames(means) <- rownames(a.samples);
  names(tt)       <- rownames(a.samples);
  
  return(new("PairComp",fc=fc,tt=tt,means=means,group=group,members=members))
}





"pairwise.comparison" <- function(x,group,members=NULL,spots=NULL,a.order=NULL,b.order=NULL,method="unlogged",logged=TRUE) {
  if(is.null(members)) {
    pd <- unique(as.character(pData(x)[,group]))
    if(is.null(pd)) {
      stop(paste("Can't find a group called",group));
    }      
    if(length(pd) != 2)  {
      stop("There must be exactly two groups for a pairwise comparison. Please specify which groups you want to compare.");
    }
    members <- pd;
  }
  if(!is.null(spots)) {
    pmac <- detection.p.val(spots);
    results <- get.fold.change.and.t.test(x,group,members,logged=logged,a.order=a.order,b.order=b.order,method=method);
    results@calls <- pmac$call;
  }
  else {
    results <- get.fold.change.and.t.test(x,group,members,logged=logged,a.order=a.order,b.order=b.order,method=method);
  }
  return(results); 
}


pairwise.filter <- function(object,x,min.exp=log2(100),min.exp.no=0,min.present.no=3,fc=1.0,tt=0.001) {

  if(class(object) != "PairComp") { stop("Can only filter an object of class 'PairComp'"); }
  if(class(x) != "exprSet") { stop("Can only filter using class 'exprSet' for parameter 'x'"); }
  pass.fc              <- (abs(fc(object)) > fc);
  pass.tt              <- (tt(object) < tt );

  samples <- exprs(get.array.subset(x,group(object),members(object)));

  no.chips             <- length(colnames(samples));

  intensity.thresh     <- array(sapply(samples[,1:no.chips],function(x) { if(x > min.exp) { 1 } else { 0 } } ),dim=dim(samples));

  min.intensity.thresh <- rowSums(intensity.thresh)

  
  pass.intensity <- (min.intensity.thresh >= min.exp.no);

  if(nrow(calls(object))>0) {
    present.thresh      <- array(sapply(calls(object)[,1:no.chips],function(x) { if(x == "P") { 1 } else { 0 } } ),dim=dim(calls(object)));
    min.present.thresh  <- rowSums(present.thresh);
    pass.present   <- (min.present.thresh >= min.present.no);

    return(object[(pass.fc & pass.tt & pass.intensity & pass.present),]);
  }
  else {
    return(object[(pass.fc & pass.tt & pass.intensity),]);
  }
}
library("methods")

#holds the results of a pairwise comparison
setClass("QCStats",representation(scale.factors="numeric",target="numeric",percent.present="numeric",average.background="numeric",minimum.background="numeric",maximum.background="numeric",bioB="character",bioC="character",bioD="character",creX="character",gapdh3="numeric",gapdhM="numeric",gapdh5="numeric",actin3="numeric",actinM="numeric",actin5="numeric"));

#accessor methods
setGeneric("sfs", function(object) standardGeneric("sfs"))
setMethod("sfs","QCStats",function(object) object@scale.factors)

setGeneric("target", function(object) standardGeneric("target"))
setMethod("target","QCStats",function(object) object@target)

setGeneric("percent.present", function(object) standardGeneric("percent.present"))
setMethod("percent.present","QCStats",function(object) object@percent.present)

setGeneric("avbg", function(object) standardGeneric("avbg"))
setMethod("avbg","QCStats",function(object) object@average.background)

setGeneric("minbg", function(object) standardGeneric("minbg"))
setMethod("minbg","QCStats",function(object) object@minimum.background)

setGeneric("maxbg", function(object) standardGeneric("maxbg"))
setMethod("maxbg","QCStats",function(object) object@maximum.background)

setGeneric("bioB", function(object) standardGeneric("bioB"))
setMethod("bioB","QCStats",function(object) object@bioB)

setGeneric("bioC", function(object) standardGeneric("bioC"))
setMethod("bioC","QCStats",function(object) object@bioC)

setGeneric("bioD", function(object) standardGeneric("bioD"))
setMethod("bioD","QCStats",function(object) object@bioD)

setGeneric("creX", function(object) standardGeneric("creX"))
setMethod("creX","QCStats",function(object) object@creX)

setGeneric("gapdh3", function(object) standardGeneric("gapdh3"))
setMethod("gapdh3","QCStats",function(object) object@gapdh3)

setGeneric("gapdhM", function(object) standardGeneric("gapdhM"))
setMethod("gapdhM","QCStats",function(object) object@gapdhM)

setGeneric("gapdh5", function(object) standardGeneric("gapdh5"))
setMethod("gapdh5","QCStats",function(object) object@gapdh5)

setGeneric("actin5", function(object) standardGeneric("actin5"))
setMethod("actin5","QCStats",function(object) object@actin5)

setGeneric("actinM", function(object) standardGeneric("actinM"))
setMethod("actinM","QCStats",function(object) object@actinM)

setGeneric("actin3", function(object) standardGeneric("actin3"))
setMethod("actin3","QCStats",function(object) object@actin3)

setGeneric("gapdh35", function(object) standardGeneric("gapdh35"))
setMethod("gapdh35","QCStats",function(object) object@gapdh3 - object@gapdh5);

setGeneric("actin35", function(object) standardGeneric("actin35"))
setMethod("actin35","QCStats",function(object) object@actin3 - object@actin5);

setGeneric("gapdh3M", function(object) standardGeneric("gapdh3M"))
setMethod("gapdh3M","QCStats",function(object) object@gapdh3 - object@gapdhM);

setGeneric("actin3M", function(object) standardGeneric("actin3M"))
setMethod("actin3M","QCStats",function(object) object@actin3 - object@actinM);

.nms <- c("atgenomecdf",	
"ath1121501cdf",	
"drosgenome1cdf",	
"hcg110cdf",		
"hgu133acdf",		
"hgu133a2cdf",	
"hgu133atagcdf",	
"hgu133bcdf",		
"hgu133plus2cdf",	
"hgu95acdf",		
"hgu95av2cdf",	
"hgu95bcdf",		
"hgu95ccdf",		
"hgu95dcdf",		
"hgu95ecdf",		
"mgu74acdf",		
"mgu74av2cdf",	
"mgu74bcdf",		
"mgu74bv2cdf",	
"mgu74ccdf",		
"mgu74cv2cdf",	
"mouse4302cdf",	
"mouse430acdf",		
"mouse430a2cdf",	
"mouse430b2cdf",		
"mu11ksubacdf",	
"mu11ksubbcdf",	
"rae230acdf",		
"rae230bcdf",		
"rgu34acdf",		
"rgu34bcdf",		
"rgu34ccdf")		

.alpha1 <- c(0.04,	
	     0.05,	
	     0.04,	
	     0.04,	
	     0.05,	
	     0.05,	
	     0.05,	
	     0.05,	
	     0.05,	
	     0.04,	
	     0.04,	
	     0.04,	
	     0.04,	
	     0.04,	
	     0.04,	
	     0.04,	
	     0.04,	
	     0.04,	
	     0.04,	
	     0.04,	
	     0.04,	
	     0.05,	
	     0.05,	
	     0.05,	
	     0.05,	
	     0.04,	
	     0.04,	
	     0.05,	
	     0.05,	
	     0.04,	
	     0.04,	
	     0.04)	

.alpha2 <- c(0.06,	
	     0.065,	
	     0.06,	
	     0.06,	
	     0.065,	
	     0.065,	
	     0.065,	
	     0.065,	
	     0.065,	
	     0.06,	
	     0.06,	
	     0.06,	
	     0.06,	
	     0.06,	
	     0.06,	
	     0.06,	
	     0.06,	
	     0.06,	
	     0.06,	
	     0.06,	
	     0.06,	
	     0.065,	
	     0.065,	
	     0.065,	
	     0.065,	
	     0.06,	
	     0.06,	
	     0.065,	
	     0.065,	
	     0.06,	
	     0.06,	
	     0.06)

.actin3 <- c("AFFX-r2-At-Actin-3_s_at",		
	     "AFFX-r2-At-Actin-3_s_at",		
	     "AFFX-Dros-ACTIN_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-HSAC07/X00351_3_at",		
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX-b-ActinMur/M12481_3_at",	
	     "AFFX_Rat_beta-actin_3_at",	
	     "AFFX_Rat_beta-actin_3_at",	
	     "AFFX_Rat_beta-actin_3_at",	
	     "AFFX_Rat_beta-actin_3_at",	
	     "AFFX_Rat_beta-actin_3_at")


.actinM <- c("AFFX-r2-At-Actin-M_s_at",		
	     "AFFX-r2-At-Actin-M_s_at",		
	     "AFFX-Dros-ACTIN_M_r_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-HSAC07/X00351_M_at",		
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX-b-ActinMur/M12481_M_at",	
	     "AFFX_Rat_beta-actin_M_at",	
	     "AFFX_Rat_beta-actin_M_at",	
	     "AFFX_Rat_beta-actin_M_at",	
	     "AFFX_Rat_beta-actin_M_at",	
	     "AFFX_Rat_beta-actin_M_at")


.actin5 <- c("AFFX-r2-At-Actin-5_s_at",		
	     "AFFX-r2-At-Actin-5_s_at",		
	     "AFFX-Dros-ACTIN_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-HSAC07/X00351_5_at",		
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX-b-ActinMur/M12481_5_at",	
	     "AFFX_Rat_beta-actin_5_at",	
	     "AFFX_Rat_beta-actin_5_at",	
	     "AFFX_Rat_beta-actin_5_at",	
	     "AFFX_Rat_beta-actin_5_at",	
	     "AFFX_Rat_beta-actin_5_at")

.gapdh3 <- c("AFFX-Athal-GAPDH_3_s_at",		
	     "AFFX-Athal-GAPDH_3_s_at",		
	     "AFFX-Dros-GAPDH_3_at",		
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-HUMGAPDH/M33197_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX-GapdhMur/M32599_3_at",	
	     "AFFX_Rat_GAPDH_3_at",		
	     "AFFX_Rat_GAPDH_3_at",		
	     "AFFX_Rat_GAPDH_3_at",		
	     "AFFX_Rat_GAPDH_3_at",		
	     "AFFX_Rat_GAPDH_3_at")		

.gapdhM <- c("AFFX-Athal-GAPDH_M_s_at",		
	     "AFFX-Athal-GAPDH_M_s_at",		
	     "AFFX-Dros-GAPDH_M_at",		
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-HUMGAPDH/M33197_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX-GapdhMur/M32599_M_at",	
	     "AFFX_Rat_GAPDH_M_at",		
	     "AFFX_Rat_GAPDH_M_at",		
	     "AFFX_Rat_GAPDH_M_at",		
	     "AFFX_Rat_GAPDH_M_at",		
	     "AFFX_Rat_GAPDH_M_at")		

.gapdh5 <- c("AFFX-Athal-GAPDH_5_s_at",		
	     "AFFX-Athal-GAPDH_5_s_at",		
	     "AFFX-Dros-GAPDH_5_at",		
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-HUMGAPDH/M33197_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX-GapdhMur/M32599_5_at",	
	     "AFFX_Rat_GAPDH_5_at",		
	     "AFFX_Rat_GAPDH_5_at",		
	     "AFFX_Rat_GAPDH_5_at",		
	     "AFFX_Rat_GAPDH_5_at",		
	     "AFFX_Rat_GAPDH_5_at")		

.biob <- c("AFFX-BioB-3_at",      
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-r2-Ec-bioB-3_at",
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at",	    
	    "AFFX-BioB-3_at")	  

.bioc <- c("AFFX-BioC-3_at",      
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-r2-Ec-bioC-3_at",
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at",	    
	    "AFFX-BioC-3_at")	    

.biod <- c("AFFX-BioD-3_at",      
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioD-3_at",	    
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-r2-Ec-bioD-3_at",
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at",	    
	    "AFFX-BioDn-3_at")	    

.crex <- c("AFFX-CreX-3_at",      
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-r2-P1-cre-3_at",
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at",	    
	    "AFFX-CreX-3_at")	    

names(.alpha1) <- .nms
names(.alpha2) <- .nms

names(.actin3) <- .nms
names(.actinM) <- .nms
names(.actin5) <- .nms

names(.gapdh3) <- .nms
names(.gapdhM) <- .nms
names(.gapdh5) <- .nms

names(.biob) <- .nms
names(.bioc) <- .nms
names(.biod) <- .nms
names(.crex) <- .nms


getTao <- function(name) {
  0.015;
}

getAlpha1 <- function(name) {
  .alpha1[name]
}

getAlpha2 <- function(name) {
  .alpha2[name]
}

getActin3 <- function(name) {
  .actin3[name]
}

getActinM <- function(name) {
  .actinM[name]
}

getActin5 <- function(name) {
  .actin5[name]
}

getGapdh3 <- function(name) {
  .gapdh3[name]
}

getGapdhM <- function(name) {
  .gapdhM[name]
}

getGapdh5 <- function(name) {
  .gapdh5[name]
}

getBioB <- function(name) {
  .biob[name]
}


getBioC <- function(name) {
  .bioc[name]
}


getBioD <- function(name) {
  .biod[name]
}


getCreX <- function(name) {
  .crex[name]
}




qc.affy <-function(unnormalised,normalised=NULL,logged=TRUE,
                   tau=getTao(cleancdfname(cdfName(unnormalised))),
                   alpha1=getAlpha1(cleancdfname(cdfName(unnormalised))),
                   alpha2=getAlpha2(cleancdfname(cdfName(unnormalised))),
	           bioB=getBioB(cleancdfname(cdfName(unnormalised))),
                   bioC=getBioC(cleancdfname(cdfName(unnormalised))),
                   bioD=getBioD(cleancdfname(cdfName(unnormalised))),
                   creX=getCreX(cleancdfname(cdfName(unnormalised))),
                   gapdh3=getGapdh3(cleancdfname(cdfName(unnormalised))),
                   gapdhM=getGapdhM(cleancdfname(cdfName(unnormalised))),
                   gapdh5=getGapdh5(cleancdfname(cdfName(unnormalised))),
                   actin3=getActin3(cleancdfname(cdfName(unnormalised))),
                   actinM=getActinM(cleancdfname(cdfName(unnormalised))),
                   actin5=getActin5(cleancdfname(cdfName(unnormalised)))) {
  if(is.null(normalised)) {
    normalised <- call.exprs(unnormalised,"mas5");
  }


  if(!is.element(cleancdfname(cdfName(unnormalised)),.nms)) {
	stop(paste("I'm sorry I do not know about chip type:",cleancdfname(cdfName(unnormalised))))
  }

  x <- exprs(normalised);

  det <- detection.p.val(unnormalised,tau=tau,alpha1=alpha1,alpha2=alpha2);

  dpv<-apply(det$call,2,function(x) { 
            x[x!="P"] <- 0;
	    x[x=="P"] <- 1;
            x<-as.numeric(x);			       
            return(100 * sum(x)/length(x));
  });

  sfs    <- normalised@description@preprocessing$sfs;
  target <- normalised@description@preprocessing$tgt;

  if(!logged) { x <- log2(x); }

  bgsts<-.bg.stats(unnormalised)$zonebg

  meanbg <- apply(bgsts,1,mean);

  minbg  <- apply(bgsts,1,min);

  maxbg  <- apply(bgsts,1,max);

  stdvbg <- sqrt(apply(bgsts,1,var));

  gapdh3 <- x[gapdh3,];
  gapdhM <- x[gapdhM,];
  gapdh5 <- x[gapdh5,];

  actin3 <- x[actin3,];
  actinM <- x[actinM,];
  actin5 <- x[actin5,];

  bioB <- det$call[bioB,];
  bioC <- det$call[bioC,];
  bioD <- det$call[bioD,];
  creX <- det$call[creX,];

  return(new("QCStats",scale.factors=sfs,target=target,percent.present=dpv,average.background=meanbg,minimum.background=minbg,maximum.background=maxbg,
              gapdh3=gapdh3,gapdhM=gapdhM,gapdh5=gapdh5,actin3=actin3,actinM=actinM,actin5=actin5,bioB=bioB,bioC=bioC,bioD=bioD,creX=creX));
}


setGeneric("qc", function(unnormalised,...) standardGeneric("qc"))
setMethod("qc","AffyBatch",function(unnormalised,...) qc.affy(unnormalised,...)) 

setGeneric("getQCParams", function(x) standardGeneric("getQCParams") )

setMethod("getQCParams","AffyBatch",function(x) {
  n <- cleancdfname(cdfName(x))

  res <- list(getGapdh3(n),getGapdhM(n),getGapdh5(n),getActin3(n),getActinM(n),getActin5(n),getBioB(n),getBioC(n),getBioD(n),getCreX(n),getAlpha1(n),getAlpha2(n),getTao(n))
  names(res) <- c("Gapdh3","GapdhM","Gapdh5","Actin3","ActinM","Actin5","BioB","BioC","BioD","CreX","Alpha1","Alpha2","Tao")
  return(res)
})


.bg.stats <- function(unnormalised, grid=c(4,4)) {
pms         <- unlist(pmindex(unnormalised))
mms         <- unlist(mmindex(unnormalised))
all         <- c(pms,mms)
intensities <- exprs(unnormalised)
rws <- nrow(unnormalised)
cls <- ncol(unnormalised)
zonebg <- c();
zonesd <- c();
for(no in 1:length(unnormalised)){
  this.array <- intensities[,no];
  result <- .C("bgmas",as.integer(as.vector(all)),as.integer(length(all)),
       as.double(as.vector(this.array)),as.integer(length(this.array)),
       as.integer(rws),
       as.integer(cls),
       as.integer(grid[1]),as.integer(grid[2]),
       zonebg=double(grid[1] * grid[2]),
       zonesd=double(grid[1] * grid[2]),corrected=double(length(this.array)),PACKAGE="simpleaffy");
  zonesd <- rbind(zonesd, result$zonesd);
  zonebg <- rbind(zonebg, result$zonebg);
  }
  colnames(zonesd) <- paste("zone",1:16,sep=".");
  colnames(zonebg) <- paste("zone",1:16,sep=".");
  rownames(zonesd) <- sampleNames(unnormalised);
  rownames(zonebg) <- sampleNames(unnormalised);
  return(list(zonebg=zonebg,zonesd=zonesd))
}



plot.qc.stats<-function(x,fc.line.col="black",sf.ok.region="light blue",chip.label.col="black",sf.thresh = 3.0,gdh.thresh = 3.0,ba.thresh = 3.0,present.thresh=10,bg.thresh=20,label=NULL,...) {

  sfs    <- log2(sfs(x))

  n      <- length(sfs)

  meansf <- mean(sfs)

  dpv <- percent.present(x)
  dpv <- (round(100*dpv))/100;

  abg <- avbg(x)
  abg <- (round(100*abg))/100;
	
  sfs <- sfs + 6.0;
  if(is.null(label)) { label <- 1:n }
  col=c("red","green");
  d1 <- 0.0;
  d2 <- 0.0;
  d3 <- 0.0;

  for(i in 1:n) {
    for(j in 1:n) { 
      d1 <- max(abs(sfs[i] - sfs[j]),d1);
      d2 <- max(abs(dpv[i] - dpv[j]),d1);
      d3 <- max(abs(abg[i] - abg[j]),d3);
    }
  }
  i<-1:101;
  arc <- 2 * pi / 100;
  x1 <- (7.5 +  meansf) * cos(i * arc);
  y1 <- (7.5 +  meansf) * sin(i * arc);
  x2 <- (4.5 +  meansf) * cos(i * arc);
  y2 <- (4.5 +  meansf) * sin(i * arc);

  plot(x1[1],y1[1],pch=".",xlim=range(-12,12),ylim=range(-12,12),xaxt="n",yaxt="n",xlab="",ylab="",col=sf.ok.region,...);

  polygon(c(x1,rev(x2),x1[1]),c(y1,rev(y2),y1[1]),col=sf.ok.region,border=sf.ok.region);

  x1 <- 3 * cos(i * arc);
  y1 <- 3 * sin(i * arc);
  points(x1,y1,pch=".",col=fc.line.col);
  x1 <- 6 * cos(i * arc);
  y1 <- 6 * sin(i * arc);
  points(x1,y1,pch=".",col=fc.line.col);
  x1 <- 9 * cos(i * arc);
  y1 <- 9 * sin(i * arc);
  points(x1,y1,pch=".",col=fc.line.col);
  text(0,3,"-3",col=fc.line.col)
  text(0,6,"0",col=fc.line.col)
  text(0,9,"+3",col=fc.line.col)
  arc <- 2 * pi / n;

  gdh <- gapdh35(x);
  ba  <- actin35(x)
  bb  <- bioB(x)

  for(i in 1:n) {
    if(d1 > sf.thresh) { col = "red" } else {col="blue"}
     x1 <- sfs[i] * cos(i * arc);
     y1 <- sfs[i] * sin(i * arc);
     x2 <- 6 * cos(i * arc);
     y2 <- 6 * sin(i * arc);
     lines(c(x2,x1),c(y2,y1),col=col);
     points(x1,y1,col=col,pch=20);
     text(x1,y1,col=chip.label.col,label=label[i],adj=0.2 * c(cos(i * arc),sin(i * arc)));
     x2 <- (6 + gdh[i]) * cos(i * arc);
     y2 <- (6 + gdh[i]) * sin(i * arc);
     if(gdh[i] > gdh.thresh) { col = "red" } else {col="blue"}	
     points(x2,y2,pch=1,col=col);
     x2 <- (6 + ba[i]) * cos(i * arc);
     y2 <- (6 + ba[i]) * sin(i * arc);
     if(ba[i] > ba.thresh) { col = "red" } else {col="blue"}	
     points(x2,y2,pch=2,col=col);

     if(d2 > present.thresh) { col = "red" } else {col="blue"}
     x2 <- (9 * cos(i * arc));
     y2 <- (9 * sin(i * arc));
     text(x2,y2,label=paste(dpv[i],"%",sep=""),col=col);

     if(d3 > bg.thresh) { col = "red" } else {col="blue"}
     x2 <- (10 * cos(i * arc));
     y2 <- (10 * sin(i * arc));
     text(x2,y2,label=abg[i],col=col);

     if(bb[i]!="P") {
       x2 <- (11 * cos(i * arc));
       y2 <- (11 * sin(i * arc));
       text(x2,y2,label="bioB",col="red");
     }
  }
  legend(-10,10,pch=1:2,c("GAPDH","Beta Actin"))
}

setMethod("plot","QCStats",function(x,y) plot.qc.stats(x,...))
"read.affy" <-
function(covdesc="covdesc",path=".",...) {
  samples <- read.phenoData( paste(path,covdesc,sep="/"));
  files.to.read <- rownames(pData(samples));
  files.to.read <- paste(path,files.to.read,sep="/")
  eset <- ReadAffy(filenames=files.to.read,...);
  newPhenoData <- cbind(pData(eset),pData(samples)[rownames(pData(eset)),]);
  colnames(newPhenoData) <- c(colnames(pData(eset)),colnames(pData(samples)));
  tmp <- as.list(colnames(newPhenoData));
  names(tmp) <- colnames(newPhenoData);
  newPhenoData <- new("phenoData",pData=newPhenoData,varLabels=tmp)

  eset@phenoData <- newPhenoData;	
  return(eset);
}
"trad.scatter.plot" <-
function(x,y,add=FALSE,fc.lines=log2(c(2,4,6,8)),draw.fc.lines=TRUE,draw.fc.line.labels=TRUE,fc.line.col="lightgrey",pch=20,...) {
   if(!add) {
     mx <- max(x,y);     
     mn <- min(x,y);
     plot(x,y,xlim=range(mn,mx),ylim=range(mn,mx),pch=pch,...);
     if(draw.fc.lines) {
       for( lne in fc.lines) {
         xc <- c(mn,mx);
         yc <- c(mn+lne,mx+lne);
         lines(xc,yc,col=fc.line.col); 
         if(draw.fc.line.labels) {text(mn+0.25,lne+mn,2^lne,col=fc.line.col); }
         xc <- c(mn, mx);
         yc <- c(mn-lne,mx-lne);
         lines(xc,yc,col=fc.line.col); 
          if(draw.fc.line.labels) {text(lne+mn,mn+0.25,2^lne,col=fc.line.col); }
       }
     }
   }
   else {
     points(x,y,pch=pch,...); 
   }
}
.First.lib <-  function(lib,pkg,where) {
  library.dynam("simpleaffy",pkg,lib);
  require(affy,quietly=TRUE);
  require(methods,quietly=TRUE);
  where <- match(paste("package:", pkg, sep=""), search());
  .initClassesAndMethods(where);
  cacheMetaData(as.environment(where));
  cat("Welcome to 'simpleaffy'                                      \n");
  cat("      Produced by The Paterson Institute for Cancer Research \n");
  cat("      and funded by CANCER RESEARCH UK.                      \n");
  cat("      http://bioinformatics.picr.man.ac.uk/simpleaffy        \n");
  cat("      mailto: microarray@picr.man.ac.uk                      \n");
}

.initClassesAndMethods <- function(where) {

}
