IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Sep 19, 2014, 4:05:27 PM (12 years ago)
Author:
eugene
Message:

merge changes from trunk

Location:
branches/eam_branches/ps2-tc3-20130727
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ps2-tc3-20130727

  • branches/eam_branches/ps2-tc3-20130727/psLib/src/fits/psFitsScale.c

    r36680 r37403  
    929929
    930930
    931 psImage *psFitsScaleForDisk(const psImage *image, const psFits *fits, double bscale, double bzero, double boffset, double bsoften,
    932                             psRandom *rng)
     931psImage *psFitsScaleForDisk(const psImage *image, const psFits *fits, double bscale, double bzero, double boffset, double bsoften, psRandom *rng)
    933932{
    934933    PS_ASSERT_IMAGE_NON_NULL(image, NULL);
    … …  
    971970    psImage *out = psImageAlloc(numCols, numRows, outType); // Output image
    972971
    973     if (!psMemIncrRefCounter(rng) && options->fuzz) {
    974         // Don't blab about which seed we're going to get --- it's not necessary for this purpose
    975         rng = psRandomAlloc(PS_RANDOM_TAUS);
    976     }
    977 
    978 #define SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, OUTTYPE) \
    979     case PS_TYPE_##OUTTYPE: { \
    980         double scale = 1.0 / bscale; \
    981         double zero = bzero; \
     972    // fprintf (stderr, "corruption check psFitsScaleForDisk 1\n");
     973    // psMemCheckCorruption (stderr, true);
     974
     975    if (options->fuzz && !psMemIncrRefCounter (rng)) {
     976      // Don't blab about which seed we're going to get --- it's not necessary for this purpose
     977      rng = psRandomAlloc(PS_RANDOM_TAUS);
     978    }
     979
     980#define SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, OUTTYPE)                  \
     981    case PS_TYPE_##OUTTYPE: {                                           \
     982        double scale = 1.0 / bscale;                                    \
     983        double zero = bzero;                                            \
    982984        /* Note: BITPIX=8 treated differently, since it uses unsigned values; the rest use signed */ \
    983         double min = bitpix == 8 ? 0 : -pow(2.0, options->bitpix - 1); \
     985        double min = bitpix == 8 ? 0 : -pow(2.0, options->bitpix - 1);  \
    984986        double max = bitpix == 8 ? 255 : (pow(2.0, options->bitpix - 1) - 1.0); \
    985         for (int y = 0; y < numRows; y++) { \
    986             for (int x = 0; x < numCols; x++) { \
    987               ps##INTYPE value; \
    988               if ((options->scaling == PS_FITS_SCALE_LOG_RANGE)||       \
    989                 (options->scaling == PS_FITS_SCALE_LOG_MANUAL)||        \
    990                   (options->scaling == PS_FITS_SCALE_LOG_STDEV_POSITIVE)|| \
    991                   (options->scaling == PS_FITS_SCALE_LOG_STDEV_NEGATIVE)|| \
    992                   (options->scaling == PS_FITS_SCALE_LOG_STDEV_BOTH)) { \
    993                 if (isfinite((IN)->data.INTYPE[y][x])) {                \
    994                   value = log10( (IN)->data.INTYPE[y][x] - boffset );   \
     987        for (int y = 0; y < numRows; y++) {                             \
     988            for (int x = 0; x < numCols; x++) {                         \
     989                ps##INTYPE value;                                       \
     990                if ((options->scaling == PS_FITS_SCALE_LOG_RANGE)||     \
     991                    (options->scaling == PS_FITS_SCALE_LOG_MANUAL)||    \
     992                    (options->scaling == PS_FITS_SCALE_LOG_STDEV_POSITIVE)|| \
     993                    (options->scaling == PS_FITS_SCALE_LOG_STDEV_NEGATIVE)|| \
     994                    (options->scaling == PS_FITS_SCALE_LOG_STDEV_BOTH)) { \
     995                    if (isfinite((IN)->data.INTYPE[y][x])) {            \
     996                        value = log10( (IN)->data.INTYPE[y][x] - boffset ); \
     997                    }                                                   \
     998                    else {                                              \
     999                        value = (IN)->data.INTYPE[y][x];                \
     1000                    }                                                   \
     1001                }                                                       \
     1002                else if ((options->scaling == PS_FITS_SCALE_ASINH_RANGE)|| \
     1003                         (options->scaling == PS_FITS_SCALE_ASINH_MANUAL)|| \
     1004                         (options->scaling == PS_FITS_SCALE_ASINH_STDEV_POSITIVE)|| \
     1005                         (options->scaling == PS_FITS_SCALE_ASINH_STDEV_NEGATIVE)|| \
     1006                         (options->scaling == PS_FITS_SCALE_ASINH_STDEV_BOTH)) { \
     1007                    if (isfinite((IN)->data.INTYPE[y][x])) {            \
     1008                        value = 1.0857362 * (asinh( ((IN)->data.INTYPE[y][x] - boffset) / (2.0 * bsoften))); \
     1009                    }                                                   \
     1010                    else {                                              \
     1011                        value = (IN)->data.INTYPE[y][x];                \
     1012                    }                                                   \
    9951013                }                                                       \
    9961014                else {                                                  \
    997                   value = (IN)->data.INTYPE[y][x];                      \
     1015                    value = (IN)->data.INTYPE[y][x];                    \
    9981016                }                                                       \
    999               }                                                         \
    1000               else if ((options->scaling == PS_FITS_SCALE_ASINH_RANGE)||        \
    1001                 (options->scaling == PS_FITS_SCALE_ASINH_MANUAL)||        \
    1002                   (options->scaling == PS_FITS_SCALE_ASINH_STDEV_POSITIVE)|| \
    1003                   (options->scaling == PS_FITS_SCALE_ASINH_STDEV_NEGATIVE)|| \
    1004                   (options->scaling == PS_FITS_SCALE_ASINH_STDEV_BOTH)) {       \
    1005                 if (isfinite((IN)->data.INTYPE[y][x])) {                \
    1006                   value = 1.0857362 * (asinh( ((IN)->data.INTYPE[y][x] - boffset) / (2.0 * bsoften))); \
     1017                if (!isfinite(value)) {                                 \
     1018                    /* This choice of "max" for non-finite pixels is mainly cosmetic --- it has to be */ \
     1019                    /* something, and "min" would produce holes in the cores of bright stars. */ \
     1020                    (OUT)->data.OUTTYPE[y][x] = max;                    \
     1021                } else {                                                \
     1022                    value = (value - zero) * scale;                     \
     1023                    if (options->fuzz && (value - (int)value != 0.0)) { \
     1024                        /* Add random factor [0.0,1.0): adds a variance of 1/12, */ \
     1025                        /* but preserves the expectation value given the floor() */ \
     1026                        value += psRandomUniform(rng) ;                 \
     1027                    }                                                   \
     1028                    /* Check for underflow and overflow; set either to max */ \
     1029                    (OUT)->data.OUTTYPE[y][x] = (value < min || value > max ? max : floor(value)); \
    10071030                }                                                       \
    1008                 else {                                                  \
    1009                   value = (IN)->data.INTYPE[y][x];                      \
    1010                 }                                                       \
    1011               }                                                         \
    1012                 else { \
    1013                   value = (IN)->data.INTYPE[y][x];      \
    1014                 }                                           \
    1015                 if (!isfinite(value)) { \
    1016                     /* This choice of "max" for non-finite pixels is mainly cosmetic --- it has to be */ \
    1017                     /* something, and "min" would produce holes in the cores of bright stars. */ \
    1018                     (OUT)->data.OUTTYPE[y][x] = max; \
    1019                 } else { \
    1020                     value = (value - zero) * scale; \
    1021                     if (options->fuzz && (value - (int)value != 0.0)) { \
    1022                        /* Add random factor [-0.5,0.5): adds a variance of 1/12, */ \
    1023                        /* but preserves the expectation value */ \
    1024                         value += psRandomUniform(rng) - 0.5; \
    1025                     } \
    1026                     /* Check for underflow and overflow; set either to max */ \
    1027                     (OUT)->data.OUTTYPE[y][x] = (value < min || value > max ? max : value); \
    1028                 } \
    1029             } \
    1030         } \
    1031         break; \
    1032     }
    1033 
    1034 #define SCALE_WRITE_IN_CASE(IN, INTYPE, OUT) \
    1035     case PS_TYPE_##INTYPE: { \
    1036         switch (outType) { \
    1037           SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, U8);;   \
    1038           SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S16);;  \
    1039           SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S32);;  \
    1040           SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S64);;  \
    1041           default: \
    1042             psAbort("Should be unreachable."); \
    1043         } \
    1044         break; \
    1045     }
    1046 
     1031            }                                                           \
     1032        }                                                               \
     1033        break;                                                          \
     1034    }
     1035   
     1036#define SCALE_WRITE_IN_CASE(IN, INTYPE, OUT)                    \
     1037    case PS_TYPE_##INTYPE: {                                    \
     1038        switch (outType) {                                      \
     1039            SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, U8);;         \
     1040            SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S16);;        \
     1041            SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S32);;        \
     1042            SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S64);;        \
     1043          default:                                              \
     1044            psAbort("Should be unreachable.");                  \
     1045        }                                                       \
     1046        break;                                                  \
     1047    }
     1048   
    10471049    switch (image->type.type) {
    1048         SCALE_WRITE_IN_CASE(image, F32, out);
    1049         SCALE_WRITE_IN_CASE(image, F64, out);
     1050        SCALE_WRITE_IN_CASE(image, F32, out);
     1051        SCALE_WRITE_IN_CASE(image, F64, out);
    10501052      default:
    10511053        psAbort("Should be unreachable.");
    … …  
    10531055
    10541056    psFree(rng);
    1055 
    10561057    return out;
    10571058}
Note: See TracChangeset for help on using the changeset viewer.