#include <stdio.h>
#include <assert.h>

#include "pslib.h"
#include "psImageFlip.h"
#include "psRegionIsBad.h"

#include "pmFPA.h"
#include "pmFPAUtils.h"
#include "pmHDU.h"
#include "pmHDUUtils.h"
#include "pmFPACopy.h"

#define MAX(x,y) ((x) > (y) ? (x) : (y))


//////////////////////////////////////////////////////////////////////////////////////////////////////////////
// File-static functions
//////////////////////////////////////////////////////////////////////////////////////////////////////////////

// Add cells in a chip to a list
static bool addCellsFromChip(psList *list, // List of cells
                             const pmChip *chip // The chip from which to add cells
                            )
{
    assert(list);
    assert(chip);

    psArray *cells = chip->cells;       // Array of cells
    bool result = true;                 // Result of adding cells
    for (int i = 0; i < cells->n; i++) {
        pmCell *cell = cells->data[i];  // A cell
        result |= psListAdd(list, PS_LIST_TAIL, cell);
    }

    return result;
}

// Add cells in an FPA to a list
static bool addCellsFromFPA(psList *list, // List of cells
                            const pmFPA *fpa // The FPA from which to add cells
                           )
{
    assert(list);
    assert(fpa);

    psArray *chips = fpa->chips;        // Array of chips
    bool result = true;                 // Result of adding cells
    for (int i = 0; i < chips->n; i++) {
        pmChip *chip = chips->data[i];  // A chip
        result |= addCellsFromChip(list, chip);
    }

    return result;
}

// Get a list of cells that share the HDU for the target cell
static bool cellList(psList *targets,   // The list of target cells
                     psList *sources,   // The list of source cells
                     pmCell *targetCell, // The target cell
                     pmCell *sourceCell // The source cell
                    )
{
    assert(targetCell);
    assert(sourceCell);

    if (targetCell->hdu) {
        if (targets) {
            psListAdd(targets, PS_LIST_TAIL, targetCell);
        }
        if (sources) {
            psListAdd(sources, PS_LIST_TAIL, sourceCell);
        }
    } else {
        pmChip *targetChip = targetCell->parent; // The target parent chip
        pmChip *sourceChip = sourceCell->parent; // The source parent chip
        if (targetChip->hdu) {
            if (targets) {
                addCellsFromChip(targets, targetChip);
            }
            if (sources) {
                addCellsFromChip(sources, sourceChip);
            }
        } else {
            pmFPA *targetFPA = targetChip->parent; // The target parent FPA
            pmFPA *sourceFPA = sourceChip->parent; // The source parent FPA
            if (targetFPA->hdu) {
                if (targets) {
                    addCellsFromFPA(targets, targetFPA);
                }
                if (sources) {
                    addCellsFromFPA(sources, sourceFPA);
                }
            } else {
                psError(PS_ERR_IO, true, "Unable to find HDU for cell to generate list!\n");
                return false;
            }
        }
    }

    return true;
}

// Get the maximum extent of the HDU from the trimsec and biassecs
static bool sizeHDU(int *xSize, int *ySize, // Size of HDU
                    psList *cells       // List of cells
                   )
{
    psListIterator *cellsIter = psListIteratorAlloc(cells, PS_LIST_HEAD, false); // Iterator for cells
    pmCell *cell = NULL;                // The cell from iteration
    bool mdok = true;                   // Status of MD lookup
    *xSize = 0;
    *ySize = 0;
    while ((cell = psListGetAndIncrement(cellsIter))) {
        psRegion *trimsec = psMetadataLookupPtr(&mdok, cell->concepts, "CELL.TRIMSEC"); // Trim section
        if (mdok && trimsec && !psRegionIsBad(*trimsec)) {
            *xSize = MAX(trimsec->x1, *xSize);
            *ySize = MAX(trimsec->y1, *ySize);
        } else {
            psFree(cellsIter);
            return false;
        }
        psList *biassecs = psMetadataLookupPtr(&mdok, cell->concepts, "CELL.BIASSEC"); // Bias sections
        if (mdok && biassecs) {
            psListIterator *biassecsIter = psListIteratorAlloc(biassecs, PS_LIST_HEAD, false); // Iterator
            psRegion *biassec = NULL;   // The bias section
            while ((biassec = psListGetAndIncrement(biassecsIter))) {
                if (!psRegionIsBad(*trimsec)) {
                    *xSize = MAX(biassec->x1, *xSize);
                    *ySize = MAX(biassec->y1, *ySize);
                } else {
                    psFree(biassecsIter);
                    psFree(cellsIter);
                    return false;
                }
            }
            psFree(biassecsIter);
        }
    }
    psFree(cellsIter);

    return (*xSize != 0 && *ySize != 0);
}


static psRegion sectionForImage(int *position, // Position on the output image, updated
                                const psImage *image, // Image containing the sizes and offsets
                                int readdir // Read direction, 1=rows, 2=cols
                               )
{
    psRegion region;
    switch (readdir) {
    case 1:                           // Read direction is rows
        region = psRegionSet(*position, *position + image->numCols, image->row0,
                             image->row0 + image->numRows);
        *position += image->numCols;
        break;
    case 2:                           // Read direction is columns
        region = psRegionSet(image->col0, image->col0 + image->numCols, *position,
                             *position + image->numRows);
        *position += image->numRows;
        break;
    default:
        psAbort(__func__, "Shouldn't ever get here!\n");
    }

    return region;
}

static bool doBiasSections(int *position, // Position on the output image, updated
                           pmCell *target, // Target cell
                           const pmCell *source // Source cell
                          )
{
    psMetadataItem *biassecItem = psMetadataLookup(target->concepts, "CELL.BIASSEC"); // Bias sections
    if (!biassecItem) {
        psLogMsg(__func__, PS_LOG_WARN, "CELL.BIASSEC has not been initialised in target cell --- "
                 "ignored.\n");
        return false;
    }
    psFree(biassecItem->data.V);        // Blow away the old list
    psList *biassecs = psListAlloc(NULL);
    biassecItem->data.V = biassecs;

    bool mdok = true;                   // Status of MD lookup
    int readdir = psMetadataLookupS32(&mdok, source->concepts, "CELL.READDIR"); // Read direction
    if (!mdok || (readdir != 1 && readdir != 2)) {
        // Probably unnecessary, but just in case....
        psLogMsg(__func__, PS_LOG_WARN, "CELL.READDIR is not set in source cell --- ignored.\n");
        return false;
    }

    pmReadout *readout = source->readouts->data[0]; // The first source readout, as representative
    psList *biases = readout->bias; // The bias images from the source readout

    psListIterator *biasIter = psListIteratorAlloc(biases, PS_LIST_HEAD, true); // Iterator for biases
    psImage *bias = NULL;       // Bias image from iteration
    while ((bias = psListGetAndIncrement(biasIter))) {
        // Construct a region
        psRegion *biassec = psAlloc(sizeof(psRegion)); // The new region; need a psMemBlock
        *biassec = sectionForImage(position, bias, readdir);
        psListAdd(biassecs, PS_LIST_TAIL, biassec);
        psFree(biassec);        // Drop reference
    }
    psFree(biasIter);

    return true;
}

// Generate CELL.TRIMSEC and CELL.BIASSEC for the target cells
static bool generateTrimBias(psList *targets, // List of target cells
                             psList *sources // List of source cells
                            )
{
    pmCell *target = NULL, *source = NULL; // Cells from iteration
    int numCells = targets->n;          // Number of cells
    int cellNum = 0;                    // The cell number
    int position = 0;                   // Position on the image
    bool mdok = true;                   // Status of MD lookup

    // First run through to do the LHS biases
    psListIterator *targetsIter = psListIteratorAlloc(targets, PS_LIST_HEAD, false); // Iterator for targets
    psListIterator *sourcesIter = psListIteratorAlloc(sources, PS_LIST_HEAD, false); // Iterator for sources
    bool done = false;                   // Done with iteration?
    while ((target = psListGetAndIncrement(targetsIter)) && (source = psListGetAndIncrement(sourcesIter)) &&
            !done) {
        if (cellNum <= numCells/2 - 1) {
            doBiasSections(&position, target, source);
            cellNum++;
        } else {
            done = true;
        }
    }

    // Second run through to do the trim sections
    psListIteratorSet(targetsIter, PS_LIST_HEAD);
    psListIteratorSet(sourcesIter, PS_LIST_HEAD);
    while ((target = psListGetAndIncrement(targetsIter)) && (source = psListGetAndIncrement(sourcesIter))) {
        psRegion *trimsec = psMetadataLookupPtr(&mdok, target->concepts, "CELL.TRIMSEC"); // Trim section
        if (!mdok || !trimsec) {
            psLogMsg(__func__, PS_LOG_WARN, "CELL.TRIMSEC has not been initialised in target cell --- "
                     "ignored.\n");
            continue;
        }

        int readdir = psMetadataLookupS32(&mdok, source->concepts, "CELL.READDIR"); // Read direction
        if (!mdok || (readdir != 1 && readdir != 2)) {
            // Probably unnecessary, but just in case....
            psLogMsg(__func__, PS_LOG_WARN, "CELL.READDIR is not set in source cell --- ignored.\n");
            continue;
        }

        pmReadout *readout = source->readouts->data[0]; // The first source readout, as representative
        psImage *image = readout->image;// The proper image
        *trimsec = sectionForImage(&position, image, readdir);
    }

    // A final run through to do the RHS biases
    psListIteratorSet(targetsIter, cellNum);
    psListIteratorSet(sourcesIter, cellNum);
    while ((target = psListGetAndIncrement(targetsIter)) && (source = psListGetAndIncrement(sourcesIter))) {
        doBiasSections(&position, target, source);
    }

    // Clean up
    psFree(targetsIter);
    psFree(sourcesIter);

    return (position > 0);
}


// Generate an HDU with the pixels
static bool generateHDU(pmCell *target,  // The target cell
                        pmCell *source, // The source cell
                        int xBin, int yBin // Binning in x and y
                       )
{
    // Get the HDU and a list of cells below it
    pmHDU *hdu = pmHDUFromCell(target); // The HDU in the target cell
    psList *targetCells = psListAlloc(NULL); // List of target cells below the target HDU
    psList *sourceCells = psListAlloc(NULL); // List of source cells below the target HDU
    if (! cellList(targetCells, sourceCells, target, source)) {
        psError(PS_ERR_IO, true, "Unable to find cells to generate HDU!\n");
        return false;
    }

    // Check the number of readouts
    int numReadouts = -1;               // Number of readouts
    {
        psListIterator *iter = psListIteratorAlloc(sourceCells, PS_LIST_HEAD, false); // Iterator for cells
        pmCell *cell = NULL;                // The cell from iteration
        while ((cell = psListGetAndIncrement(iter)))
        {
            psArray *readouts = cell->readouts;
            if (numReadouts == -1) {
                numReadouts = readouts->n;
            } else if (readouts->n != numReadouts) {
                psError(PS_ERR_IO, true, "Number of readouts doesn't match: %d vs %d\n", readouts->n,
                        numReadouts);
                return false;
            }

        }
        psFree(iter);
    }

    // Get the size of the HDU, either from existing trimsec and biassec, or generate these and try again
    int xSize = 0, ySize = 0;           // Size of HDU
    if (!sizeHDU(&xSize, &ySize, targetCells) && !(generateTrimBias(targetCells, sourceCells) &&
            sizeHDU(&xSize, &ySize, targetCells))) {
        psError(PS_ERR_IO, true, "Unable to determine size of HDU!\n");
        return false;
    }
    psFree(targetCells);
    psFree(sourceCells);

    xSize = (int)ceilf((float)xSize/(float)xBin);
    ySize = (int)ceilf((float)ySize/(float)yBin);

    hdu->images = psArrayAlloc(numReadouts);
    for (int i = 0; i < numReadouts; i++) {
        psImage *image = psImageAlloc(xSize, ySize, PS_TYPE_F32);
        psImageInit(image, 0.0);
        hdu->images->data[i] = image;
    }

    return true;
}

// Copy pixels from a target image to a source image, with flips
static bool copyPixels(psImage *target, // Target image (HDU pixels)
                       psImage *source, // Source image (from source cell)
                       psRegion region, // Region for pasting
                       bool xFlip,      // Flip in x?
                       bool yFlip       // Flip in y?
                      )
{
    psImage *overlay = psMemIncrRefCounter(source);
    if (xFlip) {
        psImage *temp = psImageFlipX(overlay);
        psFree(overlay);
        overlay = temp;
    }
    if (yFlip) {
        psImage *temp = psImageFlipY(overlay);
        psFree(overlay);
        overlay = temp;
    }
    int numPix = psImageOverlaySection(target, overlay, region.x0, region.y0, "=");
    psFree(overlay);
    return (numPix > 0);
}

// Bin a region down by specified factors in x and y
static void binRegion(psRegion *region, // Region to be binned
                      int xBin, int yBin// Binning in x and y
                     )
{
    // Want to include the lower bound: 1 binned by 4 --> 0; 3 binned by 4 --> 0; 4 binned by 4 --> 1
    region->x0 = (int)(region->x0 / xBin);
    region->y0 = (int)(region->y0 / yBin);
    // Want to exclude the upper bound: 4 binned by 4 --> 1; 5 binned by 4 --> 2; 7 binned by 4 --> 2
    region->x1 = (int)((region->x1 + xBin - 1) / xBin);
    region->y1 = (int)((region->y1 + yBin - 1) / yBin);
}


static pmHDU *findPHU(const pmCell *cell// The cell for which to find the PHU
                     )
{
    if (cell->hdu && cell->hdu->phu) {
        return cell->hdu;
    }
    pmChip *chip = cell->parent;        // The parent chip
    if (chip->hdu && chip->hdu->phu) {
        return chip->hdu;
    }
    pmFPA *fpa = chip->parent;  // The parent FPA
    if (fpa->hdu && fpa->hdu->phu) {
        return fpa->hdu;
    }

    psError(PS_ERR_IO, true, "Unable to find PHU for target cell.\n");
    return NULL;
}


static bool copyPHU(pmCell *targetCell, // The target cell
                    pmCell *sourceCell  // The source cell
                   )
{
    // Find the respective PHUs
    pmHDU *targetPHU = findPHU(targetCell); // The target PHU
    if (targetPHU->header) {
        return false;                   // No work required
    }
    // Copy the header over
    pmHDU *sourcePHU = findPHU(sourceCell); // The source PHU
    targetPHU->header = psMetadataCopy(NULL, sourcePHU->header);
    return true;
}


//////////////////////////////////////////////////////////////////////////////////////////////////////////////
// File-static engine functions --- these do all the work.  Actually, cellCopy does all the work; the others
// merely iterate on the higher-level components.
//////////////////////////////////////////////////////////////////////////////////////////////////////////////

static int cellCopy(pmCell *target,     // The target cell
                    pmCell *source,     // The source cell, to be copied
                    bool pixels,        // Copy the pixels?
                    int xBin, int yBin  // (Relative) binning factors in x and y
                   )
{
    assert(target);
    assert(source);

    psArray *sourceReadouts = source->readouts; // The source readouts
    int numReadouts = sourceReadouts->n; // Number of readouts copied

    // Copy any headers
    if (target->hdu && !target->hdu->phu) {
        pmHDU *sourceHDU = pmHDUFromCell(source);
        target->hdu->header = psMetadataCopy(target->hdu->header, sourceHDU->header);
    }

    pmHDU *hdu = pmHDUFromCell(target); // The target HDU; we need to fix this up
    if (!hdu->images) {
        generateHDU(target, source, xBin, yBin);
    }
    if (!hdu->header) {
        hdu->header = psMetadataAlloc();
    }
    // Copy the PHU over as well, if required
    copyPHU(target, source);

    // Need to check/change CELL.XPARITY and CELL.YPARITY
    bool xFlip = (psMetadataLookupS32(NULL, target->concepts, "CELL.XPARITY") !=
                  psMetadataLookupS32(NULL, source->concepts, "CELL.XPARITY")); // Switch parity in x?
    bool yFlip = (psMetadataLookupS32(NULL, target->concepts, "CELL.YPARITY") !=
                  psMetadataLookupS32(NULL, source->concepts, "CELL.YPARITY")); // Switch parity in y?
    psTrace(__func__, 3, "xFlip: %d; yFlip: %d\n", xFlip, yFlip);

    bool mdok = true;                   // Status of MD lookup
    psRegion *trimsec = psMetadataLookupPtr(&mdok, target->concepts, "CELL.TRIMSEC"); // The trim section
    if (!mdok || !trimsec || psRegionIsBad(*trimsec)) {
        psError(PS_ERR_IO, true, "CELL.TRIMSEC isn't set!\n");
        return 0;
    }
    psList *biassecs = psMetadataLookupPtr(&mdok, target->concepts, "CELL.BIASSEC"); // The bias sections
    if (!mdok || !biassecs) {
        psError(PS_ERR_IO, true, "CELL.BIASSEC isn't set!\n");
        return 0;
    }

    for (int i = 0; i < numReadouts; i++) {
        pmReadout *sourceReadout = sourceReadouts->data[i]; // The source readout
        psImage *sourceImage = sourceReadout->image; // The source image
        pmReadout *targetReadout = pmReadoutAlloc(target); // The target readout; this adds it to the cell
        if (sourceImage->numCols != trimsec->x1 - trimsec->x0 ||
                sourceImage->numRows != trimsec->y1 - trimsec->y0) {
            psString trimsecString = psRegionToString(*trimsec); // String with the trim section
            psLogMsg(__func__, PS_LOG_WARN, "Source image size (%dx%d) for readout %d doesn't match "
                     "CELL.TRIMSEC for target (%s) -- ignored.\n", sourceImage->numCols, sourceImage->numRows,
                     i, trimsecString);
            psFree(trimsecString);
        } else {
            binRegion(trimsec, xBin, yBin);
            if (pixels) {
                copyPixels(hdu->images->data[i], sourceImage, *trimsec, xFlip, yFlip);
            }
            targetReadout->image = psImageSubset(hdu->images->data[i], *trimsec);
        }

        psListIterator *biassecsIter = psListIteratorAlloc(biassecs, PS_LIST_HEAD, false); // Iterator
        psRegion *biassec = NULL;       // Bias section from iteration
        psListIterator *biasIter = psListIteratorAlloc(sourceReadout->bias, PS_LIST_HEAD, false); // Iterator
        psImage *bias = NULL;           // Bias image from iteration
        while ((biassec = psListGetAndIncrement(biassecsIter)) && (bias = psListGetAndIncrement(biasIter))) {
            if (psRegionIsBad(*biassec)) {
                psString biassecString = psRegionToString(*biassec); // String for bias section
                psLogMsg(__func__, PS_LOG_WARN, "Bias section (%s) isn't set --- ignored.\n", biassecString);
                psFree(biassecString);
                continue;
            }
            if (bias->numCols != biassec->x1 - biassec->x0 ||
                    bias->numRows != biassec->y1 - biassec->y0) {
                psString biassecString = psRegionToString(*biassec); // String with the bias section
                psLogMsg(__func__, PS_LOG_WARN, "Source image size (%dx%d) for readout %d doesn't match "
                         "CELL.BIASSEC for target (%s) -- ignored.\n", bias->numCols, bias->numRows, i,
                         biassecString);
                psFree(biassecString);
            } else {
                binRegion(biassec, xBin, yBin);
                if (pixels) {
                    copyPixels(hdu->images->data[i], bias, *biassec, xFlip, yFlip);
                }
                psImage *newBias = psImageSubset(hdu->images->data[i], *biassec);
                psListAdd(targetReadout->bias, PS_LIST_TAIL, newBias);
                psFree(newBias);        // Drop reference
            }
        }
        psFree(targetReadout);          // Drop reference
        psFree(biassecsIter);
        psFree(biasIter);
    }

    // Copy the remaining "concepts" over
    psMetadataIterator *conceptsIter = psMetadataIteratorAlloc(source->concepts, PS_LIST_HEAD, NULL);
    psMetadataItem *conceptItem = NULL; // Item from iteration
    while ((conceptItem = psMetadataGetAndIncrement(conceptsIter))) {
        psString name = conceptItem->name; // Name of concept
        if (strcmp(name, "CELL.TRIMSEC") != 0 && strcmp(name, "CELL.BIASSEC") != 0 &&
                strcmp(name, "CELL.XPARITY") != 0 && strcmp(name, "CELL.YPARITY") != 0 &&
                strcmp(name, "CELL.X0") != 0 && strcmp(name, "CELL.Y0") != 0) {
            psMetadataAddItem(target->concepts, conceptItem, PS_LIST_TAIL, PS_META_REPLACE);
        }
    }
    psFree(conceptsIter);

    // Need to update CELL.X0 and CELL.Y0 if we flipped
    if (xFlip) {
        int xZero = psMetadataLookupS32(NULL, source->concepts, "CELL.X0"); // CELL.X0 from source
        int xParity = psMetadataLookupS32(NULL, target->concepts, "CELL.XPARITY"); // Parity in x
        int xBin = psMetadataLookupS32(NULL, source->concepts, "CELL.XBIN"); // CELL.XBIN from source
        pmReadout *readout = source->readouts->data[0]; // A representative readout
        psTrace(__func__, 3, "CELL.X0: Before: %d After: %d\n", xZero,
                xZero - (readout->image->numCols - 1) * xParity * xBin);
        psTrace(__func__, 9, "(xParity: %d xBin: %d numCols: %d)\n", xParity, xBin, readout->image->numCols);
        xZero -= (readout->image->numCols - 1) * xParity * xBin; // Change the parity on the X0 position
        psMetadataItem *newItem = psMetadataLookup(target->concepts, "CELL.X0"); // CELL.X0 from target
        newItem->data.S32 = xZero;
    }
    if (yFlip) {
        int yZero = psMetadataLookupS32(NULL, source->concepts, "CELL.Y0"); // CELL.Y0 from source
        int yParity = psMetadataLookupS32(NULL, target->concepts, "CELL.YPARITY"); // Parity in y
        int yBin = psMetadataLookupS32(NULL, source->concepts, "CELL.YBIN"); // Parity in y
        pmReadout *readout = source->readouts->data[0]; // A representative readout
        psTrace(__func__, 3, "CELL.Y0: Before: %d After: %d\n", yZero,
                yZero - (readout->image->numRows - 1) * yParity * yBin);
        psTrace(__func__, 9, "(yParity: %d yBin: %d numRows: %d)\n", yParity, yBin, readout->image->numRows);
        yZero -= (readout->image->numRows - 1) * yParity * yBin; // Change the parity on the Y0 position
        psMetadataItem *newItem = psMetadataLookup(target->concepts, "CELL.Y0"); // CELL.Y0 from target
        newItem->data.S32 = yZero;
    }

    // Update the binning concepts
    psMetadataItem *binItem = psMetadataLookup(target->concepts, "CELL.XBIN");
    binItem->data.S32 *= xBin;
    binItem = psMetadataLookup(target->concepts, "CELL.YBIN");
    binItem->data.S32 *= yBin;

    return numReadouts;
}

static int chipCopy(pmChip *target,          // The target chip
                    pmChip *source,          // The source chip, to be copied
                    bool pixels,             // Copy the pixels?
                    int xBin, int yBin       // (Relative) binning factors in x and y
                   )
{
    assert(target);
    assert(source);

    psArray *targetCells = target->cells; // The target cells
    psArray *sourceCells = source->cells; // The source cells
    if (targetCells->n != sourceCells->n) {
        psError(PS_ERR_IO, true, "Number of source cells (%d) differs from the number of target cells (%d)\n",
                sourceCells->n, targetCells->n);
        return false;
    }

    // Copy any headers
    if (target->hdu && !target->hdu->phu) {
        pmHDU *sourceHDU = pmHDUFromChip(source);
        target->hdu->header = psMetadataCopy(target->hdu->header, sourceHDU->header);
    }

    int numCells = 0;                   // Number of cells copied
    for (int i = 0; i < targetCells->n; i++) {
        pmCell *targetCell = targetCells->data[i]; // The target cell
        const char *cellName = psMetadataLookupStr(NULL, targetCell->concepts, "CELL.NAME"); // Name of cell
        int cellNum = pmChipFindCell(source, cellName); // Number of cell with that name
        if (cellNum >= 0) {
            pmCell *sourceCell = sourceCells->data[cellNum]; // The source cell
            int numReadouts = cellCopy(targetCell, sourceCell, pixels, xBin, yBin); // Number of readouts
            // copied
            psTrace(__func__, 5, "Copied %d readouts for cell %s\n", numReadouts, cellName);
            numCells++;
        }
    }

    // Update the concepts
    psMetadataCopy(target->concepts, source->concepts);

    return numCells;

}

static int fpaCopy(pmFPA *target,            // The target FPA
                   pmFPA *source,            // The source FPA, to be copied
                   bool pixels,              // Copy the pixels?
                   int xBin, int yBin        // (Relative) binning factors in x and y
                  )
{
    assert(target);
    assert(source);

    psArray *targetChips = target->chips; // The target chips
    psArray *sourceChips = source->chips; // The source chips
    if (targetChips->n != sourceChips->n) {
        psError(PS_ERR_IO, true, "Number of source chips (%d) differs from the number of target chips (%d)\n",
                sourceChips->n, targetChips->n);
        return false;
    }

    // Copy any headers
    if (target->hdu && !target->hdu->phu) {
        pmHDU *sourceHDU = pmHDUFromFPA(source);
        target->hdu->header = psMetadataCopy(target->hdu->header, sourceHDU->header);
    }

    int numChips = 0;                   // Number of chips copied
    for (int i = 0; i < targetChips->n; i++) {
        pmChip *targetChip = targetChips->data[i]; // The target chip
        const char *chipName = psMetadataLookupStr(NULL, targetChip->concepts, "CHIP.NAME"); // Name of chip
        int chipNum = pmFPAFindChip(source, chipName); // Number of chip with that name
        if (chipNum >= 0) {
            pmChip *sourceChip = sourceChips->data[chipNum]; // The source chip
            int numCells = chipCopy(targetChip, sourceChip, pixels, xBin, yBin); // Number of cells copied
            psTrace(__func__, 5, "Copied %d cells for chip %s\n", numCells, chipName);
            numChips++;
        }
    }

    // Update the concepts
    psMetadataCopy(target->concepts, source->concepts);

    return numChips;
}

//////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Public functions
//////////////////////////////////////////////////////////////////////////////////////////////////////////////


int pmFPACopy(pmFPA *target,            // The target FPA
              pmFPA *source             // The source FPA, to be copied
             )
{
    return fpaCopy(target, source, true, 1, 1);
}

int pmChipCopy(pmChip *target,          // The target chip
               pmChip *source           // The source chip, to be copied
              )
{
    return chipCopy(target, source, true, 1, 1);
}

int pmCellCopy(pmCell *target,          // The target cell
               pmCell *source           // The source cell, to be copied
              )
{
    return cellCopy(target, source, true, 1, 1);
}


int pmFPACopyStructure(pmFPA *target,   // The target FPA
                       pmFPA *source,   // The source FPA, to be copied
                       int xBin, int yBin // Binning factors in x and y
                      )
{
    return fpaCopy(target, source, false, xBin, yBin);
}

int pmChipCopyStructure(pmChip *target, // The target chip
                        pmChip *source, // The source chip, to be copied
                        int xBin, int yBin // Binning factors in x and y
                       )
{
    return chipCopy(target, source, false, xBin, yBin);
}

int pmCellCopyStructure(pmCell *target, // The target cell
                        pmCell *source, // The source cell, to be copied
                        int xBin, int yBin // Binning factors in x and y
                       )
{
    return cellCopy(target, source, false, xBin, yBin);
}
