srcjar/fr/orsay/lri/varna/models/naView/NAView.java

   1 /*
   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
   6
   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.
  10
  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.
  14
  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.
  17  */
  18 package fr.orsay.lri.varna.models.naView;
  19
  20
  21 import java.util.ArrayList;
  22
  23 import fr.orsay.lri.varna.exceptions.ExceptionNAViewAlgorithm;
  24 import fr.orsay.lri.varna.interfaces.InterfaceVARNAListener;
  25 import fr.orsay.lri.varna.interfaces.InterfaceVARNAObservable;
  26
  27 public class NAView {
  28         private final double ANUM = 9999.0;
  29         private final int MAXITER = 500;
  30
  31         private ArrayList<Base> bases;
  32         private int nbase, nregion, loop_count;
  33
  34         private Loop root = new Loop();
  35         private ArrayList<Loop> loops;
  36
  37         private ArrayList<Region> regions;
  38
  39         private Radloop rlphead = new Radloop();
  40
  41         private double lencut=0.8;
  42         private final double RADIUS_REDUCTION_FACTOR = 1.4;
  43
  44         // show algorithm step by step
  45         private boolean debug = false;
  46
  47         private double angleinc;
  48
  49         private double _h;
  50
  51         private ArrayList<InterfaceVARNAListener> _listeVARNAListener = new ArrayList<InterfaceVARNAListener>();
  52
  53         private boolean noIterationFailureYet = true;
  54
  55         double HELIX_FACTOR = 0.6;
  56         double BACKBONE_DISTANCE = 27;
  57
  58         public int naview_xy_coordinates(ArrayList<Short> pair_table2,
  59                         ArrayList<Double> x, ArrayList<Double> y)
  60                         throws ExceptionNAViewAlgorithm {
  61                 if (debug)
  62                         System.out.println("naview_xy_coordinates");
  63                 if (pair_table2.size() == 0)
  64                         return 0;
  65                 int i;
  66                 ArrayList<Integer> pair_table = new ArrayList<Integer>(pair_table2
  67                                 .size() + 1);
  68                 pair_table.add(pair_table2.size());
  69
  70                 for (int j = 0; j < pair_table2.size(); j++) {
  71                         pair_table.add(pair_table2.get(j) + 1);
  72                 }
  73
  74                 if (debug) {
  75                         infoStructure(pair_table);
  76                 }
  77                 // length
  78                 nbase = pair_table.get(0);
  79                 bases = new ArrayList<Base>(nbase + 1);
  80
  81                 for (int index = 0; index < bases.size(); index++) {
  82                         bases.add(new Base());
  83                 }
  84
  85                 regions = new ArrayList<Region>();
  86                 for (int index = 0; index < nbase + 1; index++) {
  87                         regions.add(new Region());
  88                 }
  89
  90                 read_in_bases(pair_table);
  91
  92                 if (debug)
  93                         infoBasesMate();
  94
  95                 rlphead = null;
  96
  97                 find_regions();
  98
  99                 if (debug)
 100                         infoRegions();
 101
 102                 loop_count = 0;
 103                 loops = new ArrayList<Loop>(nbase + 1);
 104                 for (int index = 0; index < nbase + 1; index++) {
 105                         loops.add(new Loop());
 106                 }
 107
 108                 construct_loop(0);
 109
 110                 if (debug)
 111                         infoBasesExtracted();
 112
 113                 find_central_loop();
 114
 115                 if (debug)
 116                         infoRoot();
 117
 118                 if (debug)
 119                         dump_loops();
 120
 121                 traverse_loop(root, null);
 122
 123                 for (i = 0; i < nbase; i++) {
 124                         x.add(100 + BACKBONE_DISTANCE * bases.get(i + 1).getX());
 125                         y.add(100 + BACKBONE_DISTANCE * bases.get(i + 1).getY());
 126                 }
 127
 128                 return nbase;
 129         }
 130
 131         private void infoStructure(ArrayList<Integer> pair_table) {
 132                 System.out.println("structure:");
 133                 for (int j = 0; j < pair_table.size(); j++) {
 134                         System.out.print("#" + j + ":" + pair_table.get(j) + "\t");
 135                         if (j % 10 == 0)
 136                                 System.out.println();
 137                 }
 138                 System.out.println();
 139         }
 140
 141         private void infoBasesMate() {
 142                 System.out.println("Bases mate:");
 143                 for (int index = 0; index < bases.size(); index++) {
 144                         System.out.print("#" + index + ":" + bases.get(index).getMate()
 145                                         + "\t");
 146                         if (index % 10 == 0)
 147                                 System.out.println();
 148                 }
 149                 System.out.println();
 150         }
 151
 152         private void infoRegions() {
 153                 System.out.println("regions:");
 154                 for (int index = 0; index < regions.size(); index++) {
 155                         System.out.print("(" + regions.get(index).getStart1() + ","
 156                                         + regions.get(index).getStart2() + ";"
 157                                         + regions.get(index).getEnd1() + ","
 158                                         + regions.get(index).getEnd2() + ")\t\t");
 159                         if (index % 5 == 0)
 160                                 System.out.println();
 161                 }
 162                 System.out.println();
 163         }
 164
 165         private void infoBasesExtracted() {
 166                 System.out.println("Bases extracted:");
 167                 for (int index = 0; index < bases.size(); index++) {
 168                         System.out.print("i=" + index + ":"
 169                                         + bases.get(index).isExtracted() + "\t");
 170                         if (index % 5 == 0)
 171                                 System.out.println();
 172                 }
 173                 System.out.println();
 174         }
 175
 176         private void infoRoot() {
 177                 System.out.println("root" + root.getNconnection() + ";"
 178                                 + root.getNumber());
 179                 System.out.println("\troot : ");
 180                 System.out.println("\tdepth=" + root.getDepth());
 181                 System.out.println("\tmark=" + root.isMark());
 182                 System.out.println("\tnumber=" + root.getNumber());
 183                 System.out.println("\tradius=" + root.getRadius());
 184                 System.out.println("\tx=" + root.getX());
 185                 System.out.println("\ty=" + root.getY());
 186                 System.out.println("\tnconnection=" + root.getNconnection());
 187         }
 188
 189         private void read_in_bases(ArrayList<Integer> pair_table) {
 190                 if (debug)
 191                         System.out.println("read_in_bases");
 192
 193                 int i, npairs;
 194
 195                 // Set up an origin.
 196                 bases.add(new Base());
 197                 bases.get(0).setMate(0);
 198                 bases.get(0).setExtracted(false);
 199                 bases.get(0).setX(ANUM);
 200                 bases.get(0).setY(ANUM);
 201
 202                 for (npairs = 0, i = 1; i <= nbase; i++) {
 203                         bases.add(new Base());
 204                         bases.get(i).setExtracted(false);
 205                         bases.get(i).setX(ANUM);
 206                         bases.get(i).setY(ANUM);
 207                         bases.get(i).setMate(pair_table.get(i));
 208                         if ((int) pair_table.get(i) > i)
 209                                 npairs++;
 210                 }
 211                 // must have at least 1 pair to avoid segfault
 212                 if (npairs == 0) {
 213                         bases.get(1).setMate(nbase);
 214                         bases.get(nbase).setMate(1);
 215                 }
 216         }
 217
 218         /**
 219          * Identifies the regions in the structure.
 220          */
 221         private void find_regions()
 222
 223         {
 224                 if (debug)
 225                         System.out.println("find_regions");
 226                 int i, mate, nb1;
 227                 nb1 = nbase + 1;
 228                 ArrayList<Boolean> mark = new ArrayList<Boolean>(nb1);
 229                 for (i = 0; i < nb1; i++)
 230                         mark.add(false);
 231                 nregion = 0;
 232                 for (i = 0; i <= nbase; i++) {
 233                         if ((mate = bases.get(i).getMate()) != 0 && !mark.get(i)) {
 234                                 regions.get(nregion).setStart1(i);
 235                                 regions.get(nregion).setEnd2(mate);
 236                                 mark.set(i, true);
 237                                 mark.set(mate, true);
 238                                 bases.get(i).setRegion(regions.get(nregion));
 239                                 bases.get(mate).setRegion(regions.get(nregion));
 240                                 for (i++, mate--; i < mate && bases.get(i).getMate() == mate; i++, mate--) {
 241                                         mark.set(mate, true);
 242                                         mark.set(i, true);
 243                                         bases.get(i).setRegion(regions.get(nregion));
 244                                         bases.get(mate).setRegion(regions.get(nregion));
 245                                 }
 246                                 regions.get(nregion).setEnd1(--i);
 247                                 regions.get(nregion).setStart2(mate + 1);
 248                                 if (debug) {
 249                                         if (nregion == 0)
 250                                                 System.out.printf("\nRegions are:\n");
 251                                         System.out.printf(
 252                                                         "Region %d is %d-%d and %d-%d with gap of %d.\n",
 253                                                         nregion + 1, regions.get(nregion).getStart1(),
 254                                                         regions.get(nregion).getEnd1(), regions
 255                                                                         .get(nregion).getStart2(), regions.get(
 256                                                                         nregion).getEnd2(), regions.get(nregion)
 257                                                                         .getStart2()
 258                                                                         - regions.get(nregion).getEnd1() + 1);
 259                                 }
 260                                 nregion++;
 261                         }
 262                 }
 263         }
 264
 265         /**
 266          * Starting at residue ibase, recursively constructs the loop containing
 267          * said base and all deeper bases.
 268          *
 269          * @throws ExceptionNAViewAlgorithm
 270          */
 271         private Loop construct_loop(int ibase) throws ExceptionNAViewAlgorithm {
 272                 if (debug)
 273                         System.out.println("construct_loop");
 274                 int i, mate;
 275                 Loop retloop = new Loop(), lp = new Loop();
 276                 Connection cp = new Connection();
 277                 Region rp = new Region();
 278                 Radloop rlp = new Radloop();
 279                 retloop = loops.get(loop_count++);
 280                 retloop.setNconnection(0);
 281                 //System.out.println(""+ibase+" "+nbase);
 282                 //ArrayList<Connection> a = new ArrayList<Connection>(nbase + 1);
 283                 //retloop.setConnections(a);
 284 //              for (int index = 0; index < nbase + 1; index++)
 285 //                      retloop.getConnections().add(new Connection());
 286                 //for (int index = 0; index < nbase + 1; index++)
 287                 //      retloop.addConnection(index,new Connection());
 288                 retloop.setDepth(0);
 289                 retloop.setNumber(loop_count);
 290                 retloop.setRadius(0.0);
 291                 for (rlp = rlphead; rlp != null; rlp = rlp.getNext())
 292                         if (rlp.getLoopnumber() == loop_count)
 293                                 retloop.setRadius(rlp.getRadius());
 294                 i = ibase;
 295                 do {
 296                         if ((mate = bases.get(i).getMate()) != 0) {
 297                                 rp = bases.get(i).getRegion();
 298                                 if (!bases.get(rp.getStart1()).isExtracted()) {
 299                                         if (i == rp.getStart1()) {
 300                                                 bases.get(rp.getStart1()).setExtracted(true);
 301                                                 bases.get(rp.getEnd1()).setExtracted(true);
 302                                                 bases.get(rp.getStart2()).setExtracted(true);
 303                                                 bases.get(rp.getEnd2()).setExtracted(true);
 304                                                 lp = construct_loop(rp.getEnd1() < nbase ? rp.getEnd1() + 1
 305                                                                 : 0);
 306                                         } else if (i == rp.getStart2()) {
 307                                                 bases.get(rp.getStart2()).setExtracted(true);
 308                                                 bases.get(rp.getEnd2()).setExtracted(true);
 309                                                 bases.get(rp.getStart1()).setExtracted(true);
 310                                                 bases.get(rp.getEnd1()).setExtracted(true);
 311                                                 lp = construct_loop(rp.getEnd2() < nbase ? rp.getEnd2() + 1
 312                                                                 : 0);
 313                                         } else {
 314                                                 throw new ExceptionNAViewAlgorithm(
 315                                                                 "naview:Error detected in construct_loop. i = "
 316                                                                                 + i + " not found in region table.\n");
 317                                         }
 318                                         retloop.setNconnection(retloop.getNconnection() + 1);
 319                                         cp = new Connection();
 320                                         retloop.setConnection(retloop.getNconnection() - 1,     cp);
 321                                         retloop.setConnection(retloop.getNconnection(), null);
 322                                         cp.setLoop(lp);
 323                                         cp.setRegion(rp);
 324                                         if (i == rp.getStart1()) {
 325                                                 cp.setStart(rp.getStart1());
 326                                                 cp.setEnd(rp.getEnd2());
 327                                         } else {
 328                                                 cp.setStart(rp.getStart2());
 329                                                 cp.setEnd(rp.getEnd1());
 330                                         }
 331                                         cp.setExtruded(false);
 332                                         cp.setBroken(false);
 333                                         lp.setNconnection(lp.getNconnection() + 1);
 334                                         cp = new Connection();
 335                                         lp.setConnection(lp.getNconnection() - 1, cp);
 336                                         lp.setConnection(lp.getNconnection(), null);
 337                                         cp.setLoop(retloop);
 338                                         cp.setRegion(rp);
 339                                         if (i == rp.getStart1()) {
 340                                                 cp.setStart(rp.getStart2());
 341                                                 cp.setEnd(rp.getEnd1());
 342                                         } else {
 343                                                 cp.setStart(rp.getStart1());
 344                                                 cp.setEnd(rp.getEnd2());
 345                                         }
 346                                         cp.setExtruded(false);
 347                                         cp.setBroken(false);
 348                                 }
 349                                 i = mate;
 350                         }
 351                         if (++i > nbase)
 352                                 i = 0;
 353                 } while (i != ibase);
 354                 return retloop;
 355         }
 356
 357         /**
 358          * Displays all the loops.
 359          */
 360         private void dump_loops() {
 361                 System.out.println("dump_loops");
 362                 int il, ilp, irp;
 363                 Loop lp;
 364                 Connection cp;
 365
 366                 System.out.printf("\nRoot loop is #%d\n", loops.indexOf(root) + 1);
 367                 for (il = 0; il < loop_count; il++) {
 368                         lp = loops.get(il);
 369                         System.out.printf("Loop %d has %d connections:\n", il + 1, lp
 370                                         .getNconnection());
 371                         for (int i = 0; (cp = lp.getConnection(i)) != null; i++) {
 372                                 ilp = (loops.indexOf(cp.getLoop())) + 1;
 373                                 irp = (regions.indexOf(cp.getRegion())) + 1;
 374                                 System.out.printf("  Loop %d Region %d (%d-%d)\n", ilp, irp, cp
 375                                                 .getStart(), cp.getEnd());
 376                         }
 377                 }
 378         }
 379
 380         /**
 381          * Find node of greatest branching that is deepest.
 382          */
 383         private void find_central_loop() {
 384                 if (debug)
 385                         System.out.println("find_central_loop");
 386                 Loop lp = new Loop();
 387                 int maxconn, maxdepth, i;
 388
 389                 determine_depths();
 390                 maxconn = 0;
 391                 maxdepth = -1;
 392                 for (i = 0; i < loop_count; i++) {
 393                         lp = loops.get(i);
 394                         if (lp.getNconnection() > maxconn) {
 395                                 maxdepth = lp.getDepth();
 396                                 maxconn = lp.getNconnection();
 397                                 root = lp;
 398                         } else if (lp.getDepth() > maxdepth
 399                                         && lp.getNconnection() == maxconn) {
 400                                 maxdepth = lp.getDepth();
 401                                 root = lp;
 402                         }
 403                 }
 404         }
 405
 406         /**
 407          * Determine the depth of all loops.
 408          */
 409         private void determine_depths() {
 410                 if (debug)
 411                         System.out.println("determine_depths");
 412                 Loop lp = new Loop();
 413                 int i, j;
 414
 415                 for (i = 0; i < loop_count; i++) {
 416                         lp = loops.get(i);
 417                         for (j = 0; j < loop_count; j++)
 418                                 loops.get(j).setMark(false);
 419                         lp.setDepth(depth(lp));
 420                 }
 421         }
 422
 423         /**
 424          * Determines the depth of loop, lp. Depth is defined as the minimum
 425          * distance to a leaf loop where a leaf loop is one that has only one or no
 426          * connections.
 427          */
 428         private int depth(Loop lp) {
 429                 if (debug)
 430                         System.out.println("depth");
 431                 int count, ret, d;
 432
 433                 if (lp.getNconnection() <= 1)
 434                         return 0;
 435                 if (lp.isMark())
 436                         return -1;
 437                 lp.setMark(true);
 438                 count = 0;
 439                 ret = 0;
 440                 for (int i = 0; lp.getConnection(i) != null; i++) {
 441                         d = depth(lp.getConnection(i).getLoop());
 442                         if (d >= 0) {
 443                                 if (++count == 1)
 444                                         ret = d;
 445                                 else if (ret > d)
 446                                         ret = d;
 447                         }
 448                 }
 449                 lp.setMark(false);
 450                 return ret + 1;
 451         }
 452
 453         /**
 454          * This is the workhorse of the display program. The algorithm is recursive
 455          * based on processing individual loops. Each base pairing region is
 456          * displayed using the direction given by the circle diagram, and the
 457          * connections between the regions is drawn by equally spaced points. The
 458          * radius of the loop is set to minimize the square error for lengths
 459          * between sequential bases in the loops. The "correct" length for base
 460          * links is 1. If the least squares fitting of the radius results in loops
 461          * being less than 1/2 unit apart, then that segment is extruded.
 462          *
 463          * The variable, anchor_connection, gives the connection to the loop
 464          * processed in an previous level of recursion.
 465          *
 466          * @throws ExceptionNAViewAlgorithm
 467          */
 468         private void traverse_loop(Loop lp, Connection anchor_connection)
 469                         throws ExceptionNAViewAlgorithm {
 470                 if (debug)
 471                         System.out.println("  traverse_loop");
 472                 double xs, ys, xe, ye, xn, yn, angleinc, r;
 473                 double radius, xc, yc, xo, yo, astart, aend, a;
 474                 Connection cp, cpnext, acp, cpprev;
 475                 int i, j, n, ic;
 476                 double da, maxang;
 477                 int count, icstart, icend, icmiddle, icroot;
 478                 boolean done, done_all_connections, rooted;
 479                 int sign;
 480                 double midx, midy, nrx, nry, mx, my, vx, vy, dotmv, nmidx, nmidy;
 481                 int icstart1, icup, icdown, icnext, direction;
 482                 double dan, dx, dy, rr;
 483                 double cpx, cpy, cpnextx, cpnexty, cnx, cny, rcn, rc, lnx, lny, rl, ac, acn, sx, sy, dcp;
 484                 int imaxloop = 0;
 485
 486                 angleinc = 2 * Math.PI / (nbase + 1);
 487                 acp = null;
 488                 icroot = -1;
 489                 int indice = 0;
 490
 491                 for (ic = 0; (cp = lp.getConnection(indice)) != null; indice++, ic++) {
 492                         // xs = cos(angleinc*cp.setStart(); ys = sin(angleinc*cp.setStart();
 493                         // xe =
 494                         // cos(angleinc*cp.setEnd()); ye = sin(angleinc*cp.setEnd());
 495                         xs = -Math.sin(angleinc * cp.getStart());
 496                         ys = Math.cos(angleinc * cp.getStart());
 497                         xe = -Math.sin(angleinc * cp.getEnd());
 498                         ye = Math.cos(angleinc * cp.getEnd());
 499                         xn = ye - ys;
 500                         yn = xs - xe;
 501                         r = Math.sqrt(xn * xn + yn * yn);
 502                         cp.setXrad(xn / r);
 503                         cp.setYrad(yn / r);
 504                         cp.setAngle(Math.atan2(yn, xn));
 505                         if (cp.getAngle() < 0.0)
 506                                 cp.setAngle(cp.getAngle() + 2 * Math.PI);
 507                         if (anchor_connection != null
 508                                         && anchor_connection.getRegion() == cp.getRegion()) {
 509                                 acp = cp;
 510                                 icroot = ic;
 511                         }
 512                 }
 513                 // remplacement d'une etiquette de goto
 514                 set_radius: while (true) {
 515                         determine_radius(lp, lencut);
 516                         radius = lp.getRadius()/RADIUS_REDUCTION_FACTOR;
 517                         if (anchor_connection == null)
 518                                 xc = yc = 0.0;
 519                         else {
 520                                 xo = (bases.get(acp.getStart()).getX() + bases
 521                                                 .get(acp.getEnd()).getX()) / 2.0;
 522                                 yo = (bases.get(acp.getStart()).getY() + bases
 523                                                 .get(acp.getEnd()).getY()) / 2.0;
 524                                 xc = xo - radius * acp.getXrad();
 525                                 yc = yo - radius * acp.getYrad();
 526                         }
 527
 528                         // The construction of the connectors will proceed in blocks of
 529                         // connected connectors, where a connected connector pairs means two
 530                         // connectors that are forced out of the drawn circle because they
 531                         // are too close together in angle.
 532
 533                         // First, find the start of a block of connected connectors
 534
 535                         if (icroot == -1)
 536                                 icstart = 0;
 537                         else
 538                                 icstart = icroot;
 539                         cp = lp.getConnection(icstart);
 540                         count = 0;
 541                         if (debug)
 542                         {
 543                                 System.out.printf("Now processing loop %d\n", lp.getNumber());
 544                                 System.out.println("  "+lp);
 545                         }
 546                         done = false;
 547                         do {
 548                                 j = icstart - 1;
 549                                 if (j < 0)
 550                                         j = lp.getNconnection() - 1;
 551                                 cpprev = lp.getConnection(j);
 552                                 if (!connected_connection(cpprev, cp)) {
 553                                         done = true;
 554                                 } else {
 555                                         icstart = j;
 556                                         cp = cpprev;
 557                                 }
 558                                 if (++count > lp.getNconnection()) {
 559                                         // Here everything is connected. Break on maximum angular
 560                                         // separation between connections.
 561
 562                                         maxang = -1.0;
 563                                         for (ic = 0; ic < lp.getNconnection(); ic++) {
 564                                                 j = ic + 1;
 565                                                 if (j >= lp.getNconnection())
 566                                                         j = 0;
 567                                                 cp = lp.getConnection(ic);
 568                                                 cpnext = lp.getConnection(j);
 569                                                 ac = cpnext.getAngle() - cp.getAngle();
 570                                                 if (ac < 0.0)
 571                                                         ac += 2 * Math.PI;
 572                                                 if (ac > maxang) {
 573                                                         maxang = ac;
 574                                                         imaxloop = ic;
 575                                                 }
 576                                         }
 577                                         icend = imaxloop;
 578                                         icstart = imaxloop + 1;
 579                                         if (icstart >= lp.getNconnection())
 580                                                 icstart = 0;
 581                                         cp = lp.getConnection(icend);
 582                                         cp.setBroken(true);
 583                                         done = true;
 584                                 }
 585                         } while (!done);
 586                         done_all_connections = false;
 587                         icstart1 = icstart;
 588                         if (debug)
 589                                 System.out.printf("  Icstart1 = %d\n", icstart1);
 590                         while (!done_all_connections) {
 591                                 count = 0;
 592                                 done = false;
 593                                 icend = icstart;
 594                                 rooted = false;
 595                                 while (!done) {
 596                                         cp = lp.getConnection(icend);
 597                                         if (icend == icroot)
 598                                                 rooted = true;
 599                                         j = icend + 1;
 600                                         if (j >= lp.getNconnection()) {
 601                                                 j = 0;
 602                                         }
 603                                         cpnext = lp.getConnection(j);
 604                                         if (connected_connection(cp, cpnext)) {
 605                                                 if (++count >= lp.getNconnection())
 606                                                         break;
 607                                                 icend = j;
 608                                         } else {
 609                                                 done = true;
 610                                         }
 611                                 }
 612                                 icmiddle = find_ic_middle(icstart, icend, anchor_connection,
 613                                                 acp, lp);
 614                                 ic = icup = icdown = icmiddle;
 615                                 if (debug)
 616                                         System.out.printf("  IC start = %d  middle = %d  end = %d\n",
 617                                                         icstart, icmiddle, icend);
 618                                 done = false;
 619                                 direction = 0;
 620                                 while (!done) {
 621                                         if (direction < 0) {
 622                                                 ic = icup;
 623                                         } else if (direction == 0) {
 624                                                 ic = icmiddle;
 625                                         } else {
 626                                                 ic = icdown;
 627                                         }
 628                                         if (ic >= 0) {
 629                                                 cp = lp.getConnection(ic);
 630                                                 if (anchor_connection == null || acp != cp) {
 631                                                         if (direction == 0) {
 632                                                                 astart = cp.getAngle()
 633                                                                                 - Math.asin(1.0 / 2.0 / radius);
 634                                                                 aend = cp.getAngle()
 635                                                                                 + Math.asin(1.0 / 2.0 / radius);
 636                                                                 bases.get(cp.getStart()).setX(
 637                                                                                 xc + radius * Math.cos(astart));
 638                                                                 bases.get(cp.getStart()).setY(
 639                                                                                 yc + radius * Math.sin(astart));
 640                                                                 bases.get(cp.getEnd()).setX(
 641                                                                                 xc + radius * Math.cos(aend));
 642                                                                 bases.get(cp.getEnd()).setY(
 643                                                                                 yc + radius * Math.sin(aend));
 644                                                         } else if (direction < 0) {
 645                                                                 j = ic + 1;
 646                                                                 if (j >= lp.getNconnection())
 647                                                                         j = 0;
 648                                                                 cp = lp.getConnection(ic);
 649                                                                 cpnext = lp.getConnection(j);
 650                                                                 cpx = cp.getXrad();
 651                                                                 cpy = cp.getYrad();
 652                                                                 ac = (cp.getAngle() + cpnext.getAngle()) / 2.0;
 653                                                                 if (cp.getAngle() > cpnext.getAngle())
 654                                                                         ac -= Math.PI;
 655                                                                 cnx = Math.cos(ac);
 656                                                                 cny = Math.sin(ac);
 657                                                                 lnx = cny;
 658                                                                 lny = -cnx;
 659                                                                 da = cpnext.getAngle() - cp.getAngle();
 660                                                                 if (da < 0.0)
 661                                                                         da += 2 * Math.PI;
 662                                                                 if (cp.isExtruded()) {
 663                                                                         if (da <= Math.PI / 2)
 664                                                                                 rl = 2.0;
 665                                                                         else
 666                                                                                 rl = 1.5;
 667                                                                 } else {
 668                                                                         rl = 1.0;
 669                                                                 }
 670                                                                 bases.get(cp.getEnd()).setX(
 671                                                                                 bases.get(cpnext.getStart()).getX()
 672                                                                                                 + rl * lnx);
 673                                                                 bases.get(cp.getEnd()).setY(
 674                                                                                 bases.get(cpnext.getStart()).getY()
 675                                                                                                 + rl * lny);
 676                                                                 bases.get(cp.getStart()).setX(
 677                                                                                 bases.get(cp.getEnd()).getX() + cpy);
 678                                                                 bases.get(cp.getStart()).setY(
 679                                                                                 bases.get(cp.getEnd()).getY() - cpx);
 680                                                         } else {
 681                                                                 j = ic - 1;
 682                                                                 if (j < 0)
 683                                                                         j = lp.getNconnection() - 1;
 684                                                                 cp = lp.getConnection(j);
 685                                                                 cpnext = lp.getConnection(ic);
 686                                                                 cpnextx = cpnext.getXrad();
 687                                                                 cpnexty = cpnext.getYrad();
 688                                                                 ac = (cp.getAngle() + cpnext.getAngle()) / 2.0;
 689                                                                 if (cp.getAngle() > cpnext.getAngle())
 690                                                                         ac -= Math.PI;
 691                                                                 cnx = Math.cos(ac);
 692                                                                 cny = Math.sin(ac);
 693                                                                 lnx = -cny;
 694                                                                 lny = cnx;
 695                                                                 da = cpnext.getAngle() - cp.getAngle();
 696                                                                 if (da < 0.0)
 697                                                                         da += 2 * Math.PI;
 698                                                                 if (cp.isExtruded()) {
 699                                                                         if (da <= Math.PI / 2)
 700                                                                                 rl = 2.0;
 701                                                                         else
 702                                                                                 rl = 1.5;
 703                                                                 } else {
 704                                                                         rl = 1.0;
 705                                                                 }
 706                                                                 bases.get(cpnext.getStart()).setX(
 707                                                                                 bases.get(cp.getEnd()).getX() + rl
 708                                                                                                 * lnx);
 709                                                                 bases.get(cpnext.getStart()).setY(
 710                                                                                 bases.get(cp.getEnd()).getY() + rl
 711                                                                                                 * lny);
 712                                                                 bases.get(cpnext.getEnd()).setX(
 713                                                                                 bases.get(cpnext.getStart()).getX()
 714                                                                                                 - cpnexty);
 715                                                                 bases.get(cpnext.getEnd()).setY(
 716                                                                                 bases.get(cpnext.getStart()).getY()
 717                                                                                                 + cpnextx);
 718                                                         }
 719                                                 }
 720                                         }
 721                                         if (direction < 0) {
 722                                                 if (icdown == icend) {
 723                                                         icdown = -1;
 724                                                 } else if (icdown >= 0) {
 725                                                         if (++icdown >= lp.getNconnection()) {
 726                                                                 icdown = 0;
 727                                                         }
 728                                                 }
 729                                                 direction = 1;
 730                                         } else {
 731                                                 if (icup == icstart)
 732                                                         icup = -1;
 733                                                 else if (icup >= 0) {
 734                                                         if (--icup < 0) {
 735                                                                 icup = lp.getNconnection() - 1;
 736                                                         }
 737                                                 }
 738                                                 direction = -1;
 739                                         }
 740                                         done = icup == -1 && icdown == -1;
 741                                 }
 742                                 icnext = icend + 1;
 743                                 if (icnext >= lp.getNconnection())
 744                                         icnext = 0;
 745                                 if (icend != icstart
 746                                                 && (!(icstart == icstart1 && icnext == icstart1))) {
 747
 748                                         // Move the bases just constructed (or the radius) so that
 749                                         // the bisector of the end points is radius distance away
 750                                         // from the loop center.
 751
 752                                         cp = lp.getConnection(icstart);
 753                                         cpnext = lp.getConnection(icend);
 754                                         dx = bases.get(cpnext.getEnd()).getX()
 755                                                         - bases.get(cp.getStart()).getX();
 756                                         dy = bases.get(cpnext.getEnd()).getY()
 757                                                         - bases.get(cp.getStart()).getY();
 758                                         midx = bases.get(cp.getStart()).getX() + dx / 2.0;
 759                                         midy = bases.get(cp.getStart()).getY() + dy / 2.0;
 760                                         rr = Math.sqrt(dx * dx + dy * dy);
 761                                         mx = dx / rr;
 762                                         my = dy / rr;
 763                                         vx = xc - midx;
 764                                         vy = yc - midy;
 765                                         rr = Math.sqrt(dx * dx + dy * dy);
 766                                         vx /= rr;
 767                                         vy /= rr;
 768                                         dotmv = vx * mx + vy * my;
 769                                         nrx = dotmv * mx - vx;
 770                                         nry = dotmv * my - vy;
 771                                         rr = Math.sqrt(nrx * nrx + nry * nry);
 772                                         nrx /= rr;
 773                                         nry /= rr;
 774
 775                                         // Determine which side of the bisector the center should
 776                                         // be.
 777
 778                                         dx = bases.get(cp.getStart()).getX() - xc;
 779                                         dy = bases.get(cp.getStart()).getY() - yc;
 780                                         ac = Math.atan2(dy, dx);
 781                                         if (ac < 0.0)
 782                                                 ac += 2 * Math.PI;
 783                                         dx = bases.get(cpnext.getEnd()).getX() - xc;
 784                                         dy = bases.get(cpnext.getEnd()).getY() - yc;
 785                                         acn = Math.atan2(dy, dx);
 786                                         if (acn < 0.0)
 787                                                 acn += 2 * Math.PI;
 788                                         if (acn < ac)
 789                                                 acn += 2 * Math.PI;
 790                                         if (acn - ac > Math.PI)
 791                                                 sign = -1;
 792                                         else
 793                                                 sign = 1;
 794                                         nmidx = xc + sign * radius * nrx;
 795                                         nmidy = yc + sign * radius * nry;
 796                                         if (rooted) {
 797                                                 xc -= nmidx - midx;
 798                                                 yc -= nmidy - midy;
 799                                         } else {
 800                                                 for (ic = icstart;;) {
 801                                                         cp = lp.getConnection(ic);
 802                                                         i = cp.getStart();
 803                                                         bases.get(i).setX(
 804                                                                         bases.get(i).getX() + nmidx - midx);
 805                                                         bases.get(i).setY(
 806                                                                         bases.get(i).getY() + nmidy - midy);
 807                                                         i = cp.getEnd();
 808                                                         bases.get(i).setX(
 809                                                                         bases.get(i).getX() + nmidx - midx);
 810                                                         bases.get(i).setY(
 811                                                                         bases.get(i).getY() + nmidy - midy);
 812                                                         if (ic == icend)
 813                                                                 break;
 814                                                         if (++ic >= lp.getNconnection())
 815                                                                 ic = 0;
 816                                                 }
 817                                         }
 818                                 }
 819                                 icstart = icnext;
 820                                 done_all_connections = icstart == icstart1;
 821                         }
 822                         for (ic = 0; ic < lp.getNconnection(); ic++) {
 823                                 cp = lp.getConnection(ic);
 824                                 j = ic + 1;
 825                                 if (j >= lp.getNconnection())
 826                                         j = 0;
 827                                 cpnext = lp.getConnection(j);
 828                                 dx = bases.get(cp.getEnd()).getX() - xc;
 829                                 dy = bases.get(cp.getEnd()).getY() - yc;
 830                                 rc = Math.sqrt(dx * dx + dy * dy);
 831                                 ac = Math.atan2(dy, dx);
 832                                 if (ac < 0.0)
 833                                         ac += 2 * Math.PI;
 834                                 dx = bases.get(cpnext.getStart()).getX() - xc;
 835                                 dy = bases.get(cpnext.getStart()).getY() - yc;
 836                                 rcn = Math.sqrt(dx * dx + dy * dy);
 837                                 acn = Math.atan2(dy, dx);
 838                                 if (acn < 0.0)
 839                                         acn += 2 * Math.PI;
 840                                 if (acn < ac)
 841                                         acn += 2 * Math.PI;
 842                                 dan = acn - ac;
 843                                 dcp = cpnext.getAngle() - cp.getAngle();
 844                                 if (dcp <= 0.0)
 845                                         dcp += 2 * Math.PI;
 846                                 if (Math.abs(dan - dcp) > Math.PI) {
 847                                         if (cp.isExtruded()) {
 848                                                 warningEmition("Warning from traverse_loop. Loop "
 849                                                                 + lp.getNumber() + " has crossed regions\n");
 850                                         } else if ((cpnext.getStart() - cp.getEnd()) != 1) {
 851                                                 cp.setExtruded(true);
 852                                                 continue set_radius; // remplacement du goto
 853                                         }
 854                                 }
 855                                 if (cp.isExtruded()) {
 856                                         construct_extruded_segment(cp, cpnext);
 857                                 } else {
 858                                         n = cpnext.getStart() - cp.getEnd();
 859                                         if (n < 0)
 860                                                 n += nbase + 1;
 861                                         angleinc = dan / n;
 862                                         for (j = 1; j < n; j++) {
 863                                                 i = cp.getEnd() + j;
 864                                                 if (i > nbase)
 865                                                         i -= nbase + 1;
 866                                                 a = ac + j * angleinc;
 867                                                 rr = rc + (rcn - rc) * (a - ac) / dan;
 868                                                 bases.get(i).setX(xc + rr * Math.cos(a));
 869                                                 bases.get(i).setY(yc + rr * Math.sin(a));
 870                                         }
 871                                 }
 872                         }
 873                         break;
 874                 }
 875                 for (ic = 0; ic < lp.getNconnection(); ic++) {
 876                         if (icroot != ic) {
 877                                 cp = lp.getConnection(ic);
 878                                 generate_region(cp);
 879                                 traverse_loop(cp.getLoop(), cp);
 880                         }
 881                 }
 882                 n = 0;
 883                 sx = 0.0;
 884                 sy = 0.0;
 885                 for (ic = 0; ic < lp.getNconnection(); ic++) {
 886                         j = ic + 1;
 887                         if (j >= lp.getNconnection())
 888                                 j = 0;
 889                         cp = lp.getConnection(ic);
 890                         cpnext = lp.getConnection(j);
 891                         n += 2;
 892                         sx += bases.get(cp.getStart()).getX()
 893                                         + bases.get(cp.getEnd()).getX();
 894                         sy += bases.get(cp.getStart()).getY()
 895                                         + bases.get(cp.getEnd()).getY();
 896                         if (!cp.isExtruded()) {
 897                                 for (j = cp.getEnd() + 1; j != cpnext.getStart(); j++) {
 898                                         if (j > nbase)
 899                                                 j -= nbase + 1;
 900                                         n++;
 901                                         sx += bases.get(j).getX();
 902                                         sy += bases.get(j).getY();
 903                                 }
 904                         }
 905                 }
 906                 lp.setX(sx / n);
 907                 lp.setY(sy / n);
 908         }
 909
 910         /**
 911          * For the loop pointed to by lp, determine the radius of the loop that will
 912          * ensure that each base around the loop will have a separation of at least
 913          * lencut around the circle. If a segment joining two connectors will not
 914          * support this separation, then the flag, extruded, will be set in the
 915          * first of these two indicators. The radius is set in lp.
 916          *
 917          * The radius is selected by a least squares procedure where the sum of the
 918          * squares of the deviations of length from the ideal value of 1 is used as
 919          * the error function.
 920          */
 921         private void determine_radius(Loop lp, double lencut) {
 922                 if (debug)
 923                         System.out.println("  Determine_radius");
 924                 double mindit, ci, dt, sumn, sumd, radius, dit;
 925                 int i, j, end, start, imindit = 0;
 926                 Connection cp = new Connection(), cpnext = new Connection();
 927                 double rt2_2 = 0.7071068;
 928
 929                 do {
 930                         mindit = 1.0e10;
 931                         for (sumd = 0.0, sumn = 0.0, i = 0; i < lp.getNconnection(); i++) {
 932                                 cp = lp.getConnection(i);
 933                                 j = i + 1;
 934                                 if (j >= lp.getNconnection())
 935                                         j = 0;
 936                                 cpnext = lp.getConnection(j);
 937                                 end = cp.getEnd();
 938                                 start = cpnext.getStart();
 939                                 if (start < end)
 940                                         start += nbase + 1;
 941                                 dt = cpnext.getAngle() - cp.getAngle();
 942                                 if (dt <= 0.0)
 943                                         dt += 2 * Math.PI;
 944                                 if (!cp.isExtruded())
 945                                         ci = start - end;
 946                                 else {
 947                                         if (dt <= Math.PI / 2)
 948                                                 ci = 2.0;
 949                                         else
 950                                                 ci = 1.5;
 951                                 }
 952                                 sumn += dt * (1.0 / ci + 1.0);
 953                                 sumd += dt * dt / ci;
 954                                 dit = dt / ci;
 955                                 if (dit < mindit && !cp.isExtruded() && ci > 1.0) {
 956                                         mindit = dit;
 957                                         imindit = i;
 958                                 }
 959                         }
 960                         radius = sumn / sumd;
 961                         if (radius < rt2_2)
 962                                 radius = rt2_2;
 963                         if (mindit * radius < lencut) {
 964                                 lp.getConnection(imindit).setExtruded(true);
 965                         }
 966                 } while (mindit * radius < lencut);
 967                 if (lp.getRadius() > 0.0)
 968                         radius = lp.getRadius();
 969                 else
 970                         lp.setRadius(radius);
 971         }
 972
 973         /**
 974          * Determines if the connections cp and cpnext are connected
 975          */
 976         private boolean connected_connection(Connection cp, Connection cpnext) {
 977                 if (debug)
 978                         System.out.println("  Connected_connection");
 979                 if (cp.isExtruded()) {
 980                         return true;
 981                 } else if (cp.getEnd() + 1 == cpnext.getStart()) {
 982                         return true;
 983                 } else {
 984                         return false;
 985                 }
 986         }
 987
 988         /**
 989          * Finds the middle of a set of connected connectors. This is normally the
 990          * middle connection in the sequence except if one of the connections is the
 991          * anchor, in which case that connection will be used.
 992          *
 993          * @throws ExceptionNAViewAlgorithm
 994          */
 995         private int find_ic_middle(int icstart, int icend,
 996                         Connection anchor_connection, Connection acp, Loop lp)
 997                         throws ExceptionNAViewAlgorithm {
 998                 if (debug)
 999                         System.out.println("  Find_ic_middle");
1000                 int count, ret, ic, i;
1001                 boolean done;
1002
1003                 count = 0;
1004                 ret = -1;
1005                 ic = icstart;
1006                 done = false;
1007                 while (!done) {
1008                         if (count++ > lp.getNconnection() * 2) {
1009                                 throw new ExceptionNAViewAlgorithm(
1010                                                 "Infinite loop detected in find_ic_middle");
1011                         }
1012                         if (anchor_connection != null && lp.getConnection(ic) == acp) {
1013                                 ret = ic;
1014                         }
1015                         done = ic == icend;
1016                         if (++ic >= lp.getNconnection()) {
1017                                 ic = 0;
1018                         }
1019                 }
1020                 if (ret == -1) {
1021                         for (i = 1, ic = icstart; i < (count + 1) / 2; i++) {
1022                                 if (++ic >= lp.getNconnection())
1023                                         ic = 0;
1024                         }
1025                         ret = ic;
1026                 }
1027                 return ret;
1028         }
1029
1030         /**
1031          * Generates the coordinates for the base pairing region of a connection
1032          * given the position of the starting base pair.
1033          *
1034          * @throws ExceptionNAViewAlgorithm
1035          */
1036         private void generate_region(Connection cp) throws ExceptionNAViewAlgorithm {
1037                 if (debug)
1038                         System.out.println("  Generate_region");
1039                 int l, start, end, i, mate;
1040                 Region rp;
1041
1042                 rp = cp.getRegion();
1043                 l = 0;
1044                 if (cp.getStart() == rp.getStart1()) {
1045                         start = rp.getStart1();
1046                         end = rp.getEnd1();
1047                 } else {
1048                         start = rp.getStart2();
1049                         end = rp.getEnd2();
1050                 }
1051                 if (bases.get(cp.getStart()).getX() > ANUM - 100.0
1052                                 || bases.get(cp.getEnd()).getX() > ANUM - 100.0) {
1053                         throw new ExceptionNAViewAlgorithm(
1054                                         "Bad region passed to generate_region. Coordinates not defined.");
1055                 }
1056                 for (i = start + 1; i <= end; i++) {
1057                         l++;
1058                         bases.get(i).setX(
1059                                         bases.get(cp.getStart()).getX() + HELIX_FACTOR * l
1060                                                         * cp.getXrad());
1061                         bases.get(i).setY(
1062                                         bases.get(cp.getStart()).getY() + HELIX_FACTOR * l
1063                                                         * cp.getYrad());
1064                         mate = bases.get(i).getMate();
1065                         bases.get(mate).setX(
1066                                         bases.get(cp.getEnd()).getX() + HELIX_FACTOR * l
1067                                                         * cp.getXrad());
1068                         bases.get(mate).setY(
1069                                         bases.get(cp.getEnd()).getY() + HELIX_FACTOR * l
1070                                                         * cp.getYrad());
1071
1072                 }
1073         }
1074
1075         /**
1076          * Draws the segment of residue between the bases numbered start through
1077          * end, where start and end are presumed to be part of a base pairing
1078          * region. They are drawn as a circle which has a chord given by the ends of
1079          * two base pairing regions defined by the connections.
1080          *
1081          * @throws ExceptionNAViewAlgorithm
1082          */
1083         private void construct_circle_segment(int start, int end)
1084                         throws ExceptionNAViewAlgorithm {
1085                 if (debug)
1086                         System.out.println("  Construct_circle_segment");
1087                 double dx, dy, rr, midx, midy, xn, yn, nrx, nry, mx, my, a;
1088                 int l, j, i;
1089
1090                 dx = bases.get(end).getX() - bases.get(start).getX();
1091                 dy = bases.get(end).getY() - bases.get(start).getY();
1092                 rr = Math.sqrt(dx * dx + dy * dy);
1093                 l = end - start;
1094                 if (l < 0)
1095                         l += nbase + 1;
1096                 if (rr >= l) {
1097                         dx /= rr;
1098                         dy /= rr;
1099                         for (j = 1; j < l; j++) {
1100                                 i = start + j;
1101                                 if (i > nbase)
1102                                         i -= nbase + 1;
1103                                 bases.get(i).setX(
1104                                                 bases.get(start).getX() + dx * (double) j / (double) l);
1105                                 bases.get(i).setY(
1106                                                 bases.get(start).getY() + dy * (double) j / (double) l);
1107                         }
1108                 } else {
1109                         find_center_for_arc((l - 1), rr);
1110                         dx /= rr;
1111                         dy /= rr;
1112                         midx = bases.get(start).getX() + dx * rr / 2.0;
1113                         midy = bases.get(start).getY() + dy * rr / 2.0;
1114                         xn = dy;
1115                         yn = -dx;
1116                         nrx = midx + _h * xn;
1117                         nry = midy + _h * yn;
1118                         mx = bases.get(start).getX() - nrx;
1119                         my = bases.get(start).getY() - nry;
1120                         rr = Math.sqrt(mx * mx + my * my);
1121                         a = Math.atan2(my, mx);
1122                         for (j = 1; j < l; j++) {
1123                                 i = start + j;
1124                                 if (i > nbase)
1125                                         i -= nbase + 1;
1126                                 bases.get(i).setX(nrx + rr * Math.cos(a + j * angleinc));
1127                                 bases.get(i).setY(nry + rr * Math.sin(a + j * angleinc));
1128                         }
1129                 }
1130         }
1131
1132         /**
1133          * Constructs the segment between cp and cpnext as a circle if possible.
1134          * However, if the segment is too large, the lines are drawn between the two
1135          * connecting regions, and bases are placed there until the connecting
1136          * circle will fit.
1137          *
1138          * @throws ExceptionNAViewAlgorithm
1139          */
1140         private void construct_extruded_segment(Connection cp, Connection cpnext)
1141                         throws ExceptionNAViewAlgorithm {
1142                 if (debug)
1143                         System.out.println("  Construct_extruded_segment");
1144                 double astart, aend1, aend2, aave, dx, dy, a1, a2, ac, rr, da, dac;
1145                 int start, end, n, nstart, nend;
1146                 boolean collision;
1147
1148                 astart = cp.getAngle();
1149                 aend2 = aend1 = cpnext.getAngle();
1150                 if (aend2 < astart)
1151                         aend2 += 2 * Math.PI;
1152                 aave = (astart + aend2) / 2.0;
1153                 start = cp.getEnd();
1154                 end = cpnext.getStart();
1155                 n = end - start;
1156                 if (n < 0)
1157                         n += nbase + 1;
1158                 da = cpnext.getAngle() - cp.getAngle();
1159                 if (da < 0.0) {
1160                         da += 2 * Math.PI;
1161                 }
1162                 if (n == 2)
1163                         construct_circle_segment(start, end);
1164                 else {
1165                         dx = bases.get(end).getX() - bases.get(start).getX();
1166                         dy = bases.get(end).getY() - bases.get(start).getY();
1167                         rr = Math.sqrt(dx * dx + dy * dy);
1168                         dx /= rr;
1169                         dy /= rr;
1170                         if (rr >= 1.5 && da <= Math.PI / 2) {
1171                                 nstart = start + 1;
1172                                 if (nstart > nbase)
1173                                         nstart -= nbase + 1;
1174                                 nend = end - 1;
1175                                 if (nend < 0)
1176                                         nend += nbase + 1;
1177                                 bases.get(nstart).setX(bases.get(start).getX() + 0.5 * dx);
1178                                 bases.get(nstart).setY(bases.get(start).getY() + 0.5 * dy);
1179                                 bases.get(nend).setX(bases.get(end).getX() - 0.5 * dx);
1180                                 bases.get(nend).setY(bases.get(end).getY() - 0.5 * dy);
1181                                 start = nstart;
1182                                 end = nend;
1183                         }
1184                         do {
1185                                 collision = false;
1186                                 construct_circle_segment(start, end);
1187                                 nstart = start + 1;
1188                                 if (nstart > nbase)
1189                                         nstart -= nbase + 1;
1190                                 dx = bases.get(nstart).getX() - bases.get(start).getX();
1191                                 dy = bases.get(nstart).getY() - bases.get(start).getY();
1192                                 a1 = Math.atan2(dy, dx);
1193                                 if (a1 < 0.0)
1194                                         a1 += 2 * Math.PI;
1195                                 dac = a1 - astart;
1196                                 if (dac < 0.0)
1197                                         dac += 2 * Math.PI;
1198                                 if (dac > Math.PI)
1199                                         collision = true;
1200                                 nend = end - 1;
1201                                 if (nend < 0)
1202                                         nend += nbase + 1;
1203                                 dx = bases.get(nend).getX() - bases.get(end).getX();
1204                                 dy = bases.get(nend).getY() - bases.get(end).getY();
1205                                 a2 = Math.atan2(dy, dx);
1206                                 if (a2 < 0.0)
1207                                         a2 += 2 * Math.PI;
1208                                 dac = aend1 - a2;
1209                                 if (dac < 0.0)
1210                                         dac += 2 * Math.PI;
1211                                 if (dac > Math.PI)
1212                                         collision = true;
1213                                 if (collision) {
1214                                         ac = minf2(aave, astart + 0.5);
1215                                         bases.get(nstart).setX(
1216                                                         bases.get(start).getX() + Math.cos(ac));
1217                                         bases.get(nstart).setY(
1218                                                         bases.get(start).getY() + Math.sin(ac));
1219                                         start = nstart;
1220                                         ac = maxf2(aave, aend2 - 0.5);
1221                                         bases.get(nend).setX(bases.get(end).getX() + Math.cos(ac));
1222                                         bases.get(nend).setY(bases.get(end).getY() + Math.sin(ac));
1223                                         end = nend;
1224                                         n -= 2;
1225                                 }
1226                         } while (collision && n > 1);
1227                 }
1228         }
1229
1230         /**
1231          * Given n points to be placed equidistantly and equiangularly on a polygon
1232          * which has a chord of length, b, find the distance, h, from the midpoint
1233          * of the chord for the center of polygon. Positive values mean the center
1234          * is within the polygon and the chord, whereas negative values mean the
1235          * center is outside the chord. Also, the radial angle for each polygon side
1236          * is returned in theta.
1237          *
1238          * The procedure uses a bisection algorithm to find the correct value for
1239          * the center. Two equations are solved, the angles around the center must
1240          * add to 2*Math.PI, and the sides of the polygon excluding the chord must
1241          * have a length of 1.
1242          *
1243          * @throws ExceptionNAViewAlgorithm
1244          */
1245         private void find_center_for_arc(double n, double b)
1246                         throws ExceptionNAViewAlgorithm {
1247                 if (debug)
1248                         System.out.println("  Find_center_for_arc");
1249                 double h, hhi, hlow, r, disc, theta, e, phi;
1250                 int iter;
1251
1252                 hhi = (n + 1.0) / Math.PI;
1253                 // changed to prevent div by zero if (ih)
1254                 hlow = -hhi - b / (n + 1.000001 - b);
1255                 if (b < 1)
1256                         // otherwise we might fail below (ih)
1257                         hlow = 0;
1258                 iter = 0;
1259                 do {
1260                         h = (hhi + hlow) / 2.0;
1261                         r = Math.sqrt(h * h + b * b / 4.0);
1262                         // if (r<0.5) {r = 0.5; h = 0.5*Math.sqrt(1-b*b);}
1263                         disc = 1.0 - 0.5 / (r * r);
1264                         if (Math.abs(disc) > 1.0) {
1265                                 throw new ExceptionNAViewAlgorithm(
1266                                                 "Unexpected large magnitude discriminant = " + disc
1267                                                                 + " " + r);
1268                         }
1269                         theta = Math.acos(disc);
1270                         // theta = 2*Math.acos(Math.sqrt(1-1/(4*r*r)));
1271                         phi = Math.acos(h / r);
1272                         e = theta * (n + 1) + 2 * phi - 2 * Math.PI;
1273                         if (e > 0.0) {
1274                                 hlow = h;
1275                         } else {
1276                                 hhi = h;
1277                         }
1278                 } while (Math.abs(e) > 0.0001 && ++iter < MAXITER);
1279                 if (iter >= MAXITER) {
1280                         if (noIterationFailureYet) {
1281                                 warningEmition("Iteration failed in find_center_for_arc");
1282                                 noIterationFailureYet = false;
1283                         }
1284                         h = 0.0;
1285                         theta = 0.0;
1286                 }
1287                 _h = h;
1288                 angleinc = theta;
1289         }
1290
1291         private double minf2(double x1, double x2) {
1292                 return ((x1) < (x2)) ? (x1) : (x2);
1293         }
1294
1295         private double maxf2(double x1, double x2) {
1296                 return ((x1) > (x2)) ? (x1) : (x2);
1297         }
1298
1299         public void warningEmition(String warningMessage) throws ExceptionNAViewAlgorithm {
1300                 throw (new ExceptionNAViewAlgorithm(warningMessage));
1301         }
1302 }