PageRenderTime 5ms CodeModel.GetById 26ms app.highlight 14ms RepoModel.GetById 1ms app.codeStats 0ms

/rcdk/R/formula.R

http://github.com/rajarshi/cdkr
R | 575 lines | 326 code | 77 blank | 172 comment | 63 complexity | 414dee67dd6cc7a0b768ee0dbf67a145 MD5 | raw file
  1########################################################
  2##  set a cdkFormula function   
  3.IMolecularFormula <- "org/openscience/cdk/interfaces/IMolecularFormula"
  4
  5#' @name cdkFormula-class
  6#' @title  Class cdkFormula, ac class for handling molecular formula
  7#' @description This class handles molecular formulae. It provides extra 
  8#'  information such as the IMolecularFormula Java object, elements contained 
  9#'  and number of them.
 10#' @section Objects from the Class:
 11#'   Objects can be created using new constructor and filled with a specific 
 12#'   mass and window accuracy
 13#' @author Miguel Rojas-Cherto (\email{miguelrojasch@@yahoo.es})
 14#' @references A parallel effort to expand the Chemistry Development 
 15#'   Kit: \url{http://cdk.sourceforge.net}
 16#' @seealso  
 17#' \code{\link{get.formula}}
 18#' \code{\link{set.charge.formula}}
 19#' \code{\link{get.isotopes.pattern}}
 20#' \code{\link{isvalid.formula}}
 21#' @docType class
 22#' @aliases cdkFormula-class
 23#' @aliases show,cdkFormula-method
 24setClass("cdkFormula", representation(mass = "numeric",
 25                                      objectJ = "jobjRef",
 26                                      string = "character",
 27                                      charge = "numeric",
 28                                      isotopes = "matrix"),
 29         prototype(mass = 0,
 30                   objectJ = NULL,
 31                   string = character(0),
 32                   charge = 0,
 33                   isotopes = matrix(nrow = 0, ncol = 0))
 34         )
 35
 36#' get.formula
 37#' 
 38#' obtain molecular formula from formula string
 39#' 
 40#' @param mf Required. Molecular formula
 41#' @param charge Optional. Default \code{0}
 42#' @export
 43get.formula <- function(mf, charge=0) {
 44    
 45    manipulator <- get("mfManipulator", envir = .rcdk.GlobalEnv)
 46    if(!is.character(mf)) {
 47        stop("Must supply a Formula string");
 48    }else{
 49        dcob <- .cdkFormula.createChemObject()
 50        molecularformula <- .cdkFormula.createFormulaObject()
 51        molecularFormula <- .jcall(manipulator,
 52                                   "Lorg/openscience/cdk/interfaces/IMolecularFormula;",
 53                                   "getMolecularFormula",
 54                                   mf,
 55                                   .jcast(molecularformula,.IMolecularFormula),
 56                                   TRUE);
 57    }
 58    
 59    D <- new(J("java/lang/Integer"), as.integer(charge))
 60    .jcall(molecularFormula,"V","setCharge",D);
 61    object <- .cdkFormula.createObject(.jcast(molecularFormula,.IMolecularFormula));
 62    return(object);
 63}
 64
 65setMethod("show", "cdkFormula",
 66          function(object) {
 67              cat('cdkFormula: ',object@string,
 68                  ', mass = ',object@mass, ', charge = ',
 69                  object@charge,  '\n')
 70          })
 71#' 
 72#' get.mol2formula
 73#' 
 74#' 
 75#' @param molecule The molecule to query. Should be a `jobjRef` representing an `IAtomContainer`
 76#' @param charge Optional. Default \code{0}
 77#' @export
 78get.mol2formula <- function(molecule, charge=0) {
 79    if(attr(molecule, "jclass") != "org/openscience/cdk/interfaces/IAtomContainer") {
 80        stop("Must supply an IAtomContainerobject")
 81    }
 82    
 83    formulaJ <- .jcall('org/openscience/cdk/tools/manipulator/MolecularFormulaManipulator',
 84                       "Lorg/openscience/cdk/interfaces/IMolecularFormula;",
 85                       "getMolecularFormula",
 86                       molecule, use.true.class=FALSE);
 87    formulaJ <- .jcast(formulaJ,"org/openscience/cdk/interfaces/IMolecularFormula")
 88    
 89    ## needs that all isotopes contain the properties
 90
 91    string <- .cdkFormula.getString(formulaJ)
 92    objectF <- .cdkFormula.createFormulaObject()
 93    moleculaJT <- .jcall('org/openscience/cdk/tools/manipulator/MolecularFormulaManipulator',
 94                         "Lorg/openscience/cdk/interfaces/IMolecularFormula;",
 95                         "getMolecularFormula",string,
 96                         .jcast(objectF,"org/openscience/cdk/interfaces/IMolecularFormula"),TRUE);
 97    
 98    Do <- new(J("java/lang/Integer"), as.integer(charge))
 99    .jcall(moleculaJT,"V","setCharge",Do);	   
100
101    formula <- .cdkFormula.createObject(.jcast(moleculaJT,.IMolecularFormula))
102    return(formula);
103}
104#' 
105#' set.charge.formula
106#' 
107#' Set the charge to a cdkFormula function.
108#' 
109#' @param formula Required. Molecular formula
110#' @param charge Optional. Default \code{-1}
111#' @export
112set.charge.formula <- function(formula,charge = -1) {
113    if (class(formula) != "cdkFormula")
114        stop("Supplied object should be a cdkFormula Class")
115    
116    molecularFormula <- formula@objectJ;
117    
118    D <- new(J("java/lang/Integer"), as.integer(charge))
119    .jcall(molecularFormula,"V","setCharge",D);
120    
121    formula@objectJ <- molecularFormula;
122    formula@charge <- charge;
123    
124    return(formula)
125}
126
127#' 
128#' isvalid.formula
129#' 
130#' Validate a cdkFormula.
131#' 
132#' @param formula Required. A CDK Formula
133#' @param rule Optional. Default \code{rule=c("nitrogen","RDBE")}
134#' @export
135isvalid.formula <- function(formula,rule=c("nitrogen","RDBE")){
136    
137    if (class(formula) != "cdkFormula")
138        stop("Supplied object should be a cdkFormula Class")
139    
140    molecularFormula <- formula@objectJ;
141
142    nRule <- get("nRule", envir = .rcdk.GlobalEnv)
143    rdbeRule <- get("rdbeRule", envir = .rcdk.GlobalEnv)
144
145    for(i in 1:length(rule)){
146        ##Nitrogen Rule
147        if(rule[i] == "nitrogen"){
148            valid <- .jcall(nRule,"D","validate",molecularFormula);
149            if(valid != 1.0){
150                return (FALSE)
151            }
152        }	  
153        ##RDBE Rule
154        if(rule[i] == "RDBE"){
155            valid <- .jcall(rdbeRule,"D","validate",molecularFormula);
156            if(valid != 1.0){
157                return (FALSE)
158            }
159            else return(TRUE);
160        }
161    }
162    return(TRUE);
163}
164#' 
165#' get.isotopes.pattern
166#' 
167#' Generate the isotope pattern given a formula class
168#'
169#' @param formula Required. A CDK molecule formula
170#' @param minAbund Optional. Default \code{0.1}
171#' @export
172get.isotopes.pattern <- function(formula,minAbund=0.1){
173    
174    if (class(formula) != "cdkFormula")
175        stop("Supplied object should be a cdkFormula Class")
176    
177    molecularFormula <- formula@objectJ;
178    
179    isoGen <- .jnew("org/openscience/cdk/formula/IsotopePatternGenerator",as.double(minAbund));
180    isoPattern <- .jcall(isoGen,
181                         "Lorg/openscience/cdk/formula/IsotopePattern;",
182                         "getIsotopes",molecularFormula);
183    numIP <- .jcall(isoPattern,"I","getNumberOfIsotopes");
184    
185    ## create a matrix adding the mass and abundance of the isotope pattern
186    iso.col <- c("mass","abund");
187    
188    massVSabun <- matrix(ncol=2,nrow=numIP);
189    colnames(massVSabun)<-iso.col;
190    for (i in 1:numIP) {
191        isoContainer <- .jcall(isoPattern,"Lorg/openscience/cdk/formula/IsotopeContainer;","getIsotope",as.integer(i-1));
192        massVSabun[i,1] <- .jcall(isoContainer,"D","getMass");
193        massVSabun[i,2] <- .jcall(isoContainer,"D","getIntensity");
194    }
195    return (massVSabun);
196}
197#' 
198#' generate.formula.iter
199#' 
200#' Generate a list of possible formula objects given a mass and a mass tolerance.
201#' 
202#' @param mass Required. Mass.
203#' @param window Optional. Default \code{0.01}
204#' @param elements Optional. Default \code{
205#'                  list(c('C', 0,50), c('H', 0,50), 
206#'                  c('N', 0,50), c('O', 0,50), c('S', 0,50))}
207#' @param validation Optional. Default \code{FALSE}
208#' @param charge Optional. Default \code{FALSE}
209#' @param as.string Optional. Default \code{FALSE}
210#' @export
211generate.formula.iter <- function(mass, window = 0.01,
212                                  elements = list(
213                                      c('C', 0,50),
214                                      c('H', 0,50),
215                                      c('N', 0,50),
216                                      c('O', 0,50),
217                                      c('S', 0,50)),
218                                  validation = FALSE,
219                                  charge = 0.0,
220                                  as.string=TRUE) {
221
222    mfRange <- .jnew("org/openscience/cdk/formula/MolecularFormulaRange");
223    ifac <- .jcall("org/openscience/cdk/config/Isotopes",
224                   "Lorg/openscience/cdk/config/Isotopes;",
225                   "getInstance");
226    for (i in 1:length(elements)) {
227      
228      ## If the element list is 3, then we have sym, min, max
229      ## otherwise it should be sym, min, max, massNumber
230      if (length(elements[[i]]) == 3) {
231        isotope <- .jcall(ifac,
232                          "Lorg/openscience/cdk/interfaces/IIsotope;",
233                          "getMajorIsotope",
234                          as.character( elements[[i]][1] ), use.true.class=FALSE)
235      } else if(length(elements[[i]]) == 4) {
236        isotope <- .jcall(ifac,
237                          "Lorg/openscience/cdk/interfaces/IIsotope;",
238                          "getIsotope",
239                          as.character( elements[[i]][1] ),
240                          as.integer( elements[[i]][4] ),
241                          use.true.class=FALSE)
242        if (is.null(isotope)) stop(sprintf("Invalid mass number specified for element %s",elements[[i]][1]))
243      } else stop("Elements must be 3-tuples or 4-tuples")
244
245      .jcall(mfRange,
246             returnSig="V",
247             method="addIsotope",
248             isotope,
249             as.integer( elements[[i]][2] ),
250             as.integer( elements[[i]][3] )
251             )
252      
253    }
254
255    ## Construct range strings
256    ## rstrs <- sapply(names(elements), function(x) paste0(c(x, elements[[x]][1], elements[[x]][2]), sep='', collapse=' '))
257    ## if (length(rstrs) == 0)
258    ##     warning("The element specification resulted in a 0 length vector. This is worrisome")
259    
260    ## ## Get MF range object
261    ## mfRange <- .jcall("org/guha/rcdk/util/Misc",
262    ##                   "Lorg/openscience/cdk/formula/MolecularFormulaRange;",
263    ##                   "getMFRange",
264    ##                   .jarray(rstrs))
265    
266    ## construct generator
267    mfgen <- .jnew("org/openscience/cdk/formula/RoundRobinFormulaGenerator",
268                   get("dcob", envir = .rcdk.GlobalEnv),
269                   as.double(mass-window),
270                   as.double(mass+window),
271                   mfRange)
272
273    hasNext <- NA
274    formula <- NA
275    
276    ## hasNx <- function() {
277    ##   hasNext <<- .jcall(mfgen, "Lorg/openscience/cdk/interfaces/IMolecularFormula;", "getNextFormula")
278    ##   if (is.jnull(hasNext)) {
279    ##     ## nothing to do
280    ##     formula <<- NA
281    ##   } else formula <<- hasNext
282    ##   return(!is.jnull(hasNext))
283    ## }
284    
285    nextEl <- function() {
286        hasNext <<- NA
287        formula <- .jcall(mfgen, "Lorg/openscience/cdk/interfaces/IMolecularFormula;", "getNextFormula")
288        if (is.jnull(formula)) stop("StopIteration")
289        if (!as.string) {
290            return(formula)
291        } else {
292            return(.jcall("org/openscience/cdk/tools/manipulator/MolecularFormulaManipulator",
293                          "S", "getString", formula, FALSE, TRUE))
294        }
295    }
296
297    ##obj <- list(nextElem = nextEl, hasNext = hasNx)
298    obj <- list(nextElem = nextEl)  
299    class(obj) <- c("generate.formula2", "abstractiter", "iter")
300    return(obj)
301}
302#' generate.formula
303#' 
304#' @param mass Required. Mass.
305#' @param window Optional. Default \code{0.01}
306#' @param elements Optional. Default \code{
307#'                  list(c('C', 0,50), c('H', 0,50), 
308#'                  c('N', 0,50), c('O', 0,50), c('S', 0,50))}
309#' @param validation Optional. Default \code{FALSE}
310#' @param charge Optional. Default \code{FALSE}
311#' @export
312generate.formula <- function(mass, 
313                             window=0.01, 
314                             elements=list(c("C",0,50),c("H",0,50),c("N",0,50),c("O",0,50),c("S",0,50)), 
315                             validation=FALSE,
316                             charge=0.0){
317    
318    builder <- .cdkFormula.createChemObject();
319    mfTool <- .jnew("org/openscience/cdk/formula/MassToFormulaTool",builder);
320    ruleList <-.jcast(.jcall("org/guha/rcdk/formula/FormulaTools",
321                             "Ljava/util/List;",
322                             "createList"), "java/util/List")
323    
324    ## TOLERANCE RULE
325    toleranceRule <- .jnew("org/openscience/cdk/formula/rules/ToleranceRangeRule");
326    ruleG <- .jcast(toleranceRule, "org/openscience/cdk/formula/rules/IRule");
327    D <- new(J("java/lang/Double"), window)
328    paramsA <- .jarray(list(D,D))
329    paramsB <- .jcastToArray(paramsA)
330    .jcall(ruleG,"V","setParameters",paramsB);
331    ruleList <-.jcall("org/guha/rcdk/formula/FormulaTools", "Ljava/util/List;", "addTo",ruleList,ruleG)
332    
333    
334    ## ELEMENTS RULE
335    elementRule <- .jnew("org/openscience/cdk/formula/rules/ElementRule");
336    ruleG <- .jcast(elementRule, "org/openscience/cdk/formula/rules/IRule");
337    
338    chemObject <- .cdkFormula.createChemObject();
339    range <- .jnew("org/openscience/cdk/formula/MolecularFormulaRange");
340    ifac <- .jcall("org/openscience/cdk/config/Isotopes",
341                   "Lorg/openscience/cdk/config/Isotopes;",
342                   "getInstance");
343
344    for (i in 1:length(elements)) {
345
346        ## If the element list is 3, then we have sym, min, max
347        ## otherwise it should be sym, min, max, massNumber
348        if (length(elements[[i]]) == 3) {
349            isotope <- .jcall(ifac,
350                              "Lorg/openscience/cdk/interfaces/IIsotope;",
351                              "getMajorIsotope",
352                              as.character( elements[[i]][1] ), use.true.class=FALSE)
353        } else if(length(elements[[i]]) == 4) {
354            isotope <- .jcall(ifac,
355                              "Lorg/openscience/cdk/interfaces/IIsotope;",
356                              "getIsotope",
357                              as.character( elements[[i]][1] ),
358                              as.integer( elements[[i]][4] ),
359                              use.true.class=FALSE)
360            if (is.null(isotope)) stop(sprintf("Invalid mass number specified for element %s",elements[[i]][1]))
361        } else stop("Elements must be 3-tuples or 4-tuples")
362
363        .jcall(range,
364               returnSig="V",
365               method="addIsotope",
366               isotope,
367               as.integer( elements[[i]][2] ),
368               as.integer( elements[[i]][3] )
369               )
370
371    }
372    
373    paramsA <- .jarray(list(range))
374    paramsB <- .jcastToArray(paramsA)
375    .jcall(ruleG,"V","setParameters",paramsB);
376    
377    ruleList <-.jcall("org/guha/rcdk/formula/FormulaTools",
378                      "Ljava/util/List;", "addTo",
379                      .jcast(ruleList,"java/util/List"),
380                      ruleG)
381    
382    ## Setting the rules int FormulaTools
383    .jcall(mfTool,"V","setRestrictions",.jcast(ruleList,"java/util/List"));
384    
385    mfSet <- .jcall(mfTool,"Lorg/openscience/cdk/interfaces/IMolecularFormulaSet;",
386                    "generate",mass);
387    if (is.null(mfSet))
388        return(list())
389
390    sizeSet <- .jcall(mfSet,"I","size");
391    ecList <- list();
392    count = 1;
393    
394    for (i in 1:sizeSet) {
395        mf <- .jcall(mfSet,
396                     "Lorg/openscience/cdk/interfaces/IMolecularFormula;",
397                     "getMolecularFormula",
398                     as.integer(i-1));
399
400        .jcall(mf,"V","setCharge",new(J("java/lang/Integer"), as.integer(charge)));
401        object <- .cdkFormula.createObject(.jcast(mf,
402                                                  "org/openscience/cdk/interfaces/IMolecularFormula"));
403        
404        isValid = TRUE;
405        if(validation==TRUE)
406            isValid = isvalid.formula(object);
407        
408        if(isValid==TRUE){ ## if it's true add to the list
409            ecList[[count]] = object;
410            count = count+1;
411        }
412    }
413    ecList
414}
415
416
417#############################################################
418##  Intern functions: Creating object
419#############################################################
420
421.cdkFormula.createChemObject <- function(){
422    get("dcob", envir = .rcdk.GlobalEnv)
423                                        #  dcob <- .jcall("org/openscience/cdk/DefaultChemObjectBuilder",
424                                        #                 "Lorg/openscience/cdk/interfaces/IChemObjectBuilder;",
425                                        #                 "getInstance")
426                                        #  dcob
427}
428
429.cdkFormula.createFormulaObject <- function(){
430    ##   dcob <- .cdkFormula.createChemObject()
431    ##   klass <- J("org.openscience.cdk.interfaces.IMolecularFormula")$class
432    ##   cfob <- .jcall(dcob,
433    ##                  "Lorg/openscience/cdk/interfaces/ICDKObject;",
434    ##                  "newInstance",
435    ##                  klass
436    ##                  );  
437    ##   cfob
438    .jcast(.jnew("org/openscience/cdk/formula/MolecularFormula"),
439           "org/openscience/cdk/interfaces/IMolecularFormula")
440}
441
442#############################################################
443## extract the molecular formula string form the java object
444#############################################################
445.cdkFormula.getString <- function(molecularFormula) {
446    
447    if (attr(molecularFormula, "jclass") != 'org/openscience/cdk/interfaces/IMolecularFormula') {
448        stop("Supplied object should be a Java reference to an IMolecularFormula")
449    }
450    formula <- .jcall('org/openscience/cdk/tools/manipulator/MolecularFormulaManipulator',
451                      'S', 'getString', molecularFormula)
452}
453
454#############################################################
455## create a formula class from the molecularFormula java object
456#############################################################
457.cdkFormula.createObject <- function(molecularformula){
458    
459    object <-new("cdkFormula")
460    
461    object@objectJ <- molecularformula;
462    iterable <- .jcall(molecularformula,"Ljava/lang/Iterable;","isotopes"); 
463    isoIter <- .jcall(iterable,"Ljava/util/Iterator;","iterator");
464    size <- .jcall(molecularformula,"I","getIsotopeCount");
465    
466    isotopeList = matrix(ncol=3,nrow=size);
467    colnames(isotopeList) <- c("isoto","number","mass");
468    
469    for(i in 1:size){
470        isotope = .jcast(.jcall(isoIter,"Ljava/lang/Object;","next"), "org/openscience/cdk/interfaces/IIsotope");
471        isotopeList[i,1] <- .jcall(isotope,"S","getSymbol");
472        isotopeList[i,2] <- .jcall(molecularformula,"I","getIsotopeCount",isotope);
473        #massNum          <- .jcall(isotope,"Ljava/lang/Double;","getMassNumber");
474        massNum          <-  isotope$getMassNumber()
475        #exactMass        <- .jcall(isotope,"Ljava/lang/Double;","getExactMass");
476        exactMass        <- isotope$getExactMass()
477        #print(exactMass)
478        isotopes         <- J("org/openscience/cdk/config/Isotopes")
479        
480        if (is.null(massNum)) {
481         
482          isos         <- isotopes$getInstance()
483          majorIsotope <- isos$getMajorIsotope(isotope$getSymbol())
484          
485          if (!is.null(majorIsotope)) {
486            
487            exactMass <- majorIsotope$getExactMass()
488          }
489        } else {
490          
491          if (is.null(exactMass)) {
492            isos    <- isotopes$getInstance()
493            temp    <- isos$getIsotope(isotope$getSymbol(), massNum);
494            
495            if (!is.null(temp)) {
496              exactMass <- temp$getExactMass()
497            }
498          }
499        }
500        
501        #ch <- .jcall(isotope,"Ljava/lang/Double;","getExactMass");
502        isotopeList[i,3] <- exactMass
503        
504    }
505    
506    object@string <- .cdkFormula.getString(molecularformula);
507    manipulator <- get("mfManipulator", envir = .rcdk.GlobalEnv)
508    cMass <- .jcall(manipulator,"D","getTotalExactMass",molecularformula);
509    object@mass <- cMass;
510    chargeDO <- .jcall(molecularformula,"Ljava/lang/Integer;","getCharge");
511    charge <- .jcall(chargeDO,"D","doubleValue");
512    object <- set.charge.formula(object,charge)
513    object@isotopes <- isotopeList;
514    
515    return(object);
516}
517
518#' Construct an isotope pattern similarity calculator.
519#'
520#' A method that returns an instance of the CDK \code{IsotopePatternSimilarity}
521#' class which can be used to compute similarity scores between pairs of
522#' isotope abundance patterns.
523#'
524#' @param tol The tolerance
525#' @return A \code{jobjRef} corresponding to an instance of \code{IsotopePatternSimilarity}
526#' @seealso \code{\link{compare.isotope.pattern}}
527#' @references \url{http://cdk.github.io/cdk/1.5/docs/api/org/openscience/cdk/formula/IsotopePatternSimilarity.html}
528#' @author Miguel Rojas Cherto
529get.isotope.pattern.similarity <- function(tol = NULL) {
530    ips <- .jnew("org/openscience/cdk/formula/IsotopePatternSimilarity")
531    if (!is.null(tol)) ips$seTolerance(tol)
532    return(ips)
533}
534
535#' Construct an isotope pattern generator.
536#'
537#' Constructs an instance of the CDK \code{IsotopePatternGenerator}, with an optional
538#' minimum abundance specified. This object can be used to generate all combinatorial
539#' chemical isotopes given a structure.
540#'
541#' @param minAbundance The minimum abundance
542#' @return A \code{jobjRef} corresponding to an instance of \code{IsotopePatternGenerator}
543#' @references \url{http://cdk.github.io/cdk/1.5/docs/api/org/openscience/cdk/formula/IsotopePatternGenerator.html}
544#' @author Miguel Rojas Cherto
545get.isotope.pattern.generator <- function(minAbundance = NULL) {
546    if (is.null(minAbundance))
547        .jnew("org/openscience/cdk/formula/IsotopePatternGenerator")
548    else
549        .jnew("org/openscience/cdk/formula/IsotopePatternGenerator", as.double(minAbundance))
550}
551
552#' Compare isotope patterns.
553#'
554#' Computes a similarity score between two different isotope abundance patterns.
555#'
556#' @param iso1 The first isotope pattern, which should be a \code{jobjRef} corresponding to the \code{IsotopePattern} class
557#' @param iso2 The second isotope pattern, which should be a \code{jobjRef} corresponding to the \code{IsotopePattern} class
558#' @param ips An instance of the \code{IsotopePatternSimilarity} class. if \code{NULL} one will be constructed automatically
559#'
560#' @return A numeric value between 0 and 1 indicating the similarity between the two patterns
561#' @seealso \code{\link{get.isotope.pattern.similarity}}
562#' @export
563#' @references \url{http://cdk.github.io/cdk/2.3/docs/api/org/openscience/cdk/formula/IsotopePatternSimilarity.html}
564#' @author Miguel Rojas Cherto
565compare.isotope.pattern <- function(iso1, iso2, ips = NULL) {
566    cls <- unique(c(class(iso1), class(iso2)))
567    if (length(cls) != 1) stop("Must supply Java objects of class IsotopePattern")
568    if (cls != 'jobjRef') stop("Must supply Java objects of class IsotopePattern")
569    if(attr(iso1, "jclass") != "org/openscience/cdk/formula/IsotopePattern" ||
570       attr(iso2, "jclass") != "org/openscience/cdk/formula/IsotopePattern") {
571        stop("Must supply an IsotopePattern")
572    }
573    if (is.null(ips)) ips <- get.isotope.pattern.similarity()
574    return(ips$compare(iso1, iso2))
575}