001/*- 002 ******************************************************************************* 003 * Copyright (c) 2011, 2016 Diamond Light Source Ltd. 004 * All rights reserved. This program and the accompanying materials 005 * are made available under the terms of the Eclipse Public License v1.0 006 * which accompanies this distribution, and is available at 007 * http://www.eclipse.org/legal/epl-v10.html 008 * 009 * Contributors: 010 * Peter Chang - initial API and implementation and/or initial documentation 011 *******************************************************************************/ 012 013// GEN_COMMENT 014 015package org.eclipse.january.dataset; 016 017 018import java.util.Arrays; 019 020import org.apache.commons.math3.complex.Complex; 021 022 023/** 024 * Extend compound dataset to hold complex double values // PRIM_TYPE 025 */ 026public class ComplexDoubleDataset extends CompoundDoubleDataset { // CLASS_TYPE 027 // pin UID to base class 028 private static final long serialVersionUID = Dataset.serialVersionUID; 029 030 private static final int ISIZE = 2; // number of elements per item 031 032 @Override 033 public int getDType() { 034 return Dataset.COMPLEX128; // DATA_TYPE 035 } 036 037 /** 038 * Create a null dataset 039 */ 040 ComplexDoubleDataset() { 041 super(ISIZE); 042 } 043 044 /** 045 * Create a zero-filled dataset of given shape 046 * @param shape 047 */ 048 ComplexDoubleDataset(final int... shape) { 049 super(ISIZE, shape); 050 } 051 052 /** 053 * Create a dataset using given data (real and imaginary parts are grouped in pairs) 054 * @param data 055 * @param shape (can be null to create 1D dataset) 056 */ 057 ComplexDoubleDataset(final double[] data, final int... shape) { // PRIM_TYPE 058 super(ISIZE, data, shape); 059 } 060 061 /** 062 * Copy a dataset 063 * @param dataset 064 */ 065 ComplexDoubleDataset(final ComplexDoubleDataset dataset) { 066 super(dataset); 067 } 068 069 /** 070 * Create a dataset using given data (real and imaginary parts are given separately) 071 * @param realData 072 * @param imagData 073 * @param shape (can be null or zero-length to create 1D dataset) 074 */ 075 ComplexDoubleDataset(final double[] realData, final double[] imagData, int... shape) { // PRIM_TYPE 076 if (realData == null || imagData == null) { 077 throw new IllegalArgumentException("Data must not be null"); 078 } 079 int dsize = realData.length > imagData.length ? imagData.length : realData.length; 080 if (shape == null || shape.length == 0) { 081 shape = new int[] {dsize}; 082 } 083 isize = ISIZE; 084 size = ShapeUtils.calcSize(shape); 085 if (size != dsize) { 086 throw new IllegalArgumentException(String.format("Shape %s is not compatible with size of data array, %d", 087 Arrays.toString(shape), dsize)); 088 } 089 this.shape = size == 0 ? null : shape.clone(); 090 091 try { 092 odata = data = createArray(size); 093 } catch (Throwable t) { 094 logger.error("Could not create a dataset of shape {}", Arrays.toString(shape), t); 095 throw new IllegalArgumentException(t); 096 } 097 098 for (int i = 0, n = 0; i < size; i++) { 099 data[n++] = realData[i]; 100 data[n++] = imagData[i]; 101 } 102 } 103 104 /** 105 * Create a dataset using given data (real and imaginary parts are given separately) 106 * @param real 107 * @param imag 108 */ 109 ComplexDoubleDataset(final Dataset real, final Dataset imag) { 110 super(ISIZE, real.getShapeRef()); 111 real.checkCompatibility(imag); 112 113 IndexIterator riter = real.getIterator(); 114 IndexIterator iiter = imag.getIterator(); 115 116 for (int i = 0; riter.hasNext() && iiter.hasNext();) { 117 data[i++] = real.getElementDoubleAbs(riter.index); // ADD_CAST 118 data[i++] = imag.getElementDoubleAbs(iiter.index); // ADD_CAST 119 } 120 } 121 122 /** 123 * Copy and cast a dataset to this complex type 124 * @param dataset 125 */ 126 ComplexDoubleDataset(final Dataset dataset) { 127 super(ISIZE, dataset.getShapeRef()); 128 copyToView(dataset, this, true, false); 129 offset = 0; 130 stride = null; 131 base = null; 132 try { 133 odata = data = createArray(size); 134 } catch (Throwable t) { 135 logger.error("Could not create a dataset of shape {}", Arrays.toString(shape), t); 136 throw new IllegalArgumentException(t); 137 } 138 139 IndexIterator iter = dataset.getIterator(); 140 if (dataset.isComplex()) { 141 for (int i = 0; iter.hasNext(); i += isize) { 142 data[i] = dataset.getElementDoubleAbs(iter.index); // ADD_CAST 143 data[i+1] = dataset.getElementDoubleAbs(iter.index+1); // ADD_CAST 144 } 145 } else { 146 for (int i = 0; iter.hasNext(); i += isize) { 147 data[i] = dataset.getElementDoubleAbs(iter.index); // ADD_CAST 148 } 149 } 150 } 151 152 @Override 153 public ComplexDoubleDataset clone() { 154 return new ComplexDoubleDataset(this); 155 } 156 157 /** 158 * Create a dataset from an object which could be a Java list, array (of arrays...) 159 * or Number. Ragged sequences or arrays are padded with zeros. 160 * 161 * @param obj 162 * @return dataset with contents given by input 163 */ 164 static ComplexDoubleDataset createFromObject(final Object obj) { 165 ComplexDoubleDataset result = new ComplexDoubleDataset(); 166 167 result.shape = ShapeUtils.getShapeFromObject(obj); 168 result.size = ShapeUtils.calcSize(result.shape); 169 170 try { 171 result.odata = result.data = result.createArray(result.size); 172 } catch (Throwable t) { 173 logger.error("Could not create a dataset of shape {}", Arrays.toString(result.shape), t); 174 throw new IllegalArgumentException(t); 175 } 176 177 int[] pos = new int[result.shape.length]; 178 result.fillData(obj, 0, pos); 179 return result; 180 } 181 182 /** 183 * @param stop 184 * @return a new 1D dataset, filled with values determined by parameters 185 */ 186 static ComplexDoubleDataset createRange(final double stop) { 187 return createRange(0, stop, 1); 188 } 189 190 /** 191 * @param start 192 * @param stop 193 * @param step 194 * @return a new 1D dataset, filled with values determined by parameters 195 */ 196 static ComplexDoubleDataset createRange(final double start, final double stop, final double step) { 197 int size = calcSteps(start, stop, step); 198 ComplexDoubleDataset result = new ComplexDoubleDataset(size); 199 for (int i = 0; i < size; i ++) { 200 result.data[i*ISIZE] = (start + i*step); // ADD_CAST 201 } 202 return result; 203 } 204 205 /** 206 * @param shape 207 * @return a dataset filled with ones 208 */ 209 static ComplexDoubleDataset ones(final int... shape) { 210 return new ComplexDoubleDataset(shape).fill(1); 211 } 212 213 @Override 214 public ComplexDoubleDataset fill(final Object obj) { 215 setDirty(); 216 double vr = DTypeUtils.toReal(obj); // PRIM_TYPE // ADD_CAST 217 double vi = DTypeUtils.toImag(obj); // PRIM_TYPE // ADD_CAST 218 IndexIterator iter = getIterator(); 219 220 while (iter.hasNext()) { 221 data[iter.index] = vr; 222 data[iter.index+1] = vi; 223 } 224 225 return this; 226 } 227 228 @Override 229 public ComplexDoubleDataset getView(boolean deepCopyMetadata) { 230 ComplexDoubleDataset view = new ComplexDoubleDataset(); 231 copyToView(this, view, true, deepCopyMetadata); 232 view.data = data; 233 return view; 234 } 235 236 /** 237 * Get complex value at absolute index in the internal array. 238 * 239 * This is an internal method with no checks so can be dangerous. Use with care or ideally with an iterator. 240 * 241 * @param index absolute index 242 * @return value 243 */ 244 public Complex getComplexAbs(final int index) { 245 return new Complex(data[index], data[index+1]); 246 } 247 248 @Override 249 public Object getObjectAbs(final int index) { 250 return new Complex(data[index], data[index+1]); 251 } 252 253 @Override 254 public String getStringAbs(final int index) { 255 double di = data[index + 1]; // PRIM_TYPE 256 if (stringFormat == null) { 257 return di >= 0 ? String.format("%.8g + %.8gj", data[index], di) : // FORMAT_STRING 258 String.format("%.8g - %.8gj", data[index], -di); // FORMAT_STRING 259 } 260 StringBuilder s = new StringBuilder(); 261 s.append(stringFormat.format(data[index])); 262 if (di >= 0) { 263 s.append(" + "); 264 s.append(stringFormat.format(di)); 265 } else { 266 s.append(" - "); 267 s.append(stringFormat.format(-di)); 268 } 269 s.append('j'); 270 return s.toString(); 271 } 272 273 /** 274 * Set values at absolute index in the internal array. 275 * 276 * This is an internal method with no checks so can be dangerous. Use with care or ideally with an iterator. 277 * @param index absolute index 278 * @param val new values 279 */ 280 @SuppressWarnings("cast") 281 public void setAbs(final int index, final Complex val) { 282 setAbs(index, (double) val.getReal(), (double) val.getImaginary()); // PRIM_TYPE 283 } 284 285 @SuppressWarnings("cast") 286 @Override 287 public void setObjectAbs(final int index, final Object obj) { 288 setAbs(index, (double) DTypeUtils.toReal(obj), (double) DTypeUtils.toImag(obj)); // PRIM_TYPE 289 } 290 291 /** 292 * Set item at index to complex value given by real and imaginary parts 293 * @param index absolute index 294 * @param real 295 * @param imag 296 */ 297 public void setAbs(final int index, final double real, final double imag) { // PRIM_TYPE 298 setDirty(); 299 data[index] = real; 300 data[index+1] = imag; 301 } 302 303 /** 304 * @return item in first position 305 * @since 2.0 306 */ 307 public Complex get() { 308 int n = getFirst1DIndex(); 309 Complex z = new Complex(data[n], data[n+1]); 310 return z; 311 } 312 313 /** 314 * @param i 315 * @return item in given position 316 */ 317 public Complex get(final int i) { 318 int n = get1DIndex(i); 319 Complex z = new Complex(data[n], data[n+1]); 320 return z; 321 } 322 323 /** 324 * @param i 325 * @param j 326 * @return item in given position 327 */ 328 public Complex get(final int i, final int j) { 329 int n = get1DIndex(i, j); 330 Complex z = new Complex(data[n], data[n+1]); 331 return z; 332 } 333 334 /** 335 * @param pos 336 * @return item in given position 337 */ 338 public Complex get(final int... pos) { 339 int n = get1DIndex(pos); 340 Complex z = new Complex(data[n], data[n+1]); 341 return z; 342 } 343 344 @Override 345 public Object getObject() { 346 return get(); 347 } 348 349 @Override 350 public Object getObject(final int i) { 351 return get(i); 352 } 353 354 @Override 355 public Object getObject(final int i, final int j) { 356 return get(i, j); 357 } 358 359 @Override 360 public Object getObject(final int... pos) { 361 return getComplex(pos); 362 } 363 364 /** 365 * @return item in first position 366 * @since 2.0 367 */ 368 @SuppressWarnings("cast") 369 public double getReal() { // PRIM_TYPE 370 return (double) getFirstValue(); // PRIM_TYPE 371 } 372 373 /** 374 * @param i 375 * @return item in given position 376 */ 377 @SuppressWarnings("cast") 378 public double getReal(final int i) { // PRIM_TYPE 379 return (double) getFirstValue(i); // PRIM_TYPE 380 } 381 382 /** 383 * @param i 384 * @param j 385 * @return item in given position 386 */ 387 @SuppressWarnings("cast") 388 public double getReal(final int i, final int j) { // PRIM_TYPE 389 return (double) getFirstValue(i, j); // PRIM_TYPE 390 } 391 392 /** 393 * @param pos 394 * @return item in given position 395 */ 396 @SuppressWarnings("cast") 397 public double getReal(final int... pos) { // PRIM_TYPE 398 return (double) getFirstValue(pos); // PRIM_TYPE 399 } 400 401 /** 402 * @return item in first position 403 * @since 2.0 404 */ 405 public double getImag() { // PRIM_TYPE 406 return data[getFirst1DIndex() + 1]; 407 } 408 409 /** 410 * @param i 411 * @return item in given position 412 */ 413 public double getImag(final int i) { // PRIM_TYPE 414 return data[get1DIndex(i) + 1]; 415 } 416 417 /** 418 * @param i 419 * @param j 420 * @return item in given position 421 */ 422 public double getImag(final int i, final int j) { // PRIM_TYPE 423 return data[get1DIndex(i, j) + 1]; 424 } 425 426 /** 427 * @param pos 428 * @return item in given position 429 */ 430 public double getImag(final int... pos) { // PRIM_TYPE 431 return data[get1DIndex(pos) + 1]; 432 } 433 434 /** 435 * @return item in first position 436 * @since 2.0 437 */ 438 public Complex getComplex() { 439 return get(); 440 } 441 442 /** 443 * @param i 444 * @return item in given position 445 */ 446 public Complex getComplex(final int i) { 447 return get(i); 448 } 449 450 /** 451 * @param i 452 * @param j 453 * @return item in given position 454 */ 455 public Complex getComplex(final int i, final int j) { 456 return get(i, j); 457 } 458 459 /** 460 * @param pos 461 * @return item in given position 462 */ 463 public Complex getComplex(final int... pos) { 464 return get(pos); 465 } 466 467 @SuppressWarnings("cast") 468 @Override 469 public void set(final Object obj, final int i) { 470 setItem(new double[] {(double) DTypeUtils.toReal(obj), (double) DTypeUtils.toImag(obj)}, i); // PRIM_TYPE 471 } 472 473 @SuppressWarnings("cast") 474 @Override 475 public void set(final Object obj, final int i, final int j) { 476 setItem(new double[] {(double) DTypeUtils.toReal(obj), (double) DTypeUtils.toImag(obj)}, i, j); // PRIM_TYPE 477 } 478 479 @SuppressWarnings("cast") 480 @Override 481 public void set(final Object obj, int... pos) { 482 if (pos == null || (pos.length == 0 && shape.length > 0)) { 483 pos = new int[shape.length]; 484 } 485 486 setItem(new double[] {(double) DTypeUtils.toReal(obj), (double) DTypeUtils.toImag(obj)}, pos); // PRIM_TYPE 487 } 488 489 /** 490 * Set real and imaginary values at given position 491 * @param dr 492 * @param di 493 * @param i 494 */ 495 public void set(final double dr, final double di, final int i) { // PRIM_TYPE 496 setItem(new double[] {dr, di}, i); // PRIM_TYPE 497 } 498 499 /** 500 * Set real and imaginary values at given position 501 * @param dr 502 * @param di 503 * @param i 504 * @param j 505 */ 506 public void set(final double dr, final double di, final int i, final int j) { // PRIM_TYPE 507 setItem(new double[] {dr, di}, i, j); // PRIM_TYPE 508 } 509 510 /** 511 * Set real and imaginary values at given position 512 * @param dr 513 * @param di 514 * @param pos 515 * @since 2.0 516 */ 517 public void set(final double dr, final double di, final int... pos) { // PRIM_TYPE 518 setItem(new double[] {dr, di}, pos); // PRIM_TYPE 519 } 520 521 /** 522 * @since 2.0 523 */ 524 @Override 525 public DoubleDataset getRealPart() { // CLASS_TYPE 526 return getElements(0); 527 } 528 529 /** 530 * @since 2.0 531 */ 532 @Override 533 public DoubleDataset getRealView() { // CLASS_TYPE 534 return getElementsView(0); 535 } 536 537 /** 538 * @return imaginary part of dataset as new dataset 539 * @since 2.0 540 */ 541 public DoubleDataset getImaginaryPart() { // CLASS_TYPE 542 return getElements(1); 543 } 544 545 /** 546 * @return view of imaginary values 547 */ 548 public DoubleDataset getImaginaryView() { // CLASS_TYPE 549 return getElementsView(1); 550 } 551 552 @Override 553 public Number max(boolean... switches) { 554 throw new UnsupportedOperationException("Cannot compare complex numbers"); 555 } 556 557 @Override 558 public Number min(boolean... switches) { 559 throw new UnsupportedOperationException("Cannot compare complex numbers"); 560 } 561 562 @Override 563 public Object sum(boolean... switches) { // FIXME 564 double[] sum = (double[]) super.sum(switches); 565 return new Complex(sum[0], sum[1]); 566 } 567 568 @Override 569 public Object mean(boolean... switches) { 570 double[] mean = (double[]) super.mean(switches); 571 return new Complex(mean[0], mean[1]); 572 } 573 574 @Override 575 public int[] maxPos(boolean... switches) { 576 throw new UnsupportedOperationException("Cannot compare complex numbers"); 577 } 578 579 @Override 580 public int[] minPos(boolean... switches) { 581 throw new UnsupportedOperationException("Cannot compare complex numbers"); 582 } 583 584 @Override 585 public ComplexDoubleDataset getSlice(final SliceIterator siter) { 586 ComplexDoubleDataset result = new ComplexDoubleDataset(siter.getShape()); 587 double[] rdata = result.data; // PRIM_TYPE 588 IndexIterator riter = result.getIterator(); 589 590 while (siter.hasNext() && riter.hasNext()) { 591 rdata[riter.index] = data[siter.index]; 592 rdata[riter.index+1] = data[siter.index+1]; 593 } 594 595 result.setName(name + BLOCK_OPEN + Slice.createString(siter.shape, siter.start, siter.stop, siter.step) + BLOCK_CLOSE); 596 return result; 597 } 598 599 @Override 600 ComplexDoubleDataset setSlicedView(Dataset view, Dataset d) { 601 setDirty(); 602 final BroadcastSelfIterator it = BroadcastSelfIterator.createIterator(view, d); 603 604 if (d instanceof ComplexFloatDataset || d instanceof ComplexDoubleDataset) { 605 while (it.hasNext()) { 606 data[it.aIndex] = it.bDouble; // BCAST_WITH_CAST d.getElementDoubleAbs(it.bIndex); 607 data[it.aIndex + 1] = d.getElementDoubleAbs(it.bIndex + 1); // GET_ELEMENT_WITH_CAST 608 } 609 } else { 610 while (it.hasNext()) { 611 data[it.aIndex] = it.bDouble; // BCAST_WITH_CAST d.getElementDoubleAbs(it.bIndex); 612 data[it.aIndex + 1] = 0; 613 } 614 } 615 return this; 616 } 617 618 @Override 619 public ComplexDoubleDataset setSlice(final Object o, final IndexIterator siter) { 620 setDirty(); 621 if (o instanceof ComplexFloatDataset) { 622 ComplexFloatDataset zds = (ComplexFloatDataset) o; 623 624 if (!ShapeUtils.areShapesCompatible(siter.getShape(), zds.shape)) { 625 throw new IllegalArgumentException(String.format( 626 "Input dataset is not compatible with slice: %s cf %s", Arrays.toString(zds.shape), 627 Arrays.toString(siter.getShape()))); 628 } 629 630 IndexIterator oiter = zds.getIterator(); 631 float[] odata = zds.data; 632 633 while (siter.hasNext() && oiter.hasNext()) { 634 data[siter.index] = odata[oiter.index]; 635 data[siter.index+1] = odata[oiter.index+1]; 636 } 637 } else if (o instanceof ComplexDoubleDataset) { // IGNORE_CLASS 638 ComplexDoubleDataset zds = (ComplexDoubleDataset) o; // IGNORE_CLASS 639 640 if (!ShapeUtils.areShapesCompatible(siter.getShape(), zds.shape)) { 641 throw new IllegalArgumentException(String.format( 642 "Input dataset is not compatible with slice: %s cf %s", Arrays.toString(zds.shape), 643 Arrays.toString(siter.getShape()))); 644 } 645 646 IndexIterator oiter = zds.getIterator(); 647 double[] odata = zds.data; 648 649 while (siter.hasNext() && oiter.hasNext()) { 650 data[siter.index] = odata[oiter.index]; // PRIM_TYPE // ADD_CAST 651 data[siter.index+1] = odata[oiter.index+1]; // PRIM_TYPE // ADD_CAST 652 } 653 } else if (o instanceof IDataset) { 654 super.setSlice(o, siter); 655 } else { 656 try { 657 double vr = DTypeUtils.toReal(o); // PRIM_TYPE // ADD_CAST 658 double vi = DTypeUtils.toImag(o); // PRIM_TYPE // ADD_CAST 659 660 while (siter.hasNext()) { 661 data[siter.index] = vr; 662 data[siter.index + 1] = vi; 663 } 664 } catch (IllegalArgumentException e) { 665 throw new IllegalArgumentException("Object for setting slice is not a dataset or number"); 666 } 667 } 668 return this; 669 } 670 671 @Override 672 public ComplexDoubleDataset iadd(final Object b) { 673 setDirty(); 674 Dataset bds = b instanceof Dataset ? (Dataset) b : DatasetFactory.createFromObject(b); 675 boolean useLong = bds.getElementClass().equals(Long.class); 676 if (bds.getSize() == 1) { 677 final IndexIterator it = getIterator(); 678 final int bOffset = bds.getOffset(); 679 if (useLong) { // note no complex longs 680 final long lb = bds.getElementLongAbs(bOffset); 681 while (it.hasNext()) { 682 data[it.index] += lb; 683 } 684 } else { 685 final double db = bds.getElementDoubleAbs(bOffset); 686 if (!bds.isComplex() || bds.getElementDoubleAbs(bOffset + 1) == 0) { 687 while (it.hasNext()) { 688 data[it.index] += db; 689 } 690 } else { 691 final double vi = bds.getElementDoubleAbs(bOffset + 1); 692 while (it.hasNext()) { 693 data[it.index] += db; 694 data[it.index + 1] += vi; 695 } 696 } 697 } 698 } else { 699 final BroadcastSelfIterator it = BroadcastSelfIterator.createIterator(this, bds); 700 it.setOutputDouble(!useLong); 701 if (useLong) { // note no complex longs 702 while (it.hasNext()) { 703 data[it.aIndex] += it.bLong; 704 } 705 } else { 706 if (bds.isComplex()) { 707 while (it.hasNext()) { 708 data[it.aIndex] += it.bDouble; 709 data[it.aIndex + 1] += bds.getElementDoubleAbs(it.bIndex + 1); 710 } 711 } else { 712 while (it.hasNext()) { 713 data[it.aIndex] += it.bDouble; 714 } 715 } 716 } 717 } 718 return this; 719 } 720 721 @Override 722 public ComplexDoubleDataset isubtract(final Object b) { 723 setDirty(); 724 Dataset bds = b instanceof Dataset ? (Dataset) b : DatasetFactory.createFromObject(b); 725 boolean useLong = bds.getElementClass().equals(Long.class); 726 if (bds.getSize() == 1) { 727 final IndexIterator it = getIterator(); 728 final int bOffset = bds.getOffset(); 729 if (useLong) { // note no complex longs 730 final long lb = bds.getElementLongAbs(bOffset); 731 while (it.hasNext()) { 732 data[it.index] -= lb; 733 } 734 } else { 735 final double db = bds.getElementDoubleAbs(bOffset); 736 if (!bds.isComplex() || bds.getElementDoubleAbs(bOffset + 1) == 0) { 737 while (it.hasNext()) { 738 data[it.index] -= db; 739 } 740 } else { 741 final double vi = bds.getElementDoubleAbs(bOffset + 1); 742 while (it.hasNext()) { 743 data[it.index] -= db; 744 data[it.index + 1] -= vi; 745 } 746 } 747 } 748 } else { 749 final BroadcastSelfIterator it = BroadcastSelfIterator.createIterator(this, bds); 750 it.setOutputDouble(!useLong); 751 if (useLong) { // note no complex longs 752 while (it.hasNext()) { 753 data[it.aIndex] -= it.bLong; 754 } 755 } else { 756 if (bds.isComplex()) { 757 while (it.hasNext()) { 758 data[it.aIndex] -= it.bDouble; 759 data[it.aIndex + 1] -= bds.getElementDoubleAbs(it.bIndex + 1); 760 } 761 } else { 762 while (it.hasNext()) { 763 data[it.aIndex] -= it.bDouble; 764 } 765 } 766 } 767 } 768 return this; 769 } 770 771 @Override 772 public ComplexDoubleDataset imultiply(final Object b) { 773 setDirty(); 774 Dataset bds = b instanceof Dataset ? (Dataset) b : DatasetFactory.createFromObject(b); 775 boolean useLong = bds.getElementClass().equals(Long.class); 776 if (bds.getSize() == 1) { 777 final IndexIterator it = getIterator(); 778 final int bOffset = bds.getOffset(); 779 if (useLong) { // note no complex longs 780 final long r2 = bds.getElementLongAbs(bOffset); 781 while (it.hasNext()) { 782 data[it.index] *= r2; 783 data[it.index + 1] *= r2; 784 } 785 } else { 786 final double r2 = bds.getElementDoubleAbs(bOffset); 787 if (!bds.isComplex() || bds.getElementDoubleAbs(bOffset + 1) == 0) { 788 while (it.hasNext()) { 789 data[it.index] *= r2; 790 data[it.index + 1] *= r2; 791 } 792 } else { 793 final double i2 = bds.getElementDoubleAbs(bOffset + 1); 794 while (it.hasNext()) { 795 double r1 = data[it.index]; 796 double i1 = data[it.index + 1]; 797 data[it.index] = (r1*r2 - i1*i2); // ADD_CAST 798 data[it.index + 1] = (r1*i2 + i1*r2); // ADD_CAST 799 } 800 } 801 } 802 } else { 803 final BroadcastIterator it = BroadcastIterator.createIterator(this, bds); 804 it.setOutputDouble(!useLong); 805 if (useLong) { // note no complex longs 806 while (it.hasNext()) { 807 data[it.aIndex] *= it.bDouble; 808 data[it.aIndex + 1] *= it.bDouble; 809 } 810 } else { 811 if (bds.isComplex()) { 812 while (it.hasNext()) { 813 double r1 = it.aDouble; 814 double r2 = it.bDouble; 815 double i1 = data[it.aIndex + 1]; 816 double i2 = bds.getElementDoubleAbs(it.bIndex + 1); 817 data[it.aIndex] = (r1*r2 - i1*i2); // ADD_CAST 818 data[it.aIndex + 1] = (r1*i2 + i1*r2); // ADD_CAST 819 } 820 } else { 821 while (it.hasNext()) { 822 data[it.aIndex] *= it.bDouble; 823 data[it.aIndex + 1] *= it.bDouble; 824 } 825 } 826 } 827 } 828 return this; 829 } 830 831 @Override 832 public ComplexDoubleDataset idivide(final Object b) { 833 setDirty(); 834 Dataset bds = b instanceof Dataset ? (Dataset) b : DatasetFactory.createFromObject(b); 835 boolean useLong = bds.getElementClass().equals(Long.class); 836 if (bds.getSize() == 1) { 837 final IndexIterator it = getIterator(); 838 final int bOffset = bds.getOffset(); 839 if (useLong) { // note no complex longs 840 final long r2 = bds.getElementLongAbs(bOffset); 841 while (it.hasNext()) { 842 data[it.index] /= r2; 843 data[it.index + 1] /= r2; 844 } 845 } else { 846 final double r2 = bds.getElementDoubleAbs(bOffset); 847 if (!bds.isComplex() || bds.getElementDoubleAbs(bOffset + 1) == 0) { 848 while (it.hasNext()) { 849 data[it.index] /= r2; 850 data[it.index + 1] /= r2; 851 } 852 } else { 853 final double i2 = bds.getElementDoubleAbs(bOffset + 1); 854 if (Math.abs(r2) < Math.abs(i2)) { 855 double q = r2/i2; 856 double den = r2*q + i2; 857 while (it.hasNext()) { 858 double r1 = data[it.index]; 859 double i1 = data[it.index + 1]; 860 data[it.index] = ((r1*q + i1) / den); // ADD_CAST 861 data[it.index + 1] = ((i1*q - r1) / den); // ADD_CAST 862 } 863 } else { 864 double q = i2/r2; 865 double den = i2*q + r2; 866 if (den == 0) { 867 while (it.hasNext()) { 868 data[it.index] = Double.NaN; // CLASS_TYPE 869 data[it.index + 1] = Double.NaN; // CLASS_TYPE 870 } 871 } else { 872 while (it.hasNext()) { 873 double r1 = data[it.index]; 874 double i1 = data[it.index + 1]; 875 data[it.index] = ((i1 * q + r1) / den); // ADD_CAST 876 data[it.index + 1] = ((i1 - r1 * q) / den); // ADD_CAST 877 } 878 } 879 } 880 } 881 } 882 } else { 883 final BroadcastIterator it = BroadcastIterator.createIterator(this, bds); 884 it.setOutputDouble(!useLong); 885 if (useLong) { 886 while (it.hasNext()) { 887 data[it.aIndex] /= it.bLong; 888 data[it.aIndex + 1] /= it.bLong; 889 } 890 } else { 891 if (bds.isComplex()) { 892 while (it.hasNext()) { 893 double r1 = it.aDouble; 894 double r2 = it.bDouble; 895 double i1 = data[it.aIndex + 1]; 896 double i2 = bds.getElementDoubleAbs(it.bIndex + 1); 897 if (Math.abs(r2) < Math.abs(i2)) { 898 double q = r2/i2; 899 double den = r2*q + i2; 900 data[it.aIndex] = ((r1*q + i1) / den); // ADD_CAST 901 data[it.aIndex + 1] = ((i1*q - r1) / den); // ADD_CAST 902 } else { 903 double q = i2/r2; 904 double den = i2*q + r2; 905 if (den == 0) { 906 data[it.aIndex] = Double.NaN; // CLASS_TYPE 907 data[it.aIndex + 1] = Double.NaN; // CLASS_TYPE 908 } else { 909 data[it.aIndex] = ((i1 * q + r1) / den); // ADD_CAST 910 data[it.aIndex + 1] = ((i1 - r1 * q) / den); // ADD_CAST 911 } 912 } 913 } 914 } else { 915 while (it.hasNext()) { 916 data[it.aIndex] /= it.bDouble; 917 data[it.aIndex + 1] /= it.bDouble; 918 } 919 } 920 } 921 } 922 return this; 923 } 924 925 @Override 926 public ComplexDoubleDataset iremainder(final Object b) { 927 throw new UnsupportedOperationException("Unsupported method for class"); 928 } 929 930 @Override 931 public ComplexDoubleDataset ipower(final Object b) { 932 setDirty(); 933 Dataset bds = b instanceof Dataset ? (Dataset) b : DatasetFactory.createFromObject(b); 934 if (bds.getSize() == 1) { 935 final IndexIterator it = getIterator(); 936 final int bOffset = bds.getOffset(); 937 final double r2 = bds.getElementDoubleAbs(bOffset); 938 if (!bds.isComplex() || bds.getElementDoubleAbs(bOffset + 1) == 0) { 939 while (it.hasNext()) { 940 final Complex zd = new Complex(data[it.index], data[it.index + 1]).pow(r2); 941 data[it.index] = zd.getReal(); // ADD_CAST 942 data[it.index + 1] = zd.getImaginary(); // ADD_CAST 943 } 944 } else { 945 final Complex zv = new Complex(r2, bds.getElementDoubleAbs(bOffset + 1)); 946 while (it.hasNext()) { 947 final Complex zd = new Complex(data[it.index], data[it.index + 1]).pow(zv); 948 data[it.index] = zd.getReal(); // ADD_CAST 949 data[it.index + 1] = zd.getImaginary(); // ADD_CAST 950 } 951 } 952 } else { 953 final BroadcastIterator it = BroadcastIterator.createIterator(this, bds); 954 it.setOutputDouble(true); 955 if (bds.isComplex()) { 956 while (it.hasNext()) { 957 final Complex zv = new Complex(it.bDouble, bds.getElementDoubleAbs(it.bIndex + 1)); 958 final Complex zd = new Complex(it.aDouble, data[it.aIndex + 1]).pow(zv); 959 data[it.aIndex] = zd.getReal(); // ADD_CAST 960 data[it.aIndex + 1] = zd.getImaginary(); // ADD_CAST 961 } 962 } else { 963 while (it.hasNext()) { 964 final Complex zd = new Complex(it.aDouble, data[it.aIndex + 1]).pow(it.bDouble); 965 data[it.aIndex] = zd.getReal(); // ADD_CAST 966 data[it.aIndex + 1] = zd.getImaginary(); // ADD_CAST 967 } 968 } 969 } 970 return this; 971 } 972 973 @Override 974 public double residual(final Object b, Dataset w, boolean ignoreNaNs) { 975 Dataset bds = b instanceof Dataset ? (Dataset) b : DatasetFactory.createFromObject(b); 976 final BroadcastIterator it = BroadcastIterator.createIterator(this, bds); 977 it.setOutputDouble(true); 978 double sum = 0; 979 double comp = 0; 980 final int bis = bds.getElementsPerItem(); 981 982 if (bis == 1) { 983 if (w == null) { 984 while (it.hasNext()) { 985 double diffr = it.aDouble - it.bDouble; 986 double diffi = data[it.aIndex + 1]; 987 if (ignoreNaNs && (Double.isNaN(diffr) || Double.isNaN(diffi))) { 988 continue; 989 } 990 double err = diffr * diffr - comp; 991 double temp = sum + err; 992 comp = (temp - sum) - err; 993 sum = temp; 994 995 err = diffi * diffi - comp; 996 temp = sum + err; 997 comp = (temp - sum) - err; 998 sum = temp; 999 } 1000 } else { 1001 IndexIterator itw = w.getIterator(); 1002 while (it.hasNext() && itw.hasNext()) { 1003 final double dw = w.getElementDoubleAbs(itw.index); 1004 double diffr = it.aDouble - it.bDouble; 1005 double diffi = data[it.aIndex + 1]; 1006 if (ignoreNaNs && (Double.isNaN(diffr) || Double.isNaN(diffi))) { 1007 continue; 1008 } 1009 double err = diffr * diffr * dw - comp; 1010 double temp = sum + err; 1011 comp = (temp - sum) - err; 1012 sum = temp; 1013 1014 err = diffi * diffi * dw - comp; 1015 temp = sum + err; 1016 comp = (temp - sum) - err; 1017 sum = temp; 1018 } 1019 } 1020 } else { 1021 if (w == null) { 1022 while (it.hasNext()) { 1023 double diffr = it.aDouble - it.bDouble; 1024 double diffi = data[it.aIndex] - bds.getElementDoubleAbs(it.bIndex + 1); 1025 if (ignoreNaNs && (Double.isNaN(diffr) || Double.isNaN(diffi))) { 1026 continue; 1027 } 1028 double err = diffr * diffr - comp; 1029 double temp = sum + err; 1030 comp = (temp - sum) - err; 1031 sum = temp; 1032 1033 err = diffi * diffi - comp; 1034 temp = sum + err; 1035 comp = (temp - sum) - err; 1036 sum = temp; 1037 } 1038 } else { 1039 IndexIterator itw = w.getIterator(); 1040 while (it.hasNext() && itw.hasNext()) { 1041 final double dw = w.getElementDoubleAbs(itw.index); 1042 double diffr = it.aDouble - it.bDouble; 1043 double diffi = data[it.aIndex] - bds.getElementDoubleAbs(it.bIndex + 1); 1044 if (ignoreNaNs && (Double.isNaN(diffr) || Double.isNaN(diffi))) { 1045 continue; 1046 } 1047 double err = diffr * diffr * dw - comp; 1048 double temp = sum + err; 1049 comp = (temp - sum) - err; 1050 sum = temp; 1051 1052 err = diffi * diffi * dw - comp; 1053 temp = sum + err; 1054 comp = (temp - sum) - err; 1055 sum = temp; 1056 } 1057 } 1058 } 1059 return sum; 1060 } 1061}