IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 8480


Ignore:
Timestamp:
Aug 22, 2006, 2:48:04 PM (20 years ago)
Author:
Paul Price
Message:

Adding GPL stuff, for distribution

Location:
trunk/pois
Files:
1 added
17 edited

Legend:

Unmodified
Added
Removed

  • trunk/pois/src/errorCodes-skl.c

    r3809 r8480  
    2424    for (int i = 0; i < nerror; i++) {
    2525       psErrorDescription *tmp = psAlloc(sizeof(psErrorDescription));
    26         p_psMemSetPersistent(tmp, true);
     26        p_psMemSetPersistent(tmp, true);
    2727       *tmp = errors[i];
    2828       psErrorRegister(tmp, 1);
    29        psFree(tmp);                     /* it's on the internal list */
     29       psFree(tmp);                     /* it's on the internal list */
    3030    }
    31     nerror = 0;                                 // don't register more than once
     31    nerror = 0;                                 // don't register more than once
    3232}

  • trunk/pois/src/pois.c

    r8077 r8480  
    1 /*
    2  * POIS: Pan-STARRS (or Pricey's) Optimal Image Subtraction
    3  *
    4  * Attempts to fit for a general convolution kernel that will match the PSFs between two images.
    5  *
    6  */
     1// POIS: Pan-STARRS (or Pricey's) Optimal Image Subtraction
     2//
     3// Attempts to fit for a general convolution kernel that will match the PSFs between two images.
     4//
     5// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     6//
     7// This program is free software; you can redistribute it and/or modify
     8// it under the terms of the GNU General Public License as published by
     9// the Free Software Foundation; either version 2 of the License, or
     10// (at your option) any later version.
     11//
     12// This program is distributed in the hope that it will be useful,
     13// but WITHOUT ANY WARRANTY; without even the implied warranty of
     14// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     15// GNU General Public License for more details.
     16//
     17// You should have received a copy of the GNU General Public License
     18// along with this program; if not, write to the Free Software
     19// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
    720
    821#include <stdio.h>

  • trunk/pois/src/poisAddPenalty.c

    r3792 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include <assert.h>
    … …  
    420#include "pois.h"
    521
    6 void poisAddPenalty(psImage *matrix,    // The matrix
    7                     psVector *vector,   // The vector
    8                     float penalty,      // Penalty value
    9                     const psArray *kernels, // Kernel basis functions
    10                     const poisConfig *config // Configuration
     22void poisAddPenalty(psImage *matrix,    // The matrix
     23                    psVector *vector,   // The vector
     24                    float penalty,      // Penalty value
     25                    const psArray *kernels, // Kernel basis functions
     26                    const poisConfig *config // Configuration
    1127    )
    1228{
    … …  
    2541
    2642    for (int j = 0; j < numKernels; j++) {
    27         poisKernelBasis *kernel = kernels->data[j];
    28         if (kernel->xOrder == 0 && kernel->yOrder == 0) {
    29             for (int k = j + 1; k < numKernels; k++) {
    30                 poisKernelBasis *compare = kernels->data[k];
    31                 if (compare->xOrder == 0 && compare->yOrder == 0 &&
    32                     (compare->u == kernel->u + 1 || compare->u == kernel->u - 1 ||
    33                      compare->v == kernel->v + 1 || compare->v == kernel->v - 1)) {
    34                     matrix->data.F64[j][k] -= penalty;
    35                     matrix->data.F64[k][j] -= penalty;
    36                 }
    37             }
    38         }
    39            
    40         if (abs(kernel->u) == config->xKernel) {
    41             if (abs(kernel->v) == config->yKernel) {
    42                 vector->data.F64[j] -= 3.0 * penalty;
    43             } else {
    44                 vector->data.F64[j] -= 5.0 * penalty;
    45             }
    46         } else {
    47             if (abs(kernel->v) == config->yKernel) {
    48                 vector->data.F64[j] -= 5.0 * penalty;
    49             } else {
    50                 vector->data.F64[j] -= 8.0 * penalty;
    51             }
    52         }
    53        
     43        poisKernelBasis *kernel = kernels->data[j];
     44        if (kernel->xOrder == 0 && kernel->yOrder == 0) {
     45            for (int k = j + 1; k < numKernels; k++) {
     46                poisKernelBasis *compare = kernels->data[k];
     47                if (compare->xOrder == 0 && compare->yOrder == 0 &&
     48                    (compare->u == kernel->u + 1 || compare->u == kernel->u - 1 ||
     49                     compare->v == kernel->v + 1 || compare->v == kernel->v - 1)) {
     50                    matrix->data.F64[j][k] -= penalty;
     51                    matrix->data.F64[k][j] -= penalty;
     52                }
     53            }
     54        }
     55
     56        if (abs(kernel->u) == config->xKernel) {
     57            if (abs(kernel->v) == config->yKernel) {
     58                vector->data.F64[j] -= 3.0 * penalty;
     59            } else {
     60                vector->data.F64[j] -= 5.0 * penalty;
     61            }
     62        } else {
     63            if (abs(kernel->v) == config->yKernel) {
     64                vector->data.F64[j] -= 5.0 * penalty;
     65            } else {
     66                vector->data.F64[j] -= 8.0 * penalty;
     67            }
     68        }
     69
    5470    }
    5571

  • trunk/pois/src/poisCalculateDeviations.c

    r7075 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include <assert.h>

  • trunk/pois/src/poisCalculateEquation.c

    r6891 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include <math.h>

  • trunk/pois/src/poisCheckKernel.c

    r6891 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include "pslib.h"

  • trunk/pois/src/poisConvolveImage.c

    r6452 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include <assert.h>

  • trunk/pois/src/poisExtractKernel.c

    r6452 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include <assert.h>

  • trunk/pois/src/poisFindStamps.c

    r6891 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include <assert.h>

  • trunk/pois/src/poisKernelBasisFunctions.c

    r6891 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include <assert.h>

  • trunk/pois/src/poisMakeMask.c

    r5717 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include <assert.h>
    … …  
    723#define MIN(a,b) ((a) < (b) ? (a) : (b))
    824
    9 psImage *poisMakeMask(psImage *mask,    // Mask to update, or NULL
    10                       const psImage *refImage, // Reference image
    11                       const psImage *inImage, // Input image
    12                       const poisConfig *config // Configuration
     25psImage *poisMakeMask(psImage *mask,    // Mask to update, or NULL
     26                      const psImage *refImage, // Reference image
     27                      const psImage *inImage, // Input image
     28                      const poisConfig *config // Configuration
    1329    )
    1430{
    … …  
    2036    assert(config);
    2137
    22     int nx = config->xImage;            // Size in x
    23     int ny = config->yImage;            // Size in y
     38    int nx = config->xImage;            // Size in x
     39    int ny = config->yImage;            // Size in y
    2440
    2541    psTrace("pois.makeMask", 5, "Creating mask: %dx%d\n", nx, ny);
    … …  
    2743    // Allocate the mask image, if necessary
    2844    if (mask == NULL) {
    29         mask = psImageAlloc(nx, ny, PS_TYPE_U8);
    30         // Initialise to zero
    31         for (int j = 0; j < ny; j++) {
    32             for (int i = 0; i < nx; i++) {
    33                 mask->data.U8[j][i] = POIS_MASK_OK;
    34             }
    35         }
     45        mask = psImageAlloc(nx, ny, PS_TYPE_U8);
     46        // Initialise to zero
     47        for (int j = 0; j < ny; j++) {
     48            for (int i = 0; i < nx; i++) {
     49                mask->data.U8[j][i] = POIS_MASK_OK;
     50            }
     51        }
    3652    }
    3753
    3854    // Iterate over the images
    3955    for (int y = config->yKernel; y < ny - config->yKernel; y++) {
    40         for (int x = config->xKernel; x < nx - config->xKernel; x++) {
    41             // Mask all around something that's saturated in the reference, since it gets convolved
    42             // But only mask single pixel in the input image, since it doesn't get convolved
    43             if (refImage->data.F32[y][x] >= config->refSat || refImage->data.F32[y][x] <= config->refBad ||
    44                 isnan(refImage->data.F32[y][x])) {
    45                 mask->data.U8[y][x] |= POIS_MASK_BAD;
    46                 int xlow = MAX(0, x - config->xKernel);
    47                 int xhigh = MIN(nx, x + config->xKernel);
    48                 int ylow = MAX(0, y - config->yKernel);
    49                 int yhigh = MIN(ny, y + config->yKernel);
    50                 for (int v = ylow; v <= yhigh; v++) {
    51                     for (int u = xlow; u <= xhigh; u++) {
    52                         mask->data.U8[v][u] |= POIS_MASK_NEAR_BAD;
    53                     }
    54                 }
    55             } else if (inImage->data.F32[y][x] >= config->inSat || inImage->data.F32[y][x] <= config->inBad ||
    56                 isnan(inImage->data.F32[y][x])) {
    57                 mask->data.U8[y][x] |= POIS_MASK_BAD;
    58             }
    59         }
     56        for (int x = config->xKernel; x < nx - config->xKernel; x++) {
     57            // Mask all around something that's saturated in the reference, since it gets convolved
     58            // But only mask single pixel in the input image, since it doesn't get convolved
     59            if (refImage->data.F32[y][x] >= config->refSat || refImage->data.F32[y][x] <= config->refBad ||
     60                isnan(refImage->data.F32[y][x])) {
     61                mask->data.U8[y][x] |= POIS_MASK_BAD;
     62                int xlow = MAX(0, x - config->xKernel);
     63                int xhigh = MIN(nx, x + config->xKernel);
     64                int ylow = MAX(0, y - config->yKernel);
     65                int yhigh = MIN(ny, y + config->yKernel);
     66                for (int v = ylow; v <= yhigh; v++) {
     67                    for (int u = xlow; u <= xhigh; u++) {
     68                        mask->data.U8[v][u] |= POIS_MASK_NEAR_BAD;
     69                    }
     70                }
     71            } else if (inImage->data.F32[y][x] >= config->inSat || inImage->data.F32[y][x] <= config->inBad ||
     72                isnan(inImage->data.F32[y][x])) {
     73                mask->data.U8[y][x] |= POIS_MASK_BAD;
     74            }
     75        }
    6076    }
    6177

  • trunk/pois/src/poisMaskOps.c

    r3806 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117/*
    218 * Some utility routines that should be in psLib
    … …  
    824#include "pois.h"
    925
    10 psImage *poisImageSetVal(psImage *restrict in,       // Input image
    11                          const psC64 val             // set to this value
     26psImage *poisImageSetVal(psImage *restrict in,       // Input image
     27                         const psC64 val             // set to this value
    1228    )
    1329{
    … …  
    1632    switch (in->type.type) {
    1733      case PS_TYPE_U8:
    18         {
    19             const int numCols = in->numCols;
    20             psU8 *row0 = in->data.U8[0];
    21             for (int c = 0; c < numCols; c++) {
    22                 row0[c] = val;
    23             }
    24             for (int r = 1; r < in->numRows; r++) {
    25                 memcpy(in->data.U8[r], row0, numCols*sizeof(psU8));
    26             }
    27         }
    28         break;
     34        {
     35            const int numCols = in->numCols;
     36            psU8 *row0 = in->data.U8[0];
     37            for (int c = 0; c < numCols; c++) {
     38                row0[c] = val;
     39            }
     40            for (int r = 1; r < in->numRows; r++) {
     41                memcpy(in->data.U8[r], row0, numCols*sizeof(psU8));
     42            }
     43        }
     44        break;
    2945      case PS_TYPE_F32:
    30         {
    31             const int numCols = in->numCols;
    32             psF32 *row0 = in->data.F32[0];
    33             for (int c = 0; c < numCols; c++) {
    34                 row0[c] = val;
    35             }
    36             for (int r = 1; r < in->numRows; r++) {
    37                 memcpy(in->data.F32[r], row0, numCols*sizeof(psF32));
    38             }
    39         }
    40         break;
     46        {
     47            const int numCols = in->numCols;
     48            psF32 *row0 = in->data.F32[0];
     49            for (int c = 0; c < numCols; c++) {
     50                row0[c] = val;
     51            }
     52            for (int r = 1; r < in->numRows; r++) {
     53                memcpy(in->data.F32[r], row0, numCols*sizeof(psF32));
     54            }
     55        }
     56        break;
    4157      default:
    42         psError(POIS_ERR_UNSUPPORTED_TYPE, true,
    43                 "Type %d is not supported", in->type.type);
    44         return NULL;
     58        psError(POIS_ERR_UNSUPPORTED_TYPE, true,
     59                "Type %d is not supported", in->type.type);
     60        return NULL;
    4561    }
    4662
    … …  
    4864}
    4965
    50 psImage *poisImageSetValInMask(psImage *out,                    // Image to update, or NULL
    51                                const psImage *restrict in,      // Input image
    52                                const psImage *mask,             // mask for image
    53                                const psC64 val,                 // set to this value
    54                                const unsigned long bits         // set pixels where mask & bits != 0
     66psImage *poisImageSetValInMask(psImage *out,                    // Image to update, or NULL
     67                               const psImage *restrict in,      // Input image
     68                               const psImage *mask,             // mask for image
     69                               const psC64 val,                 // set to this value
     70                               const unsigned long bits         // set pixels where mask & bits != 0
    5571    )
    5672{
    … …  
    6581    // Allocate the output image, if necessary
    6682    if (out == NULL) {
    67         out = psImageCopy(NULL, in, in->type.type);
     83        out = psImageCopy(NULL, in, in->type.type);
    6884    }
    6985
    7086    switch (in->type.type) {
    7187      case PS_TYPE_F32:
    72         {
    73             const int numCols = in->numCols;
    74             for (int r = 0; r < in->numRows; r++) {
    75                 psF32 *row = out->data.F32[r]; // n.b. a copy of in
    76                 psU8 *mrow = mask->data.U8[r];
    77                 for (int c = 0; c < numCols; c++) {
    78                     if (mrow[c] & bits) {
    79                         row[c] = val;
    80                     }
    81                 }
    82             }
    83         }
    84         break;
     88        {
     89            const int numCols = in->numCols;
     90            for (int r = 0; r < in->numRows; r++) {
     91                psF32 *row = out->data.F32[r]; // n.b. a copy of in
     92                psU8 *mrow = mask->data.U8[r];
     93                for (int c = 0; c < numCols; c++) {
     94                    if (mrow[c] & bits) {
     95                        row[c] = val;
     96                    }
     97                }
     98            }
     99        }
     100        break;
    85101      default:
    86         psError(POIS_ERR_UNSUPPORTED_TYPE, true,
    87                 "Type %d is not supported", in->type.type);
    88         return NULL;
     102        psError(POIS_ERR_UNSUPPORTED_TYPE, true,
     103                "Type %d is not supported", in->type.type);
     104        return NULL;
    89105    }
    90106

  • trunk/pois/src/poisParseConfig.c

    r5717 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include <unistd.h>
    … …  
    723{
    824    fprintf (stderr, "POIS: Pan-STARRS (or Pricey's) Optimal Image Subtraction\n"
    9              "Usage: pois [-h] [-v] [-k X Y] [-t T] [-s S1 S2] [-b B1 B2] [-B BG] [-o N] [-n NSTAMPS] [-S STAMPS] [-i ITER] [-r REJ] [-m] REF IN OUT\n"
    10              "where\n"
    11              "\t-h         Help (this info)\n"
    12              "\t-v         Increase verbosity level\n"
    13              "\t-k X Y     Kernel half-size in x and y (5)\n"
    14              "\t-t T       Threshold for stamps on reference image (0)\n"
    15              "\t-s S1 S2   Saturation level for reference (S1) and input (S2) (50000)\n"
    16              "\t-b B1 B2   Bad level for both images (0)\n"
    17              "\t-B BG      Level to add to background to get above zero (0)\n"
    18              "\t-o ORDER   Order of spatial variations in the kernel (0)\n"
    19              "\t-n NSTAMPS Number of stamps in each dimension (5)\n"
    20              "\t-S STAMPS  File from which to read stamps\n"
    21              "\t-R         Reverse the sense of the subtraction (REF - IN)\n"
    22              "\t-i ITER    Number of rejection iterations (10)\n"
    23              "\t-r REJ     Rejection level in standard deviations (2.5)\n"
    24              "\t-m         Output mask frame (OUT.mask)\n"
    25              "\tREF        Reference image\n"
    26              "\tIN         Input image (to be matched to REF)\n"
    27              "\tOUT        Output image\n"
    28         );
     25             "Usage: pois [-h] [-v] [-k X Y] [-t T] [-s S1 S2] [-b B1 B2] [-B BG] [-o N] [-n NSTAMPS] [-S STAMPS] [-i ITER] [-r REJ] [-m] REF IN OUT\n"
     26             "where\n"
     27             "\t-h         Help (this info)\n"
     28             "\t-v         Increase verbosity level\n"
     29             "\t-k X Y     Kernel half-size in x and y (5)\n"
     30             "\t-t T       Threshold for stamps on reference image (0)\n"
     31             "\t-s S1 S2   Saturation level for reference (S1) and input (S2) (50000)\n"
     32             "\t-b B1 B2   Bad level for both images (0)\n"
     33             "\t-B BG      Level to add to background to get above zero (0)\n"
     34             "\t-o ORDER   Order of spatial variations in the kernel (0)\n"
     35             "\t-n NSTAMPS Number of stamps in each dimension (5)\n"
     36             "\t-S STAMPS  File from which to read stamps\n"
     37             "\t-R         Reverse the sense of the subtraction (REF - IN)\n"
     38             "\t-i ITER    Number of rejection iterations (10)\n"
     39             "\t-r REJ     Rejection level in standard deviations (2.5)\n"
     40             "\t-m         Output mask frame (OUT.mask)\n"
     41             "\tREF        Reference image\n"
     42             "\tIN         Input image (to be matched to REF)\n"
     43             "\tOUT        Output image\n"
     44        );
    2945}
    3046
    3147poisConfig *poisParseConfig(int argc, // Number of command-line arguments
    32                             char **argv // Command-line arguments
     48                            char **argv // Command-line arguments
    3349    )
    3450{
    35     poisConfig *config;                 // Configuration values
     51    poisConfig *config;                 // Configuration values
    3652
    3753    /* Variables for getopt */
    … …  
    6884    while ((opt = getopt(argc, argv, "hvRmk:f:q:t:s:o:b:B:n:i:r:p:S:")) != -1) {
    6985        switch (opt) {
    70           case 'h':
    71             help();
    72             exit(EXIT_SUCCESS);
    73           case 'v':
     86          case 'h':
     87            help();
     88            exit(EXIT_SUCCESS);
     89          case 'v':
    7490            config->verbose++;
    7591            break;
    76           case 'k':
     92          case 'k':
    7793            if ((sscanf(argv[optind-1], "%d", &config->xKernel) != 1) ||
    7894                (sscanf(argv[optind++], "%d", &config->yKernel) != 1)) {
    … …  
    84100                exit(EXIT_FAILURE);
    85101            }
    86             break;
    87           case 'f':
    88             if (sscanf(optarg, "%d", &config->footprint) != 1) {
    89                 help();
    90                 exit(EXIT_FAILURE);
    91             }
    92             break;
    93           case 't':
    94             if (sscanf(optarg, "%f", &config->threshold) != 1) {
    95                 help();
    96                 exit(EXIT_FAILURE);
    97             }
    98             break;
    99           case 's':
     102            break;
     103          case 'f':
     104            if (sscanf(optarg, "%d", &config->footprint) != 1) {
     105                help();
     106                exit(EXIT_FAILURE);
     107            }
     108            break;
     109          case 't':
     110            if (sscanf(optarg, "%f", &config->threshold) != 1) {
     111                help();
     112                exit(EXIT_FAILURE);
     113            }
     114            break;
     115          case 's':
    100116            if ((sscanf(argv[optind-1], "%f", &config->refSat) != 1) ||
    101117                (sscanf(argv[optind++], "%f", &config->inSat) != 1)) {
    … …  
    107123                exit(EXIT_FAILURE);
    108124            }
    109             break;
    110           case 'o':
    111             if (sscanf(optarg, "%d", &config->spatialOrder) != 1) {
    112                 help();
    113                 exit(EXIT_FAILURE);
    114             }
    115             break;
    116           case 'b':
     125            break;
     126          case 'o':
     127            if (sscanf(optarg, "%d", &config->spatialOrder) != 1) {
     128                help();
     129                exit(EXIT_FAILURE);
     130            }
     131            break;
     132          case 'b':
    117133            if ((sscanf(argv[optind-1], "%f", &config->refBad) != 1) ||
    118134                (sscanf(argv[optind++], "%f", &config->inBad) != 1)) {
    … …  
    124140                exit(EXIT_FAILURE);
    125141            }
    126             break;
    127           case 'B':
    128             if (sscanf(optarg, "%f", &config->background) != 1) {
    129                 help();
    130                 exit(EXIT_FAILURE);
    131             }
    132             break;
    133           case 'n':
     142            break;
     143          case 'B':
     144            if (sscanf(optarg, "%f", &config->background) != 1) {
     145                help();
     146                exit(EXIT_FAILURE);
     147            }
     148            break;
     149          case 'n':
    134150            if (sscanf(argv[optind-1], "%d", &config->nsx) != 1) {
    135151                help();
    136152                exit(EXIT_FAILURE);
    137153            }
    138             config->nsy = config->nsx;
    139             break;
    140           case 'i':
     154            config->nsy = config->nsx;
     155            break;
     156          case 'i':
    141157            if (sscanf(argv[optind-1], "%d", &config->numIter) != 1) {
    142158                help();
    143159                exit(EXIT_FAILURE);
    144160            }
    145             break;
    146           case 'r':
     161            break;
     162          case 'r':
    147163            if (sscanf(argv[optind-1], "%f", &config->sigmaRej) != 1) {
    148164                help();
    149165                exit(EXIT_FAILURE);
    150166            }
    151             break;
    152           case 'p':
     167            break;
     168          case 'p':
    153169            if (sscanf(argv[optind-1], "%f", &config->penalty) != 1) {
    154170                help();
    155171                exit(EXIT_FAILURE);
    156172            }
    157             break;
    158           case 'S':
    159             config->stampFile = argv[optind-1];
    160             config->threshold = -INFINITY;
    161             break;
    162           case 'R':
    163             config->reverse = true;
    164             break;
    165           case 'm':
     173            break;
     174          case 'S':
     175            config->stampFile = argv[optind-1];
     176            config->threshold = -INFINITY;
     177            break;
     178          case 'R':
     179            config->reverse = true;
     180            break;
     181          case 'm':
    166182            config->mask = true;
    167183            break;
    168           default:
    169             help();
    170         }
     184          default:
     185            help();
     186        }
    171187    }
    172188

  • trunk/pois/src/poisReadStamps.c

    r6891 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include "pslib.h"

  • trunk/pois/src/poisRejectStamps.c

    r5717 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include <assert.h>

  • trunk/pois/src/poisSolveEquation.c

    r6891 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include <assert.h>

  • trunk/pois/src/poisStamp.c

    r5717 r8480  
     1// Copyright (C) 2006  Paul A. Price (price@ifa.hawaii.edu)
     2//
     3// This program is free software; you can redistribute it and/or modify
     4// it under the terms of the GNU General Public License as published by
     5// the Free Software Foundation; either version 2 of the License, or
     6// (at your option) any later version.
     7//
     8// This program is distributed in the hope that it will be useful,
     9// but WITHOUT ANY WARRANTY; without even the implied warranty of
     10// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     11// GNU General Public License for more details.
     12//
     13// You should have received a copy of the GNU General Public License
     14// along with this program; if not, write to the Free Software
     15// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
     16
    117#include <stdio.h>
    218#include "pslib.h"
    … …  
    1935{
    2036    if (stamp->matrix) {
    21         psFree(stamp->matrix);
     37        psFree(stamp->matrix);
    2238    }
    2339    if (stamp->vector) {
    24         psFree(stamp->vector);
     40        psFree(stamp->vector);
    2541    }
    2642    psFree(stamp);
    … …  
    3046// Return true if the stamp is OK
    3147bool poisCheckStamp(const psImage *image, // Image to check for threshold
    32                     const psImage *mask, // Mask to check for bad pixels
    33                     int x, int y,       // Pixel coordinates
    34                     const poisConfig *config // Configuration
     48                    const psImage *mask, // Mask to check for bad pixels
     49                    int x, int y,       // Pixel coordinates
     50                    const poisConfig *config // Configuration
    3551    )
    3652{
    3753    if (x < config->footprint + config->xKernel ||
    38         y < config->footprint + config->yKernel ||
    39         x > config->xImage - config->footprint - config->xKernel ||
    40         y > config->yImage - config->footprint - config->yKernel ||
    41         image->data.F32[y][x] < config->threshold) {
    42         psTrace("pois.checkStamp", 9, "Rejecting stamp at %d,%d (border/threshold)\n", x, y);
    43         return false;
     54        y < config->footprint + config->yKernel ||
     55        x > config->xImage - config->footprint - config->xKernel ||
     56        y > config->yImage - config->footprint - config->yKernel ||
     57        image->data.F32[y][x] < config->threshold) {
     58        psTrace("pois.checkStamp", 9, "Rejecting stamp at %d,%d (border/threshold)\n", x, y);
     59        return false;
    4460    }
    4561
    4662    // Check the footprint
    47     bool ok = true;                     // Is the footprint OK?
    48     int footprint = config->footprint;  // Footprint size
     63    bool ok = true;                     // Is the footprint OK?
     64    int footprint = config->footprint;  // Footprint size
    4965    for (int v = y - footprint; v <= y + footprint && ok; v++) {
    50         for (int u = x - footprint; u <= x + footprint && ok; u++) {
    51             if (mask->data.U8[v][u] &
    52                 (POIS_MASK_BAD | POIS_MASK_NEAR_BAD | POIS_MASK_STAMP)) {
    53                 psTrace("pois.checkStamp", 9, "Rejecting stamp at %d,%d (bad footprint)\n", x, y);
    54                 ok = false;
    55             }
    56         }
     66        for (int u = x - footprint; u <= x + footprint && ok; u++) {
     67            if (mask->data.U8[v][u] &
     68                (POIS_MASK_BAD | POIS_MASK_NEAR_BAD | POIS_MASK_STAMP)) {
     69                psTrace("pois.checkStamp", 9, "Rejecting stamp at %d,%d (bad footprint)\n", x, y);
     70                ok = false;
     71            }
     72        }
    5773    }
    58    
     74
    5975    return ok;
    6076}
Note: See TracChangeset for help on using the changeset viewer.