"fixup.nom" <-
function(nom, fact, cvrt, str.cross=FALSE)
{
  # nom:  a nomogram object as returned by 'nomogram' in Harrell's rms package 
  #
  # fact: (character, currently length 1) use this argument to specify the name of
  #       a factor upon which there is partial interaction (some but not all of the levels)
  #       at present we handle only two way partial interactions betweeen a factor and
  #       ordinal covariate(s).
  #
  # cvrt: (a named list of character vectors) use this argument to specify each covariate that
  #       the above mentioned factor has partial interaction with (the name of each component)
  #       and the levels upon which there is interaction (the character strings in each component)
  #       currently interaction covariates are assumed to be categorical, with a referent category
  #       equal to zero.
  #
  # example call
  #
  # nom.new <- fixup.nom(nom=nom.old, fact="treatmt.", cvrt=list(stage=c("horm","surg"), HIcomorbs=c("horm","rad")))
  #
  "%,%" <- function(x,y)paste(x,y,sep="")

  #------------------------------------------------
  # first, get rid of the 'HIcomorbs (arm=1)' scale
  # this is specific to my problem --I also have an interaction with
  # the stratification variable, 'arm' (i.e. separate baseline hazards
  # per arm) and I want to do all predictions for the control arm.
  # and then store this original nomogram and its attributes 
  # note: see below where I need to subtract out the arm=1 
  

  nom.sv.nms <- names(nom)
  nom.sv <- nom
  nom.nms <- names(nom)
  if(str.cross)
  {
    ind.str.dr <- which(nom.sv.nms=="HIcomorbs (arm=arm=1)")
    ind.str.kp <- which(nom.sv.nms=="HIcomorbs (arm=arm=0)")
    nom[[ind.str.dr]]<- NULL
    nom.nms <- nom.nms[-ind.str.dr]
    nom.nms[ind.str.kp] <- "HIcomorbs"
  }
  names(nom) <- nom.nms.sv.r <- nom.nms
  nom.sv.r <- nom
  n.nom.sv.r <- length(nom.sv.r)
  atr.nom.sv.r <- attributes(nom.sv.r)
 
  # initialize
  fact.nm <- fact
  fact.lvls <- nom[[fact.nm]][[fact.nm]]
  cross.with <- names(cvrt)
  n.cross.with <- length(cross.with)
  cross.with.lvls <- list()
  n.cross.with.lvls <- rep(0, n.cross.with)
  eg.call <- list()
  eg.call[[1]] <- as.name("expand.grid")
  
  for(k in 1:n.cross.with)
  {
    cross.with.lvls[[k]] <- nom[[cross.with[k]]][[cross.with[k]]]
    eg.call[[1+k]] <- cross.with.lvls[[k]]
    n.cross.with.lvls[k] <- length(cross.with.lvls[[k]])
  }
  names(cross.with.lvls) <- names(n.cross.with.lvls) <- cross.with
  eg.call <- as.call(eg.call)
  eg <- eval(eg.call)
  names(eg) <- cross.with
  
  n.axes <- prod(n.cross.with.lvls)
  
  # initialize
  #
  # original code:
  #---------------------------------------------------------------------
  # treatmt.00 <- treatmt.01 <- treatmt.10 <- treatmt.11 <- nom$treatmt.
  #----------------------------------------------------------------------
  
  # make a new scale for the factor variable for each level of the partial interaction.  First add
  # interaction category specific intercept. Then add the partial interaction specific component
  #
  # original code:
  #------------------------------------------------------------------
  # treatmt.01$Xbeta <- nom$treatmt.$Xbeta + nom$stage$Xbeta[which(nom$stage$Xbeta!=0)]
  #  treatmt.01$Xbeta[treatmt.01$treatmt.=="horm"] <-
  #    treatmt.01$Xbeta[treatmt.01$treatmt.=="horm"] + nom$`horm.x.stage`$Xbeta[which(nom$`horm.x.stage`$Xbeta!=0)]

  #  treatmt.10$Xbeta <- nom$treatmt.$Xbeta + nom$HIcomorbs$Xbeta[which(nom$HIcomorbs$Xbeta!=0)]
  #  treatmt.10$Xbeta[treatmt.10$treatmt.=="horm"] <-
  #    treatmt.10$Xbeta[treatmt.10$treatmt.=="horm"] + nom$`horm.x.HIcomorbs`$Xbeta[which(nom$`horm.x.HIcomorbs`$Xbeta!=0)]

  #  treatmt.11$Xbeta <- nom$treatmt.$Xbeta + nom$stage$Xbeta[which(nom$stage$Xbeta!=0)] +
  #                        nom$HIcomorbs$Xbeta[which(nom$HIcomorbs$Xbeta!=0)]
  #  treatmt.11$Xbeta[treatmt.11$treatmt.=="horm"] <-
  #   treatmt.11$Xbeta[treatmt.11$treatmt.=="horm"] + nom$`horm.x.stage`$Xbeta[which(nom$`horm.x.stage`$Xbeta!=0)]
  #      nom$`horm.x.HIcomorbs`$Xbeta[which(nom$`horm.x.HIcomorbs`$Xbeta!=0)]

  fact.axes <- list()
  nms.fact.axes <- rep("", n.axes)
  for(k in 1:n.axes)
  {
    nms.fact.axes[k] <- fact.nm %,% "("
    
    # e.g. cw.vals.k <- c(stage=1, HIcomorbs=1)
    cw.vals.k <- eg[k,]
    # initial assignment takes care of 'fact' main effects
    fact.axes[[k]] <- nom[[fact.nm]]

    .sep. <- c(rep(", ", n.cross.with-1), ")")
    old.cross.nms <- NULL
    for(j in 1:n.cross.with)
    {
      # e.g. cvrt.j <- "stage"
      cvrt.j <- cross.with[j]

      if(cw.vals.k[j]!=0)
      {
        # now we add the main effects of the 'cvrt's with which there is interaction
        fact.axes[[k]]$Xbeta <- fact.axes[[k]]$Xbeta + nom[[cvrt.j]]$Xbeta[which(nom[[cvrt.j]][[cvrt.j]]==cw.vals.k[j])]

        # now we look for interaction terms
        # e.g. cross.nms.j <- c("stage.x.horm", "stage.x.surg")
        cross.nms.j <- cvrt.j %,% ".x." %,% cvrt[[cvrt.j]]
        n.cross.nms.j <- length(cross.nms.j)
        for(l in 1:n.cross.nms.j)
        {
          # e.g. cross.nm.jl <- "stage.x.horm"
          cross.nm.jl <- cross.nms.j[l]
          old.cross.nms <- c(old.cross.nms, cross.nm.jl)
          # e.g. fact.lvl.jl <- "horm"
          fact.lvl.jl <- cvrt[[cvrt.j]][l]
          fact.axes[[k]]$Xbeta[which(fact.lvls==fact.lvl.jl)] <- fact.axes[[k]]$Xbeta[which(fact.lvls==fact.lvl.jl)] +
            nom[[cross.nm.jl]]$Xbeta[which(fact.lvls==fact.lvl.jl)]
        }
      }
      nms.fact.axes[k] <- nms.fact.axes[k] %,% cvrt.j %,% "=" %,% cw.vals.k[j] %,% .sep.[j]
    }
  }

  # Original Code:
  #---------------------------------
  # remove the main effect treatmt. axes and all of the variables it interacts with
  # and then in the following, replace with interaction specific treatmt. axes
  #  dr.nms <- c("treatmt.", "stage", "HIcomorbs", "horm.x.stage", "horm.x.HIcomorbs")
  #  min.dr.idx <- min(sapply(dr.nms, FUN=function(x, y)which(y==x), y=nms.nom))
  #
  #  n.nom <- length(nom)
  #  tmpl <- list()
  #  repl <- list(treatmt.00, treatmt.01, treatmt.10, treatmt.11)
  #  l <- 1
  #  for(k in 1:n.nom)
  #  {
  #    if(k<min.dr.idx) tmpl[[k]] <- nom[[k]]
  #    if(k>=min.dr.idx & k <=min.dr.idx+3) tmpl[[k]] <- repl[[k-min.dr.idx+1]]
  #    if(k>=min.dr.idx+4 & k<=min.dr.idx+5) tmpl[[k]] <- nom[[n.nom-1 + k-min.dr.idx-4]]
  # }

  # Now, do the same deletion and replacement for the names component of the attributes
  # and then replace 'nom' by 'tmpl' (can't replace until the names attribute is of the
  # right length
  #
  #  pos.start <- which(nom.atr.nms=="treatmt.")-1
  #  pos.tag <- which(nom.atr.nms=="total.points")
  #  tmpn <- nom.atr.nms[1:pos.start]
  #  tmpn <- c(tmpn, "treatmt.(HIcomorbs=0, stage=0)", "treatmt.(HIcomorbs=0, stage=1)",
  #                "treatmt.(HIcomorbs=1, stage=0)", "treatmt.(HIcomorbs=1, stage=1)")
  #  tmpn <- c(tmpn, nom.atr.nms[pos.tag:(pos.tag+1)])
  #  nom.atr.nms <- tmpn
  #
  #  attr(tmpl, "names") <- nom.atr.nms
  #  nom <- tmpl

  #  attr(nom, "info") <- atr.sv$info
  dr.nms <- c(fact.nm, cross.with, unique(old.cross.nms))
  n.nom <- length(nom)
  nom.axes.nms <- nom.nms[-c(which(nom.nms[1:(n.nom-2)] %in% dr.nms), (n.nom-1):n.nom)]
  nom.last2.nms <- nom.nms[(n.nom-1):n.nom]
  tmpl <- list()
  grouped <- NULL
  # repl <-> fact.axes
  # n.repl <-> n.axes
  l <- 1
  for(k in 1:(n.nom-2))
  {
    if(!(nom.nms[k] %in% dr.nms))
    {
      tmpl[[l]] <- nom[[k]]
      l <- l+1
    }
  }
  for(k in 1:n.axes)
  {
    tmpl[[l]] <- fact.axes[[k]]
    grouped <- c(grouped, l)
    l <- l+1
  }
  for(k in (n.nom-1):n.nom)
  {
    tmpl[[l]] <- nom[[k]]
    l <- l+1
  }
  new.nom.nms <- c(nom.axes.nms, nms.fact.axes, nom.last2.nms)

  names(tmpl) <- new.nom.nms
  nom <- tmpl


  # now we just make sure that all components of the 'info' attributes component have appropriate
  # entries corresponding to our four new axes:
  # 1. info$discrete

  attr(nom, "info") <- atr.nom.sv.r$info
  attr(nom, "info")$discrete[grouped] <- TRUE

  # locate the '0' and '100' of the 'points' scale
  # i.e. the smallest and largest value of Xbeta
  # and then use this to generate the new points for each scale  

  supp.xb.new <- sapply(nom, FUN=function(x)if(!is.null(x$Xbeta))range(x$Xbeta) else c(-Inf, Inf))
  supp.xb.new <- supp.xb.new[,-which(supp.xb.new[1,]==-Inf)]
  attr(nom, "info")$R <- supp.xb.new
  
  lower.limit <- supp.xb.new[1,]
  lower.limit[grouped] <- rep(min(lower.limit[grouped], length(grouped)))
  min.xb <- min(supp.xb.new[1, ])
  max.xb <- max(supp.xb.new[2, ])
  supp.xb <- c(min.xb, max.xb)
  rng.xb <- diff(supp.xb)

  supp.pts.old <- NULL
  n.old <-  sum(sapply(nom.sv.r, FUN=function(x)any(names(x)=="Xbeta")))
  for(k in 1:n.old)
  {
    supp.pts.old <- cbind(supp.pts.old, range(nom.sv.r[[k]]$points))
  }
  dimnames(supp.pts.old)[[2]] <- x.old.nms <- nom.nms.sv.r[1:(n.nom.sv.r-2)]
                                       
  n <- sum(sapply(nom, FUN=function(x)any(names(x)=="Xbeta")))
  supp.pts.new <- NULL
  for(k in 1:n)
  {
    nom[[k]]$points <- 100*(nom[[k]]$Xbeta - lower.limit[k])/rng.xb
    supp.pts.new <- cbind(supp.pts.new, range(nom[[k]]$points))
  }
  dimnames(supp.pts.new)[[2]] <- x.new.nms <- nom.nms[1:n]
  
  cmn <- intersect(x.old.nms, x.new.nms)
  cmn <- cmn[-which(cmn %in% c(cross.with, fact.nm))]

  .slope. <- unique(supp.pts.new[2,cmn]/supp.pts.old[2,cmn])[1]
  .int. <- supp.pts.new[1,cross.with[1]] - .slope. * supp.pts.old[1,cross.with[1]]
  nom.sv[["HIcomorbs (arm=arm=1)"]]$points[1]
  
  len.nom <- length(nom)

  x.new <- .slope.*nom[[len.nom]]$x + .int.
  ind <- order(x.new)
  nom[[len.nom]]$x <- x.new[ind]
  nom[[len.nom]]$x.real <- nom[[len.nom]]$x.real[ind]
  nom[[len.nom]]$fat <- nom[[len.nom]]$fat[ind]
  len.x <- length(nom[[len.nom]]$x)
  nom[[len.nom-1]]$x <- pretty(c(0, nom[[len.nom]]$x), n=len.x)
  class(nom) <- "nomogram"
  nom
}