SUBROUTINE EVBIN( STATUS ) *+ * Name: * EVBIN * Purpose: * Bins N lists into an N dimensional object. [N < = 7] * Language: * Starlink Fortran * Type of Module: * ASTERIX task * Invocation: * CALL EVBIN( STATUS ) * Arguments: * STATUS = INTEGER (Given and Returned) * The global status. * Description: * An event dataset is binned up into a binned dataset. Each event list * selected (except QUALITY) will be binned into an output dimension. * * The binning of each axis may be either regularly or irregularly spaced. * If regularly spaced, then bins may be specified either by bin * width, or by the number of bins required. Also, each axis may have * either increasing or decreasing axis values. This is not controlled by * the user, but by the decreasing flag in each list. * * All bins are inclusive of their lower bound and exclusive of their * upper bound. Note that if the axis values are decreasing the lower * bound will have a larger axis value than the upper bound. * * The way QUALITY lists are handled is controlled by 2 parameters QVAL, * and QKEEP. Values present in the quality list > QVAL are treated as bad * quality values. If QKEEP is true, then bad events are written to the * output data array, and the corresponding element of the output quality * array is set to bad (i.e.1). If QKEEP is false, then all bad events are * simply ignored, and no output quality array is produced. * Usage: * evbin {parameter_usage} * Environment Parameters: * INP - Name of input EVDS (UNIV) * LISTS - Index No(s) of input list(s) (_CHAR) * OPT1 - Allow irregular binning. (_LOGICAL) * OPT2 - Select how each regular axis is to be binned (_LOGICAL) * OPT3 - TRUE = specify all regular axes by No of bins(_LOGICAL) * FALSE = specify all regular axes by bin width * QVAL - Event quality > QVAL = bad (_INTEGER) * QKEEP - Produce output QUALITY array? (_LOGICAL) * REG1..7 - Is this axis to be regularly spaced? (_LOGICAL) * USEBINSIZE1..7 - Is this axis to be specified by bin width? (_LOGICAL) * BASE1..7 - Base value if regular axis (_REAL) * BINSIZE1..7 - Actual axis bin width (_REAL) * NBINS1..7 - Actual number of bins required (_INTEGER) * RANGES1..7 - Actual bin limits for irregular bins (_CHAR) * OUT - Name of output dataset (UNIV) * {parameter_name}[pdims] = {parameter_type} ({parameter_access_mode}) * {parameter_description} * Examples: * {routine_example_text} * {routine_example_description} * Pitfalls: * {pitfall_description}... * Notes: * {routine_notes}... * Prior Requirements: * {routine_prior_requirements}... * Side Effects: * {routine_side_effects}... * Algorithm: * {algorithm_description}... * Accuracy: * {routine_accuracy} * Timing: * {routine_timing} * Implementation Status: * {routine_implementation_status} * External Routines Used: * {name_of_facility_or_package}: * {routine_used}... * Implementation Deficiencies: * Limited to 1000 irregular bins per axis. Limit will be less than this * if, on average, more than 15 characters are used to specify each bin - * including delimiters. * Irregular bins must not overlap. * References: * {task_references}... * Keywords: * evbin, usage:public * Copyright: * Copyright (C) University of Birmingham, 1995 * Authors: * PLA: Phillip Andrews (ROSAT, University of Birmingham) * DJA: David J. Allan (Jet-X, University of Birmingham) * {enter_new_authors_here} * History: * 29 Feb 1988 (PLA): * TWODBIN based on IMAGEBIN by Jim Peden. Changes: Counts lists * deleted, bin_upper deleted, code simplified, HDS type conversion * used on mapping, new calc of default bin sizes, new binning * routine. ASTERIX88 structures. * 13 Apr 1988 (PLA): * TWODBIN converted to ASTERIX88 * 20 Apr 1988 (PLA): * Original EVENTBIN started * 15 Apr 1988 (PLA): * Extra HISTORY info added * 14 Oct 1988 (PLA): * X_ axis bin inversion removed * 4 Apr 1989 (PLA): * Rewrite to accomodate irregular bins * 25 Sep 1989 V1.0-7 (PLA): * Altered to allow increasing or decreasing axes * 11 Dec 1989 V1.0-8 (DJA): * Put structure definition in EVBIN_STR. Re-written * to avoid excessive paging * 18 Dec 1989 V1.0-9 (DJA): * Rationalised EVBIN_BIN7Q to avoid naming list arrays explicitly * 8 Jan 1990 V1.0-10 (DJA): * Assumes 0.0 for QUANTUM if not present in list when * the DATA_ARRAY is of floating type, otherwise 1.0 * 17 Jan 1990 V1.0-11 (DJA): * Informs user of number of events after binning * 24 May 1990 V1.2-0 (DJA): * Handles datasets with only one event. Explicitly resets axis * normalisation. * 29 Aug 1990 V1.3-0 (DJA): * BASE<1..7> parameters added. Quality handling improved * 20 Sep 1990 V1.3-1 (DJA): * BASE option allowed when OPT3 true or false * 10 Sep 1991 V1.5-0 (DJA): * Final re-write. History improved, free all objects * and works twice as fast. * 5 Jun 1992 V1.6-0 (DJA): * Fixed bug if irregular axis wasn't first axis * 25 Feb 1994 V1.7-0 (DJA): * Updated quality handling * 24 Nov 1994 V1.8-0 (DJA): * Now use USI for user interface * 20 Apr 1995 V1.8-1 (DJA): * New data interface for output * 19 Sep 1995 V2.0-0 (DJA): * Full ADI port. Removed Fortran structures * 7 Nov 1995 V2.0-1 (DJA): * Allow field extrema to be absent, but warn if so. * {enter_changes_here} * Bugs: * {note_any_bugs_here} *- * Type Definitions: IMPLICIT NONE ! No implicit typing * Global Constants: INCLUDE 'SAE_PAR' ! Standard SAE constants INCLUDE 'ADI_PAR' INCLUDE 'QUAL_PAR' * Status: INTEGER STATUS ! Global status * Local Constants: INTEGER SP_SIZE ! Default # spatial bins PARAMETER ( SP_SIZE = 512 ) INTEGER MX_HTEXT ! Max lines of history text PARAMETER ( MX_HTEXT = ADI__MXDIM + 7 ) INTEGER MXRANGE ! 2 x max No of irregular bins PARAMETER ( MXRANGE = 2000 ) INTEGER MAXNAX ! Max # axes + 1 for PARAMETER ( MAXNAX = ADI__MXDIM+1)! quality CHARACTER*30 VERSION PARAMETER ( VERSION = 'EVBIN Version V2.0-0' ) * Local Variables: CHARACTER*80 ACTION(MX_HTEXT) ! History text CHARACTER*20 LNAME ! List name string CHARACTER*20 OCLASS ! Output object class CHARACTER*20 O_NAME(MAXNAX) ! Output axis name CHARACTER*40 O_UNITS(MAXNAX) ! Output axis units CHARACTER*12 PAR ! Parameter name CHARACTER*80 TEMP ! Temporary store for LABEL ! & UNITS info REAL ARANGE ! Axis range, diff of RHS,LHS REAL O_BSIZE(MAXNAX) ! Regular bin widths REAL O_LHS(MAXNAX) ! Lower field extreme REAL O_RHS(MAXNAX) ! Upper field extreme REAL O_QUANT(MAXNAX) ! Quanta REAL RANGES(MXRANGE,ADI__MXDIM)! Irregular bin ranges. REAL SPARR(2) ! Spaced array values INTEGER AXPTR ! Output axis array INTEGER ALEN ! Action length INTEGER ANUM ! List # from EDI_QFND INTEGER BADQUAL ! Bad event quality ! threshold INTEGER BID ! Output BinDS object INTEGER I, J ! Loop counters INTEGER LBIN(ADI__MXDIM) ! Lists to bin INTEGER IFID ! Input dataset id INTEGER NACT ! # used history lines INTEGER NEVENT ! # input events INTEGER NINDEX ! # lists selected INTEGER INLIST ! # lists in input INTEGER LID ! List identifier INTEGER O_LID(MAXNAX) ! List id per axis INTEGER O_PTR(MAXNAX) ! Mapped list per axis INTEGER OFID ! Output dataset id INTEGER ONDIM ! Output dimensionality INTEGER ODIMS(ADI__MXDIM) ! Output dimensions INTEGER ONELM ! # output elements INTEGER ODPTR ! Output data array INTEGER OQPTR ! Output quality INTEGER WPTR ! Output axis widths LOGICAL OK ! General validity test LOGICAL USEBINSIZE ! Users supplies binsize ? ! (or No bins) LOGICAL IRREG ! Allow irregular binning? LOGICAL IDVSEL ! Select how to bin each ! regular axis individualy. LOGICAL ALLNUM ! If TRUE specify regular axis ! by No of bins. Else by width. LOGICAL O_DECR(MAXNAX) ! Axis decreasing? LOGICAL O_REG(MAXNAX) ! Axis regular? LOGICAL O_QOK(MAXNAX) ! Quantum present? LOGICAL QKEEP ! Produce o/p quality? LOGICAL QUALITY ! QUALITY list present? *. * Check inherited global status. IF ( STATUS .NE. SAI__OK ) RETURN * Version id CALL MSG_PRNT( VERSION ) * Initialise ASTERIX CALL AST_INIT() * Get event dataset CALL USI_ASSOC( 'INP', 'EventDS', 'READ', IFID, STATUS ) IF ( STATUS .NE. SAI__OK ) GOTO 99 * Get number of events and number of lists CALL EDI_GETNS( IFID, NEVENT, INLIST, STATUS ) IF ( INLIST .EQ. 0 ) THEN CALL MSG_PRNT( 'Dataset does not contain any valid lists!' ) GOTO 99 END IF * Display list of list names to user CALL MSG_PRNT( 'The available lists are:' ) CALL EDI_DISP( IFID, STATUS ) IF ( STATUS .NE. SAI__OK ) GOTO 99 * Select the lists CALL MSG_BLNK() CALL MSG_PRNT( 'Select lists to be binned by entering the '// : 'index numbers, eg. 1 2 3') CALL EDI_SELCT( 'LISTS', INLIST, 1, ADI__MXDIM, LBIN, NINDEX, : STATUS ) IF ( STATUS .NE. SAI__OK ) GOTO 99 * Write chosen lists into structure QUALITY = .FALSE. ONDIM = 0 DO I = 1, NINDEX * Extract index element J = LBIN(I) * Index the list CALL EDI_IDX( IFID, J, LID, STATUS ) * Get list name CALL ADI_NAME( LID, LNAME, STATUS ) * Is it the quality list? IF ( LNAME .EQ. 'QUALITY' ) THEN * We don't bin this list, but store at the end of the binnable lists QUALITY = .TRUE. O_NAME(NINDEX) = LNAME O_LID(NINDEX) = LID * Map list DATA_ARRAY CALL EDI_MAPI( IFID, LNAME, 'READ', 0, 0, O_PTR(NINDEX), : STATUS ) ELSE ONDIM = ONDIM + 1 O_NAME(ONDIM) = LNAME O_LID(ONDIM) = LID * Map list DATA_ARRAY CALL EDI_MAPR( IFID, LNAME, 'READ', 0, 0, O_PTR(ONDIM), : STATUS ) * See if UNITS is present CALL ADI_THERE( LID, 'Units', OK, STATUS ) IF ( OK ) THEN CALL ADI_CGET0C( LID, 'Units', O_UNITS(ONDIM), STATUS ) ELSE O_UNITS(ONDIM) = ' ' END IF * See if DECREASING is present CALL ADI_THERE( LID, 'Decreasing', OK, STATUS ) IF ( OK ) THEN CALL ADI_CGET0L( LID, 'Decreasing', O_DECR(ONDIM), STATUS ) ELSE O_DECR(ONDIM) = .FALSE. END IF END IF * Check status IF ( STATUS .NE. SAI__OK ) GOTO 99 END DO * See if ONDIM is OK IF ( ONDIM .LT. 1 ) THEN STATUS = SAI__ERROR CALL ERR_REP( ' ', 'FATAL ERROR: No LIST to bin!', STATUS ) END IF * Obtain user input CALL USI_GET0L( 'OPT1', IRREG, STATUS ) CALL USI_GET0L( 'OPT2', IDVSEL, STATUS ) CALL USI_GET0L( 'OPT3', ALLNUM, STATUS ) IF ( QUALITY ) THEN CALL USI_GET0I( 'QVAL', BADQUAL, STATUS ) CALL USI_GET0L( 'QKEEP', QKEEP, STATUS ) END IF IF ( STATUS .NE. SAI__OK ) GOTO 99 * Get binning info from user DO I = 1, ONDIM * Title for the axis CALL MSG_BLNK() CALL MSG_SETC( 'NAME', O_NAME(I) ) CALL MSG_PRNT( 'Axis ^NAME :' ) * Regular bins for this axis? IF ( IRREG ) THEN * Find out if this axis is to be irregularly binned WRITE( PAR, '(A3,I1)' ) 'REG', I CALL USI_GET0L( PAR, O_REG(I), STATUS ) CALL USI_CANCL( PAR, STATUS ) IF ( STATUS .NE. SAI__OK ) GOTO 99 ELSE O_REG(I) = .TRUE. END IF * Get field range CALL EVBIN_GETRANGE( O_LID(I), O_NAME(I), O_DECR(I), O_QOK(I), : O_QUANT(I), NEVENT, %VAL(O_PTR(I)), : O_LHS(I), O_RHS(I), STATUS ) * Tell user the data range CALL MSG_SETR( 'LHS', O_LHS(I) ) CALL MSG_SETR( 'RHS', O_RHS(I) ) CALL MSG_SETC( 'UNITS', O_UNITS(I) ) CALL MSG_PRNT( ' Data range is ^LHS to ^RHS ^UNITS' ) * and the intrinsic width if present IF ( O_QOK(I) ) THEN CALL MSG_SETR( 'QUANTUM', O_QUANT(I) ) CALL MSG_SETC( 'UNITS', O_UNITS(I) ) CALL MSG_PRNT( ' Intrinsic width is ^QUANTUM ^UNITS' ) END IF IF ( O_REG(I) ) THEN * Find out if BINSIZE or number of bins are to be specified IF ( IDVSEL ) THEN WRITE( PAR, '(A10,I1)' ) 'USEBINSIZE', I CALL USI_GET0L( PAR, USEBINSIZE, STATUS ) CALL USI_CANCL( PAR, STATUS ) IF ( STATUS .NE. SAI__OK ) GOTO 99 ELSE IF ( ALLNUM ) THEN USEBINSIZE = .FALSE. ELSE USEBINSIZE = .TRUE. END IF * Is user going to override bin base? WRITE( PAR, '(A4,I1)' ) 'BASE', I CALL USI_DEF0R( PAR, O_LHS(I), STATUS ) CALL USI_GET0R( PAR, O_LHS(I), STATUS ) IF ( STATUS .NE. SAI__OK ) GOTO 99 * Set range based on left and right bounds ARANGE = ABS( O_RHS(I) - O_LHS(I) ) IF ( USEBINSIZE ) THEN * Anyone think of a more sensible default? O_BSIZE(I) = ARANGE / SP_SIZE IF ( O_QOK(I) .AND. (O_BSIZE(I) .LT. O_QUANT(I)) ) THEN O_BSIZE(I) = O_QUANT(I) END IF * Get binsize from user WRITE( PAR, '(A7,I1)' ) 'BINSIZE', I CALL USI_DEF0R( PAR, O_BSIZE(I), STATUS ) CALL USI_GET0R( PAR, O_BSIZE(I), STATUS ) CALL USI_CANCL( PAR, STATUS ) IF ( STATUS .NE. SAI__OK ) GOTO 99 IF ( O_QOK(I) .AND. (O_BSIZE(I) .LT. O_QUANT(I)) ) THEN CALL MSG_PRNT( 'WARNING: Bin size less than QUANTUM' ) END IF ODIMS(I) = NINT( ARANGE / O_BSIZE(I)) IF ( ODIMS(I) .GT. 1 ) THEN CALL MSG_SETI( 'NBIN', ODIMS(I) ) CALL MSG_SETC( 'NAME', O_NAME(I) ) CALL MSG_PRNT( 'This will give ^NBIN ^NAME bins' ) ELSE CALL MSG_PRNT( 'WARNING: You have chosen only one bin!' ) ODIMS(I) = 1 END IF ELSE WRITE( PAR, '(A5,I1)' ) 'NBINS', I ODIMS(I) = 512 O_BSIZE(I) = ARANGE / REAL(ODIMS(I)) IF ( O_QOK(I) .AND. (O_BSIZE(I) .LT. O_QUANT(I)) ) THEN ODIMS(I) = INT(ARANGE / O_QUANT(I) ) END IF CALL USI_DEF0I( PAR, ODIMS(I), STATUS ) CALL USI_GET0I( PAR, ODIMS(I), STATUS ) CALL USI_CANCL( PAR, STATUS ) IF ( STATUS .NE. SAI__OK ) GOTO 99 IF ( ODIMS(I) .LE. 1 ) THEN CALL MSG_PRNT( 'WARNING: You have chosen only one bin !') ODIMS(I) = 1 END IF O_BSIZE(I) = ARANGE / REAL(ODIMS(I)) IF ( O_QOK(I) .AND. (O_BSIZE(I) .LT. O_QUANT(I)) ) THEN CALL MSG_PRNT( 'WARNING: Bin size less than QUANTUM' ) END IF CALL MSG_SETR( 'BSIZE', O_BSIZE(I) ) CALL MSG_SETC( 'UNITS', O_UNITS(I) ) CALL MSG_PRNT( 'This will give a bin width of ^BSIZE '// : '^UNITS' ) END IF * Set sign of bin size IF ( O_DECR(I) ) O_BSIZE(I) = - O_BSIZE(I) * Irregular axis ELSE IF ( O_DECR(I) ) THEN CALL MSG_PRNT( 'You must specify DECREASING ranges e.g. '// : '30:20:10') ELSE CALL MSG_PRNT( 'You must specify INCREASING ranges e.g. '// : '10:20:30') END IF * Get bin boundaries WRITE( PAR, '(A6,I1)' ) 'RANGES', I CALL PRS_GETRANGES( PAR, MXRANGE, 1, O_RHS(I), O_LHS(I), : RANGES(1,I), ODIMS(I), STATUS ) IF ( STATUS .NE. SAI__OK ) GOTO 99 END IF END DO * Check status IF ( STATUS .NE. SAI__OK ) GOTO 99 * Try to guess o/p dataset class my looking at the lists which have been * chosen OCLASS = 'BinDS' IF ( ONDIM .EQ. 1 ) THEN CALL EDI_QFND( IFID, 'T', LNAME, ANUM, STATUS ) IF ( ANUM .EQ. LBIN(1) ) THEN OCLASS = 'TimeSeries' ELSE CALL EDI_QFND( IFID, 'E', LNAME, ANUM, STATUS ) IF ( ANUM .EQ. LBIN(1) ) THEN OCLASS = 'Spectrum' END IF END IF ELSE IF ( ONDIM .EQ. 2 ) THEN CALL EDI_QFND( IFID, 'X', LNAME, ANUM, STATUS ) IF ( ANUM .EQ. LBIN(1) ) THEN CALL EDI_QFND( IFID, 'Y', LNAME, ANUM, STATUS ) IF ( ANUM .EQ. LBIN(2) ) THEN OCLASS = 'XYimage' END IF END IF END IF CALL USI_DEF0C( 'OUTFORM', OCLASS, STATUS ) CALL USI_GET0C( 'OUTFORM', OCLASS, STATUS ) IF ( STATUS .NE. SAI__OK ) GOTO 99 * Create output dataset CALL MSG_BLNK() CALL BDI_NEW( OCLASS, ONDIM, ODIMS, 'REAL', BID, STATUS ) CALL USI_CREAT( 'OUT', BID, OFID, STATUS ) * Find total number of output elements CALL ARR_SUMDIM( ONDIM, ODIMS, ONELM ) * Map and initialise data CALL BDI_MAPR( OFID, 'Data', 'WRITE/ZERO', ODPTR, STATUS ) * Map quality IF ( QUALITY .AND. QKEEP ) THEN CALL BDI_MAPUB( OFID, 'Quality', 'WRITE/QMISSING', OQPTR, : STATUS ) CALL BDI_PUT0UB( OFID, 'QualityMask', QUAL__MASK, STATUS ) END IF IF ( STATUS .NE. SAI__OK ) GOTO 99 * Loop over AXIS structure writing the values. DO I = 1, ONDIM * Write text CALL BDI_AXPUT0C( OFID, I, 'Label', O_NAME(I), STATUS ) CALL BDI_AXPUT0C( OFID, I, 'Units', O_UNITS(I), STATUS ) CALL BDI_AXPUT0L( OFID, I, 'Normalised', .FALSE., STATUS ) * Write bin characteristics IF ( O_REG(I) ) THEN * Find leftmost bin centre and signed bin width SPARR(1) = O_LHS(I) + O_BSIZE(I)/2.0 SPARR(2) = O_BSIZE(I) CALL BDI_AXPUT1R( OFID, I, 'SpacedData', 2, SPARR, STATUS ) ELSE * Centres and widths if irregular CALL BDI_AXMAPR( OFID, I, 'Data', 'WRITE', AXPTR, STATUS ) CALL BDI_AXMAPR( OFID, I, 'Width', 'WRITE', WPTR, STATUS ) CALL ARR_BND2CWR( ODIMS(I), RANGES, %VAL(AXPTR), : %VAL(WPTR), STATUS ) CALL BDI_AXUNMAP( OFID, I, 'Data', AXPTR, STATUS ) CALL BDI_AXUNMAP( OFID, I, 'Width', WPTR, STATUS ) END IF END DO * Copy header info CALL ADI_CGET0C( IFID, 'Title', TEMP, STATUS ) IF ( STATUS .NE. SAI__OK ) THEN CALL ERR_ANNUL( STATUS ) TEMP = ' ' END IF IF ( TEMP .GT. ' ' ) THEN CALL BDI_PUT0C( OFID, 'Title', TEMP, STATUS ) END IF CALL BDI_PUT0C( OFID, 'Label', 'Events', STATUS ) * Copy the ancillaries CALL UDI_COPANC( IFID, ' ', OFID, STATUS ) IF ( STATUS .NE. SAI__OK ) GOTO 99 * Pad dimensions with ones to simulate 7D data CALL AR7_PAD( ONDIM, ODIMS, STATUS ) * Bin the data CALL EVBIN_INT( ONDIM, O_REG, O_LHS, O_BSIZE, O_DECR, O_PTR, : ODIMS, ODIMS(1), ODIMS(2), ODIMS(3), : ODIMS(4), ODIMS(5), ODIMS(6), ODIMS(7), : NEVENT, RANGES, %VAL(O_PTR(MAXNAX)), BADQUAL, : (QKEEP.AND.QUALITY), %VAL(OQPTR), : %VAL(ODPTR), STATUS ) * Copy history from input CALL HSI_COPY( IFID, OFID, STATUS ) CALL HSI_ADD( OFID, VERSION, STATUS ) * Write essential data to history * ...Input filename ACTION(1) = 'Input {INP}' NACT = MX_HTEXT CALL USI_TEXT( 1, ACTION, NACT, STATUS ) * ...List names NACT = NACT + 1 ACTION(NACT) = ' ' NACT = NACT + 1 IF ( ONDIM .GT. 1 ) THEN ACTION(NACT) = 'Dataset created from the lists,' DO I = 1, ONDIM NACT = NACT + 1 WRITE( ACTION(NACT), '(15X,A)') O_NAME(I) END DO ELSE ACTION(NACT) = 'Dataset created from the list '//O_NAME(1) END IF * ...Number of events NACT = NACT + 1 ACTION(NACT) = ' ' NACT = NACT + 1 CALL MSG_SETI( 'NEV', NEVENT ) CALL MSG_MAKE( 'Containing data on ^NEV events', ACTION(NACT), : ALEN ) * ...Quality processing IF ( QUALITY ) THEN NACT = NACT + 1 IF ( QKEEP ) THEN ACTION(NACT) = 'QUALITY array created' ELSE ACTION(NACT) = 'Bad quality events excluded' END IF END IF * Write the text CALL HSI_PTXT( OFID, NACT, ACTION, STATUS ) * Free input lists DO I = 1, NINDEX CALL EDI_UNMAP( IFID, O_NAME(I), STATUS ) CALL ADI_ERASE( O_LID(I), STATUS ) END DO * Free output CALL USI_ANNUL( 'OUT', STATUS ) * Tidy up 99 CALL AST_CLOSE CALL AST_ERR( STATUS ) END *+ EVBIN_GETRANGE - Obtains min, max & range of data SUBROUTINE EVBIN_GETRANGE( LID, NAME, DECR, QOK, QUANT, LLEN, : LDAT, LHS, RHS, STATUS ) * Description : * Obtains FIELD_MIN/MAX & range for the list pointed to by AXN. * QUANTUM is taken into account if present. * Method : * Obtains field_min, field_max from dataset & adjusts using * QUANTUM. * Deficiencies : * Only uses QUANTUM if it is SCALAR. * Bugs : * Authors : * Phil Andrews (BHVAD::PLA) * History : * * 23-MAR-1989 : Original (PLA) * * Type Definitions : * IMPLICIT NONE * * Global constants : * INCLUDE 'SAE_PAR' * * Import : * INTEGER LID ! List identifier CHARACTER*(*) NAME ! List name LOGICAL DECR INTEGER LLEN REAL LDAT(*) * * Export : * LOGICAL QOK REAL QUANT REAL LHS, RHS ! Field extrema * * Status : * INTEGER STATUS * * Local variables : * CHARACTER*8 TYPE ! Basic list type REAL DIR ! Direction vector REAL FMIN, FMAX ! Field extrema LOGICAL QVEC ! Vector quantum? *- * Check inherited global status IF ( STATUS .NE. SAI__OK ) RETURN * Try to get field extrema CALL ADI_CGET0R( LID, 'Min', FMIN, STATUS ) CALL ADI_CGET0R( LID, 'Max', FMAX, STATUS ) IF ( STATUS .NE. SAI__OK ) THEN CALL ERR_ANNUL( STATUS ) CALL MSG_SETC( 'NAME', NAME ) CALL MSG_PRNT( 'WARNING : Error reading field min/max for '/ : /'list ^NAME, data range will be used instead' ) CALL ARR_RANG1R( LLEN, LDAT, FMIN, FMAX, STATUS ) END IF * Does quantum exist? CALL ADI_THERE( LID, 'Quantum', QOK, STATUS ) IF ( QOK ) THEN * Make sure its a scalar quantum CALL ADI_CGET0L( LID, 'VectorQuantum', QVEC, STATUS ) IF ( QVEC ) THEN QOK = .FALSE. ELSE CALL ADI_CGET0R( LID, 'Quantum', QUANT, STATUS ) END IF ELSE CALL ADI_CGET0C( LID, 'TYPE', TYPE, STATUS ) IF ( (TYPE(1:4) .EQ. 'REAL') .OR. : (TYPE(1:6) .EQ. 'DOUBLE') ) THEN QUANT = 0.0 ELSE QUANT = 1.0 QOK = .TRUE. END IF CALL MSG_SETC( 'LIST', NAME ) CALL MSG_SETR( 'VAL', QUANT ) CALL MSG_PRNT( 'WARNING : No QUANTUM component in list '/ : /'^LIST, have assumed ^VAL' ) END IF * Set left and right centres of primitive bins IF ( DECR ) THEN LHS = FMAX RHS = FMIN ELSE LHS = FMIN RHS = FMAX END IF * Convert from centres to edges if quantum present IF ( QOK ) THEN IF ( DECR ) THEN DIR = -1.0 ELSE DIR = 1.0 END IF LHS = LHS - DIR*QUANT/2.0 RHS = RHS + DIR*QUANT/2.0 END IF * Exit 99 IF ( STATUS .NE. SAI__OK ) THEN CALL AST_REXIT( 'EVBIN_GETRANGE', STATUS ) END IF END *+ EVBIN_REGIDX - Return output pixel SUBROUTINE EVBIN_REGIDX( LHS, BSIZE, NBIN, IN, VALID, INDEX ) * Description : * Calculates the output pixel value * Authors : * Phil Andrews (BHVAD::PLA) * History : * * 23 Mar 89 : Original (PLA) * * Type Definitions : * IMPLICIT NONE * * Global constants : * INCLUDE 'SAE_PAR' * * Import : * REAL LHS, BSIZE INTEGER NBIN REAL IN ! Input value * * Export : * LOGICAL VALID ! Is output pixel valid INTEGER INDEX ! Index *- INDEX = 1 + INT((IN - LHS) / BSIZE) VALID = (INDEX .LE. NBIN .AND. INDEX .GT. 0) END *+ EVBIN_IRREGIDX - Returns output pixel for irregularly binned axis SUBROUTINE EVBIN_IRREGIDX( DECR, START, STOP, INC1, RANGES, IN, : VALID, INDEX ) * Description : * Authors : * Phil Andrews (BHVAD::PLA) * History : * * 23 Mar 89 : Original (PLA) * * Type Definitions : * IMPLICIT NONE * * Global constants : * INCLUDE 'SAE_PAR' * * Import : * LOGICAL DECR ! Values decreasing? INTEGER START ! Initial value for J INTEGER STOP ! Last element of RANGES INTEGER INC1 ! initial value for INC REAL RANGES(2000) ! Irregular bins REAL IN ! Event value * Export : LOGICAL VALID ! Is output pixel valid INTEGER INDEX ! Index * Local variables : INTEGER J ! Current search position INTEGER INC ! Increment for J *- IF (.NOT. DECR ) THEN IF (RANGES(1) .LE. IN .AND. IN .LT. RANGES(STOP)) THEN VALID = .TRUE. INDEX = 0 J = START INC = INC1 * 2 DO WHILE (INDEX .EQ. 0) IF (IN .LT. RANGES(j)) THEN J = J - INC IF (INC .GT. 2) INC = INC / 2 ELSE IF (IN .GE. RANGES(J+1)) THEN J = J + INC IF (INC .GT. 2) INC = INC / 2 ELSE INDEX = 1 + ((J - 1) / 2) END IF END DO ELSE VALID = .FALSE. END IF ELSE IF (RANGES(1) .GE. IN .AND. IN .GT. RANGES(STOP)) THEN VALID = .TRUE. INDEX = 0 J = START INC = INC1 * 2 DO WHILE (INDEX .EQ. 0) IF (IN .GT. RANGES(j)) THEN J = J - INC IF (INC .GT. 2) INC = INC / 2 ELSE IF (IN .LE. RANGES(J+1)) THEN J = J + INC IF (INC .GT. 2) INC = INC / 2 ELSE INDEX = 1 + ((J - 1) / 2) END IF END DO ELSE VALID = .FALSE. END IF END IF END *+ EVBIN_INT - Performs binning operation for 7D output with QUALITY. SUBROUTINE EVBIN_INT( NDIM, O_REG, O_LHS, O_BSIZE, O_DECR, O_PTR, : ODIMS, L1, L2, L3, L4, L5, L6, L7, : EVENTS, RANGES, QIN, BADQUAL, QKEEP, : QUAL, OUT, STATUS ) * Description : * Bins the event information in IN and writes it to OUT * Method : * Deficiencies : * Bugs : * Authors : * * David J. Allan ( BHVAD::DJA ) * * History : * * 11 Dec 89 : Original (DJA) * 6 Jul 92 : removed intrinsic fn so that it compiles on SUN (RDS) * * Type Definitions : * IMPLICIT NONE * * Global constants : * INCLUDE 'SAE_PAR' INCLUDE 'ADI_PAR' INCLUDE 'PRM_PAR' INCLUDE 'QUAL_PAR' * * Import : * INTEGER NDIM ! Output dimensionality LOGICAL O_REG(*) ! Axis is regular? REAL O_LHS(*) ! 1st bin lower bounds REAL O_BSIZE(*) ! Regular bin widths LOGICAL O_DECR(*) INTEGER O_PTR(*) INTEGER ODIMS(NDIM) ! Dimensions array INTEGER L1,L2,L3,L4,L5,L6,L7 INTEGER EVENTS ! Number of events in input lists INTEGER BADQUAL ! Ignore events with quality ! > this value. INTEGER QIN(EVENTS) ! QUALITY list data REAL RANGES(2000,ADI__MXDIM) ! Irregular bins LOGICAL QKEEP ! Keep bad quality events? * * Export : * BYTE QUAL(L1,L2,L3,L4,L5,L6,L7) REAL OUT (L1,L2,L3,L4,L5,L6,L7) * Status: INTEGER STATUS * Local Variables : INTEGER I,J ! loop counters INTEGER BN(ADI__MXDIM) ! Index to OUT INTEGER START(ADI__MXDIM)! Start bin for search if irreg INTEGER INC1(ADI__MXDIM) ! Initial increment for search INTEGER STOP(ADI__MXDIM) ! No of elements in RANGES(x,n) INTEGER TBIN ! Total bin counter LOGICAL ALL_REG ! All axes regular? LOGICAL VALID ! Has a valid bin been found? * * Inline functions : * INTEGER EVD *- * Check inherited global status IF ( STATUS .NE. SAI__OK ) RETURN * Boundaries for irregular axes ALL_REG = .TRUE. DO I = 1, NDIM IF ( .NOT. O_REG(I) ) THEN START(I) = 1 + (2 * ((ODIMS(I) / 2) - 1)) STOP(I) = 2 * ODIMS(I) INC1(I) = START (I) / 2 ALL_REG = .FALSE. END IF END DO * Initialise BN array for dummy dimensions CALL AR7_PAD( NDIM, BN, STATUS ) * Initialise total bin counter TBIN = 0 * Loop over the event lists IF ( QKEEP ) THEN IF ( ALL_REG ) THEN DO I = 1, EVENTS * Find bin on each dimension DO J = 1, NDIM EVD = O_PTR(J) + (I-1)*VAL__NBR CALL EVBIN_REGIDX( O_LHS(J), O_BSIZE(J), ODIMS(J), : %VAL(EVD), VALID, BN(J) ) IF ( .NOT. VALID ) GOTO 10 END DO * Bump up bin counter OUT(BN(1),BN(2),BN(3),BN(4),BN(5),BN(6),BN(7)) = : OUT(BN(1),BN(2),BN(3),BN(4),BN(5),BN(6),BN(7))+1.0 TBIN = TBIN + 1 IF ( QIN(I) .GT. BADQUAL ) THEN QUAL(BN(1),BN(2),BN(3),BN(4),BN(5),BN(6),BN(7)) = : QUAL__BAD END IF 10 CONTINUE END DO ELSE DO I = 1, EVENTS * Find bin on each dimension DO J = 1, NDIM EVD = O_PTR(J) + (I-1)*VAL__NBR IF ( O_REG(J) ) THEN CALL EVBIN_REGIDX( O_LHS(J), O_BSIZE(J), ODIMS(J), : %VAL(EVD), VALID, BN(J) ) ELSE CALL EVBIN_IRREGIDX( O_DECR(J), START(J),STOP(J), : INC1(J),RANGES(1,J), %VAL(EVD), VALID, BN(J) ) END IF IF ( .NOT. VALID ) GOTO 20 END DO * Bump up bin counter OUT(BN(1),BN(2),BN(3),BN(4),BN(5),BN(6),BN(7)) = : OUT(BN(1),BN(2),BN(3),BN(4),BN(5),BN(6),BN(7))+1.0 QUAL(BN(1),BN(2),BN(3),BN(4),BN(5),BN(6),BN(7)) = : QUAL__GOOD TBIN = TBIN + 1 20 CONTINUE END DO END IF ELSE IF ( ALL_REG ) THEN DO I = 1, EVENTS * Find bin on each dimension DO J = 1, NDIM EVD = O_PTR(J) + (I-1)*VAL__NBR CALL EVBIN_REGIDX( O_LHS(J), O_BSIZE(J), ODIMS(J), : %VAL(EVD), VALID, BN(J) ) IF ( .NOT. VALID ) GOTO 30 END DO * Bump up bin counter OUT( BN(1),BN(2),BN(3),BN(4),BN(5),BN(6),BN(7) ) = : OUT( BN(1),BN(2),BN(3),BN(4),BN(5),BN(6),BN(7) )+1.0 TBIN = TBIN + 1 30 CONTINUE END DO ELSE DO I = 1, EVENTS * Find bin on each dimension DO J = 1, NDIM EVD = O_PTR(J) + (I-1)*VAL__NBR IF ( O_REG(J) ) THEN CALL EVBIN_REGIDX( O_LHS(J), O_BSIZE(J), ODIMS(J), : %VAL(EVD), VALID, BN(J) ) ELSE CALL EVBIN_IRREGIDX( O_DECR(J), START(J),STOP(J), : INC1(J),RANGES(1,J), %VAL(EVD), VALID, BN(J) ) END IF IF ( .NOT. VALID ) GOTO 40 END DO * Bump up bin counter OUT( BN(1),BN(2),BN(3),BN(4),BN(5),BN(6),BN(7) ) = : OUT( BN(1),BN(2),BN(3),BN(4),BN(5),BN(6),BN(7) )+1.0 TBIN = TBIN + 1 40 CONTINUE END DO END IF END IF IF ( TBIN .NE. 0 ) THEN CALL MSG_SETI( 'N', TBIN ) IF ( TBIN .NE. EVENTS ) THEN CALL MSG_SETI( 'NP', EVENTS ) CALL MSG_PRNT( 'A total of ^N events were binned out'/ : /' of ^NP input' ) ELSE CALL MSG_PRNT( 'A total of ^N events were binned' ) END IF END IF END