2 VARNA is a tool for the automated drawing, visualization and annotation of the secondary structure of RNA, designed as a companion software for web servers and databases.
3 Copyright (C) 2008 Kevin Darty, Alain Denise and Yann Ponty.
4 electronic mail : Yann.Ponty@lri.fr
5 paper mail : LRI, bat 490 Universit� Paris-Sud 91405 Orsay Cedex France
7 This file is part of VARNA version 3.1.
8 VARNA version 3.1 is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License
9 as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
11 VARNA version 3.1 is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY;
12 without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
13 See the GNU General Public License for more details.
15 You should have received a copy of the GNU General Public License along with VARNA version 3.1.
16 If not, see http://www.gnu.org/licenses.
18 package fr.orsay.lri.varna.controlers;
20 import java.awt.event.ActionEvent;
21 import java.awt.event.ActionListener;
22 import java.awt.geom.Point2D;
23 import java.util.ArrayList;
24 import java.util.Arrays;
25 import java.util.LinkedList;
26 import java.util.Vector;
28 import javax.swing.Timer;
30 import fr.orsay.lri.varna.VARNAPanel;
31 import fr.orsay.lri.varna.exceptions.MappingException;
32 import fr.orsay.lri.varna.models.VARNAConfig;
33 import fr.orsay.lri.varna.models.rna.Mapping;
34 import fr.orsay.lri.varna.models.rna.ModeleBase;
35 import fr.orsay.lri.varna.models.rna.RNA;
37 public class ControleurInterpolator extends Thread implements ActionListener {
39 private static final int START = 0;
40 private static final int LOOP = 1;
41 private static final int END = 2;
44 protected int _numSteps = 25; // BH SwingJS
45 private long _timeDelay = /** @j2sNative 2 || */15;
46 private boolean _running = false;
47 Targets _d = new Targets();
48 protected int _step; // BH SwingJS
49 protected boolean _firstHalf; // BH SwingJS
51 public ControleurInterpolator(VARNAPanel vpn) {
55 public synchronized void addTarget(RNA target, Mapping mapping) {
56 addTarget(target, null, mapping);
60 public synchronized void addTarget(RNA target, VARNAConfig conf, Mapping mapping) {
61 _d.add(new TargetsHolder(target,conf,mapping));
64 public synchronized void addTarget(RNA target) {
65 addTarget(target, null, Mapping.DefaultOutermostMapping(_vpn.getRNA()
66 .get_listeBases().size(), target.get_listeBases().size()));
69 public boolean isInterpolationInProgress()
74 private Point2D.Double computeDestination(Point2D pli, Point2D pri,
75 Point2D pi, int i, int n, Point2D plf, Point2D prf) {
76 Point2D.Double plm = new Point2D.Double(
77 (pli.getX() + plf.getX()) / 2.0,
78 (pli.getY() + plf.getY()) / 2.0);
79 Point2D.Double prm = new Point2D.Double(
80 (pri.getX() + prf.getX()) / 2.0,
81 (pri.getY() + prf.getY()) / 2.0);
82 Point2D.Double pm = new Point2D.Double(((n - i) * plm.getX() + i
84 / n, ((n - i) * plm.getY() + i * prm.getY()) / n);
85 Point2D.Double v = new Point2D.Double(pm.getX() - pi.getX(), pm.getY()
87 Point2D.Double pf = new Point2D.Double(pi.getX() + 2.0 * v.getX(), pi
93 private Vector<Vector<Integer>> clusterIndices(int numIndices,
94 int[] mappedIndices) throws MappingException {
95 int[] indices = new int[numIndices];
96 for (int i = 0; i < numIndices; i++) {
99 return clusterIndices(indices, mappedIndices);
103 * Builds and returns an interval array, alternating unmatched regions and
104 * matched indices. Namely, returns a vector of vector such that the vectors
105 * found at odd positions contain those indices that are NOT associated with
106 * any other base in the current mapping. On the other hand, vectors found
107 * at even positions contain only one mapped index. <br/>
108 * Ex: If indices=<code>[1,2,3,4,5]</code> and mappedIndices=
109 * <code>[2,5]</code> then the function will return
110 * <code>[[1],[2],[3,4],[5],[]]</code>.
113 * The total list of indices
114 * @param mappedIndices
115 * Matched indices, should be a subset of <code>indices</code>
116 * @return A clustered array
117 * @throws MappingException
118 * If one of the parameters is an empty array
120 private Vector<Vector<Integer>> clusterIndices(int[] indices,
121 int[] mappedIndices) throws MappingException {
122 if ((mappedIndices.length == 0) || (indices.length == 0)) {
123 throw new MappingException(
124 "Mapping Error: Cannot cluster indices in an empty mapping");
126 Vector<Vector<Integer>> res = new Vector<Vector<Integer>>();
128 Arrays.sort(indices);
129 Arrays.sort(mappedIndices);
131 Vector<Integer> tmp = new Vector<Integer>();
132 for (i = 0; (i < indices.length) && (j < mappedIndices.length); i++) {
133 if (indices[i] == mappedIndices[j]) {
135 tmp = new Vector<Integer>();
138 tmp = new Vector<Integer>();
145 for (i = k; (i < indices.length); i++) {
157 TargetsHolder d = _d.get();
160 nextTarget(d.target,d.conf,d.mapping);
164 System.err.println(e);
178 * Compute the centroid of the RNA bases that have their indexes
179 * in the given array.
181 private static Point2D.Double computeCentroid(ArrayList<ModeleBase> rnaBases, int[] indexes) {
182 double centroidX = 0, centroidY = 0;
183 for (int i=0; i<indexes.length; i++) {
184 int index = indexes[i];
185 Point2D coords = rnaBases.get(index).getCoords();
186 centroidX += coords.getX();
187 centroidY += coords.getY();
189 centroidX /= indexes.length;
190 centroidY /= indexes.length;
192 return new Point2D.Double(centroidX, centroidY);
196 * The purpose of this class is to compute the rotation that minimizes the
197 * RMSD between the bases of the first RNA and the bases of the second RNA.
199 * @author Raphael Champeimont
201 private static class MinimizeRMSD {
202 private double[] X1, X2, Y1, Y2;
204 public MinimizeRMSD(double[] X1, double[] Y1, double[] X2, double[] Y2) {
212 * A function such that minimizing it is equivalent to
213 * minimize the RMSD between between the two arrays of vectors,
214 * supposing we rotate the points in [X2,Y2] by theta.
216 private double f(double theta) {
217 double cos_theta = Math.cos(theta);
218 double sin_theta = Math.sin(theta);
221 for (int i=0; i<n; i++) {
222 double dsx = X2[i]*cos_theta - Y2[i]*sin_theta - X1[i];
223 double dsy = X2[i]*sin_theta + Y2[i]*cos_theta - Y1[i];
224 sum = sum + dsx*dsx + dsy*dsy;
232 private double fprime(double theta) {
233 double cos_theta = Math.cos(theta);
234 double sin_theta = Math.sin(theta);
237 for (int i=0; i<n; i++) {
239 + (X1[i]*X2[i] + Y1[i]*Y2[i]) * sin_theta
240 + (X1[i]*Y2[i] - X2[i]*Y1[i]) * cos_theta;
246 * d^2/dtheta^2 f(theta)
248 private double fsecond(double theta) {
249 double cos_theta = Math.cos(theta);
250 double sin_theta = Math.sin(theta);
253 for (int i=0; i<n; i++) {
255 + (X1[i]*X2[i] + Y1[i]*Y2[i]) * cos_theta
256 + (X2[i]*Y1[i] - X1[i]*Y2[i]) * sin_theta;
263 * Find the theta that minimizes f(theta).
264 * We use Newton's method.
266 public double computeOptimalTheta() {
267 final double epsilon = 1E-5;
270 // In practice the number of steps to reach precision 1E-5 is smaller that
271 // 10 most of the time.
272 final int maxnumsteps = 100;
276 numsteps = numsteps + 1;
277 double d = fsecond(x_n);
279 // if f''(x_n) is 0 we cannot divide by it,
280 // so we move a little.
281 x_n_plus_1 = x_n + epsilon * (Math.random() - 0.5);
283 x_n_plus_1 = x_n - fprime(x_n)/fsecond(x_n);
284 if (Math.abs(x_n_plus_1 - x_n) < epsilon) {
289 if (numsteps >= maxnumsteps) {
290 // If things go bad (for example f''(x) keeps being 0)
291 // we need to give up after what we consider to be too many steps.
292 // In practice this can happen only in pathological cases
293 // like f being constant, which is very unlikely.
300 // We now have either found the min or the max at x = result.
301 // If we have the max at x we know the min is at x+pi.
302 if (f(result + Math.PI) < f(result)) {
303 result = result + Math.PI;
314 * We suppose we have two lists of points. The coordinates of the first
315 * list of points are X1 and Y1, and X2 and Y2 for the second list.
316 * We suppose that the centroids of [X1,Y1] and [X2,Y2] are both (0,0).
317 * This function computes the rotation to apply to the second set of
318 * points such that the RMSD (Root mean square deviation) between them
319 * is minimum. The returned angle is in radians.
321 private static double minimizeRotateRMSD(double[] X1, double[] Y1, double[] X2, double[] Y2) {
322 MinimizeRMSD minimizer = new MinimizeRMSD(X1, Y1, X2, Y2);
323 return minimizer.computeOptimalTheta();
328 * Move rna2 using a rotation so that it would move as little as possible
329 * (in the sense that we minimize the RMSD) when transforming
330 * it to rna1 using the given mapping.
332 public static void moveNearOtherRNA(RNA rna1, RNA rna2, Mapping mapping) {
333 int[] rna1MappedElems = mapping.getSourceElems();
334 int[] rna2MappedElems = mapping.getTargetElems();
335 ArrayList<ModeleBase> rna1Bases = rna1.get_listeBases();
336 ArrayList<ModeleBase> rna2Bases = rna2.get_listeBases();
337 int n = rna1MappedElems.length;
339 // If there is less than 2 points, it is useless to rotate the RNA.
342 // We can now assume that n >= 2.
344 // Compute the centroids of both RNAs
345 Point2D.Double rna1MappedElemsCentroid = computeCentroid(rna1Bases, rna1MappedElems);
346 Point2D.Double rna2MappedElemsCentroid = computeCentroid(rna2Bases, rna2MappedElems);
348 // Compute the optimal rotation
349 // We first compute coordinates for both RNAs, changing the origins
350 // to be the centroids.
351 double[] X1 = new double[rna1MappedElems.length];
352 double[] Y1 = new double[rna1MappedElems.length];
353 double[] X2 = new double[rna2MappedElems.length];
354 double[] Y2 = new double[rna2MappedElems.length];
355 for (int i=0; i<rna1MappedElems.length; i++) {
356 int base1Index = rna1MappedElems[i];
357 Point2D.Double coords1 = rna1Bases.get(base1Index).getCoords();
358 X1[i] = coords1.x - rna1MappedElemsCentroid.x;
359 Y1[i] = coords1.y - rna1MappedElemsCentroid.y;
360 Point2D.Double coords2 = rna2Bases.get(mapping.getPartner(base1Index)).getCoords();
361 X2[i] = coords2.x - rna2MappedElemsCentroid.x;
362 Y2[i] = coords2.y - rna2MappedElemsCentroid.y;
365 // Compute the optimal rotation angle theta
366 double theta = minimizeRotateRMSD(X1, Y1, X2, Y2);
369 rna2.globalRotation(theta * 180.0 / Math.PI);
372 public void nextTarget(RNA _target, VARNAConfig _conf, Mapping _mapping) {
373 nextTarget(_target, _conf, _mapping, false);
376 Runnable _loop, _end;
379 * The argument moveTarget specifies whether the RNA _target should
380 * be rotated so that bases move as little as possible when switching
381 * from the current RNA to _target using the animation.
382 * Note that this will modify the _target object directly.
384 public void nextTarget(final RNA _target, final VARNAConfig _conf, Mapping _mapping, boolean moveTarget)
386 _end = new Runnable() {
390 _vpn.showRNA(_target);
397 final RNA source = _vpn.getRNA();
399 if (moveTarget) moveNearOtherRNA(source, _target, _mapping);
401 if (source.getSize()!=0&&_target.getSize()!=0)
403 ArrayList<ModeleBase> currBases = source.get_listeBases();
404 ArrayList<ModeleBase> destBases = _target.get_listeBases();
405 Vector<Vector<Integer>> intArrSource = new Vector<Vector<Integer>>();
406 Vector<Vector<Integer>> intArrTarget = new Vector<Vector<Integer>>();
407 // Building interval arrays
408 intArrSource = clusterIndices(currBases.size(), _mapping
410 intArrTarget = clusterIndices(destBases.size(), _mapping
413 // Duplicating source and target coordinates
414 final Point2D.Double[] initPosSource = new Point2D.Double[currBases
416 final Point2D.Double[] finalPosTarget = new Point2D.Double[destBases
419 for (int i = 0; i < currBases.size(); i++) {
420 Point2D tmp = currBases.get(i).getCoords();
421 initPosSource[i] = new Point2D.Double(tmp.getX(), tmp.getY());
423 for (int i = 0; i < destBases.size(); i++) {
424 Point2D tmp = destBases.get(i).getCoords();
425 finalPosTarget[i] = new Point2D.Double(tmp.getX(), tmp.getY());
429 * Assigning final (Source) and initial (Target) coordinates
431 final Point2D.Double[] finalPosSource = new Point2D.Double[initPosSource.length];
432 final Point2D.Double[] initPosTarget = new Point2D.Double[finalPosTarget.length];
433 // Final position of source model
434 for (int i = 0; i < finalPosSource.length; i++) {
435 if (_mapping.getPartner(i) != Mapping.UNKNOWN) {
437 dest = finalPosTarget[_mapping.getPartner(i)];
438 finalPosSource[i] = new Point2D.Double(dest.getX(), dest
443 for (int i = 0; i < intArrSource.size(); i += 2) {
444 int matchedNeighborLeft, matchedNeighborRight;
446 matchedNeighborLeft = intArrSource.get(1).get(0);
447 matchedNeighborRight = intArrSource.get(1).get(0);
448 } else if (i == intArrSource.size() - 1) {
449 matchedNeighborLeft = intArrSource.get(
450 intArrSource.size() - 2).get(0);
451 matchedNeighborRight = intArrSource.get(
452 intArrSource.size() - 2).get(0);
454 matchedNeighborLeft = intArrSource.get(i - 1).get(0);
455 matchedNeighborRight = intArrSource.get(i + 1).get(0);
457 Vector<Integer> v = intArrSource.get(i);
458 for (int j = 0; j < v.size(); j++) {
459 int index = v.get(j);
460 finalPosSource[index] = computeDestination(
461 initPosSource[matchedNeighborLeft],
462 initPosSource[matchedNeighborRight],
463 initPosSource[index], j + 1, v.size() + 1,
464 finalPosSource[matchedNeighborLeft],
465 finalPosSource[matchedNeighborRight]);
468 for (int i = 0; i < initPosTarget.length; i++) {
469 if (_mapping.getAncestor(i) != Mapping.UNKNOWN) {
471 dest = initPosSource[_mapping.getAncestor(i)];
472 initPosTarget[i] = new Point2D.Double(dest.getX(), dest.getY());
475 for (int i = 0; i < intArrTarget.size(); i += 2) {
476 int matchedNeighborLeft, matchedNeighborRight;
478 matchedNeighborLeft = intArrTarget.get(1).get(0);
479 matchedNeighborRight = intArrTarget.get(1).get(0);
480 } else if (i == intArrTarget.size() - 1) {
481 matchedNeighborLeft = intArrTarget.get(
482 intArrTarget.size() - 2).get(0);
483 matchedNeighborRight = intArrTarget.get(
484 intArrTarget.size() - 2).get(0);
486 matchedNeighborLeft = intArrTarget.get(i - 1).get(0);
487 matchedNeighborRight = intArrTarget.get(i + 1).get(0);
489 Vector<Integer> v = intArrTarget.get(i);
490 for (int j = 0; j < v.size(); j++) {
491 int index = v.get(j);
492 initPosTarget[index] = computeDestination(
493 finalPosTarget[matchedNeighborLeft],
494 finalPosTarget[matchedNeighborRight],
495 finalPosTarget[index], j + 1, v.size() + 1,
496 initPosTarget[matchedNeighborLeft],
497 initPosTarget[matchedNeighborRight]);
502 _loop = new Runnable(){
507 RNA current = (_firstHalf ? source : _target);
508 if (i == _numSteps / 2) {
509 _vpn.showRNA(_target);
513 {_vpn.setConfig(_conf);}
515 for (int j = 0; j < initPosSource.length; j++) {
516 source.setCoord(j, initPosSource[j]);
519 ArrayList<ModeleBase> currBases = current.get_listeBases();
520 for (int j = 0; j < currBases.size(); j++) {
521 ModeleBase m = currBases.get(j);
524 mpc = initPosSource[j];
525 mnc = finalPosSource[j];
527 mpc = initPosTarget[j];
528 mnc = finalPosTarget[j];
530 m.setCoords(new Point2D.Double(((_numSteps - 1 - i)
531 * mpc.getX() + (i) * mnc.getX())
532 / (_numSteps - 1), ((_numSteps - 1 - i)
533 * mpc.getY() + i * mnc.getY())
540 actionPerformed(null);
544 } catch (MappingException e) {
547 }catch (Exception e) {
555 private class TargetsHolder
558 public VARNAConfig conf;
559 public Mapping mapping;
560 public TargetsHolder(RNA t, VARNAConfig c, Mapping m)
568 private class Targets
570 LinkedList<TargetsHolder> _d = new LinkedList<TargetsHolder>();
572 // BH j2s SwingJS added only to remove Eclipse warning
574 public synchronized void add(TargetsHolder d)
578 @SuppressWarnings("unused")
579 // BH SwingJS no notify()
580 Runnable interpolator = ControleurInterpolator.this;
581 if (/** @j2sNative 1?true:*/false)
587 public synchronized TargetsHolder get() {
590 * BH SwingJS no wait()
597 while (_d.size() == 0) {
600 } catch (InterruptedException e) {
605 TargetsHolder x = _d.getLast();
613 public void actionPerformed(ActionEvent e) {
619 private void runAnimation() {
627 if (_step < _numSteps) {
638 Timer t = new Timer((int) _timeDelay, this);
642 // for (int i = 0; i < _numSteps; i++) {
646 // sleep(_timeDelay);
649 // } catch (InterruptedException e) {
650 // e.printStackTrace();