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