
<<HistogramsImpl.Mesa>>
    <<Last tweaked by Mike Spreitzer on March 16, 1989 4:26:26 pm PST>>

    DIRECTORY Histograms, HistogramsExtras, HistogramsPrivate, IO, Real, RealFns, Rope;

    HistogramsImpl: CEDAR MONITOR
        LOCKS h USING h: Histogram
        IMPORTS HistogramsPrivate, IO, Real, RealFns, Rope
        EXPORTS Histograms, HistogramsExtras, HistogramsPrivate =

    BEGIN OPEN Histograms, HistogramsExtras, HistogramsPrivate;

    Error: PUBLIC ERROR [msg: ROPE] = CODE;

    Histogram: TYPE = REF HistogramRep;
    HistogramRep: PUBLIC TYPE = HistogramsPrivate.HistogramRep;

    BinNameProc: TYPE ~ Histograms.BinNameProc;

    Create1D: PUBLIC PROC
        [
            factor: REAL _ 1.0,
            offset: REAL _ 0.0,
            logarithmic: BOOL _ FALSE,
            BinNamer: BinNameProc _ NIL,
            clientData: REF ANY _ NIL]
        RETURNS [h: Histogram] =
        BEGIN
        d: DataRef _ NEW[Data[0]];
        h _ NEW [HistogramRep _ [
            dimensionality: 1,
            BinNamer: BinNamer,
            clientData: clientData,
            data: d,
            iFactor: factor, iOffset: offset, jFactor: 1, jOffset: 0,
            xChg: factor, x0: offset, yChg: 1, y0: 0,
            log: [logarithmic, FALSE]
            ]];
        IF h.log[X] THEN {h.iOffset _ RealFns.Ln[offset]; h.iFactor _ RealFns.Ln[factor]};
        END;

    Create2D: PUBLIC PROC
        [    --sorry, no labelling yet
            iFactor, jFactor: REAL _ 1.0,
            iOffset, jOffset: REAL _ 0.0,
            logI, logJ: BOOL _ FALSE,
            BinNamer: BinNameProc _ NIL,
            clientData: REF ANY _ NIL]
        RETURNS [h: Histogram] ~ {
        h _ NEW [HistogramRep _ [
            dimensionality: 2,
            BinNamer: BinNamer,
            clientData: clientData,
            data: NIL,
            nI: 1, nJ: 1,
            iFactor: iFactor, iOffset: iOffset, jFactor: jFactor, jOffset: jOffset,
            xChg: iFactor, x0: iOffset, yChg: jFactor, y0: jOffset,
            log: [logI, logJ]
            ]];
        h.data _ NEW [Data[h.nI * h.nJ]];
        FOR k: NAT IN [0 .. h.data.length) DO h.data[k] _ 0 ENDLOOP;
        IF h.log[X] THEN {h.iOffset_RealFns.Ln[iOffset]; h.iFactor_RealFns.Ln[iFactor]};
        IF h.log[Y] THEN {h.jOffset_RealFns.Ln[jOffset]; h.jFactor_RealFns.Ln[jFactor]};
        RETURN};

    Ensure1: PUBLIC PROC [old: DataRef, n: INT] RETURNS [new: DataRef] = {
        oldLen: NAT = IF old # NIL THEN old.length ELSE 0;
        n _ MAX[n, 0];
        new _ NEW[Data[MIN[n+n/2+1, NAT.LAST]]];
        FOR i: NAT _ 0, i+1 WHILE i<oldLen DO
            new.counts[i] _ old.counts[i];
            ENDLOOP;
        FOR i: NAT _ oldLen, i+1 WHILE i<new.length DO
            new.counts[i] _ 0;
            ENDLOOP;
        RETURN};

    Ensure2: PROC [h: Histogram, i, j: INT] RETURNS [errMsg: ROPE] ~ {
        IF i < h.nI AND j < h.nJ THEN RETURN [NIL];
        IF i>=NAT.LAST OR j>=NAT.LAST THEN RETURN ["too much data"];
        {oldLen: IntPair ~ [h.nI, h.nJ];
        newLen: IntPair ~ [
            i: MIN[MAX[i+i/2+1, h.nI], NAT.LAST],
            j: MIN[MAX[j+j/2+1, h.nJ], NAT.LAST]];
        newArea: INT ~ INT[newLen.i]*newLen.j;
        IF newArea>=NAT.LAST THEN RETURN ["too much data"];
        {old: DataRef ~ h.data;
        new: DataRef ~ h.data _ NEW[Data[newArea]];
        nk, ok: NAT _ 0;
        FOR i: NAT IN [0 .. oldLen.i) DO
            FOR j: NAT IN [0 .. oldLen.j) DO
                new[nk] _ old[ok];
                nk _ nk + 1;
                ok _ ok + 1;
                ENDLOOP;
            FOR j: NAT IN [oldLen.j .. newLen.j) DO
                new[nk] _ 0;
                nk _ nk + 1;
                ENDLOOP;
            ENDLOOP;
        FOR i: NAT IN [oldLen.i .. newLen.i) DO
            FOR j: NAT IN [0 .. newLen.j) DO
                new[nk] _ 0;
                nk _ nk + 1;
                ENDLOOP;
            ENDLOOP;
        [h.nI, h.nJ] _ newLen;
        RETURN}}};

    AddNoter: PUBLIC ENTRY PROC [h: Histogram, n: Noter] ~ {
        ENABLE UNWIND => NULL;
        h.noters _ CONS[n, h.noters];
        RETURN};

    RemNoter: PUBLIC ENTRY PROC [h: Histogram, n: Noter] ~ {
        ENABLE UNWIND => NULL;
        last: NoterList _ NIL;
        FOR nl: NoterList _ h.noters, nl.rest WHILE nl # NIL DO
            IF nl.first # n THEN last _ nl
            ELSE IF last=NIL THEN h.noters _ nl.rest
            ELSE last.rest _ nl.rest;
            ENDLOOP;
        last _ last; RETURN};

    BroadcastChange: INTERNAL PROC [h: Histogram, c: ChangeNote] ~ {
        h _ h;
        FOR nl: NoterList _ h.noters, nl.rest WHILE nl#NIL DO
            nl.first.Note[nl.first.data, h, c];
            ENDLOOP;
        RETURN};

    Change: PUBLIC ENTRY PROC [h: Histogram, i: NAT, howMuch: INTEGER] = {
        ENABLE UNWIND => NULL;
        IF h.dimensionality # 1 THEN RETURN WITH ERROR Error["This is a 1-D proc, stupid"];
        {errMsg: ROPE = ChangeInt[h, i, howMuch];
        IF errMsg # NIL THEN RETURN WITH ERROR Error[errMsg];
        }};

    ChangeInt: INTERNAL PROC [h: Histogram, i: NAT, howMuch: INTEGER] RETURNS [errMsg: ROPE _ NIL] = {
        IF i < h.data.length THEN NULL
        ELSE IF i >= NAT.LAST THEN RETURN["too much data"]
        ELSE h.data _ Ensure1[h.data, i];
        IF (IF howMuch < 0 THEN h.data.counts[i] < Count[-howMuch] ELSE h.data.counts[i] > Count.LAST-Count[howMuch]) THEN RETURN["counter 
        tried to go out of bounds"];
        h.data.counts[i] _ h.data.counts[i] + howMuch;
        IF h.maxValid THEN SELECT howMuch FROM
            <  0 => h.maxValid _ FALSE;
            >= 0 => h.maxCount _ MAX[h.maxCount, h.data.counts[i]];
            ENDCASE => ERROR;
        BroadcastChange[h, [range: [[i, i], [0, 0]], change: add[howMuch]]];
        };

    Change2: INTERNAL PROC [h: Histogram, i, j: INT, delta: INTEGER] RETURNS [errMsg: ROPE _ NIL] = {
        IF (errMsg _ Ensure2[h, i, j])#NIL THEN RETURN;
        {index: INT ~ i*h.nJ + j;
        oc: Count ~ h.data.counts[index];
        IF (IF delta < 0 THEN oc < Count[-delta] ELSE oc > Count.LAST-Count[delta]) THEN RETURN["counter tried to go out of bounds"];
        {nc: Count ~ oc + delta;
        h.data.counts[index] _ nc;
        IF h.maxValid THEN SELECT delta FROM
            <  0 => h.maxValid _ FALSE;
            >= 0 => h.maxCount _ MAX[h.maxCount, nc];
            ENDCASE => ERROR;
        BroadcastChange[h, [[[i, i], [j, j]], add[delta]]];
        RETURN}}};

    IncrementTransformed: PUBLIC ENTRY PROC [h: Histogram, xmin, xmax, x: REAL] = {
        ENABLE UNWIND => NULL;
        i: INTEGER ~ InvertI[h, MIN[xmax, MAX[xmin, x]]];
        IF i<0 THEN RETURN WITH ERROR Error["That goes to a negative bin, you fool"];
        {errMsg: ROPE = ChangeInt[h, i, 1];
        IF errMsg # NIL THEN RETURN WITH ERROR Error[errMsg];
        }};

    Change1DTransformed: PUBLIC ENTRY PROC [h: Histogram, x: REAL, delta: INTEGER _ 1] ~ {
        i: INTEGER ~ InvertI[h, x];
        IF h.dimensionality#1 THEN RETURN WITH ERROR Error["not 1d"];
        IF i<0 THEN RETURN WITH ERROR Error["sample out of range"];
        {errMsg: ROPE ~ ChangeInt[h, i, delta];
        IF errMsg # NIL THEN RETURN WITH ERROR Error[errMsg];
        RETURN}};

    ChangeTransformed: PUBLIC ENTRY PROC [h: Histogram, x: REAL, y: REAL _ 0, delta: INTEGER _ 1] = {
        ENABLE UNWIND => NULL;
        ip: IntPair ~ IF h.dimensionality=2 THEN Invert2[h, x, y] ELSE [InvertI[h, x], 0];
        IF ip.i<0 OR ip.j<0 THEN RETURN WITH ERROR Error["sample point out of range"];
        {errMsg: ROPE = SELECT h.dimensionality FROM
            1 => ChangeInt[h, ip.i, delta],
            2 => Change2[h, ip.i, ip.j, delta],
            ENDCASE => ERROR;
        IF errMsg # NIL THEN RETURN WITH ERROR Error[errMsg];
        }};

    ClearAll: PUBLIC ENTRY PROC [h: Histogram] = {
        ENABLE UNWIND => NULL;
        FOR i: NAT IN [0 .. h.data.length) DO h.data[i] _ 0 ENDLOOP;
        h.maxCount _ 0; h.maxValid _ TRUE;
        SELECT h.dimensionality FROM
            1 => BroadcastChange[h, [[[0, h.data.length-1], [0, 0]], set[0]]];
            2 => BroadcastChange[h, [[[0, h.nI-1], [0, h.nJ-1]], set[0]]];
            ENDCASE => ERROR;
        };

    ReadSpecs: PUBLIC PROC [from: IO.STREAM] RETURNS [create: CreateData _ NIL, show: ShowData _ NIL] ~ {
        ENABLE IO.Error => IF ec=SyntaxError THEN Error["Syntax"];
        dim: NATURAL _ 0;
        DO
            [] _ from.SkipWhitespace[];
            SELECT from.PeekChar[] FROM
                '1 => {dim _ 1; EXIT};
                '2 => {dim _ 2; EXIT};
                's => {
                    show _ NEW [ShowDataPrivate _ []];
                    IF NOT from.GetID[].Equal["show"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal["["] THEN Error["Syntax"];
                    IF NOT from.GetID[].Equal["format"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal[":"] THEN Error["Syntax"];
                    show.format _ from.GetRopeLiteral[];
                    IF NOT from.GetCedarTokenRope[].token.Equal[","] THEN Error["Syntax"];
                    IF NOT from.GetID[].Equal["width"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal[":"] THEN Error["Syntax"];
                    show.width _ from.GetInt[];
                    IF NOT from.GetCedarTokenRope[].token.Equal[","] THEN Error["Syntax"];
                    IF NOT from.GetID[].Equal["base"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal[":"] THEN Error["Syntax"];
                    show.base _ from.GetReal[];
                    IF NOT from.GetCedarTokenRope[].token.Equal["]"] THEN Error["Syntax"];
                    };
                'o => {
                    cd: REF CreateDataPrivate[oneD] ~ NEW [CreateDataPrivate[oneD] _ [oneD[]]];
                    IF NOT from.GetID[].Equal["oneD"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal["["] THEN Error["Syntax"];
                    IF NOT from.GetID[].Equal["factor"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal[":"] THEN Error["Syntax"];
                    cd.factor _ from.GetReal[];
                    IF NOT from.GetCedarTokenRope[].token.Equal[","] THEN Error["Syntax"];
                    IF NOT from.GetID[].Equal["offset"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal[":"] THEN Error["Syntax"];
                    cd.offset _ from.GetReal[];
                    IF NOT from.GetCedarTokenRope[].token.Equal[","] THEN Error["Syntax"];
                    IF NOT from.GetID[].Equal["log"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal[":"] THEN Error["Syntax"];
                    cd.log _ from.GetBool[];
                    IF NOT from.GetCedarTokenRope[].token.Equal["]"] THEN Error["Syntax"];
                    create _ cd};
                't => {
                    cd: REF CreateDataPrivate[twoD] ~ NEW [CreateDataPrivate[twoD] _ [twoD[]]];
                    IF NOT from.GetID[].Equal["twoD"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal["["] THEN Error["Syntax"];
                    IF NOT from.GetID[].Equal["factors"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal[":"] THEN Error["Syntax"];
                    cd.iFactor _ from.GetReal[];
                    cd.jFactor _ from.GetReal[];
                    IF NOT from.GetCedarTokenRope[].token.Equal[","] THEN Error["Syntax"];
                    IF NOT from.GetID[].Equal["offsets"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal[":"] THEN Error["Syntax"];
                    cd.iOffset _ from.GetReal[];
                    cd.jOffset _ from.GetReal[];
                    IF NOT from.GetCedarTokenRope[].token.Equal[","] THEN Error["Syntax"];
                    IF NOT from.GetID[].Equal["log"] THEN Error["Syntax"];
                    IF NOT from.GetCedarTokenRope[].token.Equal[":"] THEN Error["Syntax"];
                    cd.log[X] _ from.GetBool[];
                    cd.log[Y] _ from.GetBool[];
                    IF NOT from.GetCedarTokenRope[].token.Equal["]"] THEN Error["Syntax"];
                    create _ cd};
                ENDCASE => Error["Syntax"];
            ENDLOOP;
        IF show=NIL THEN show _ NEW [ShowDataPrivate _ [base: 2]];
        IF create=NIL THEN SELECT dim FROM
            1 => create _ NEW [CreateDataPrivate _ [oneD[]]];
            2 => create _ NEW [CreateDataPrivate _ [twoD[]]];
            ENDCASE => ERROR;
        RETURN};

    AddFrom: PUBLIC PROC [h: Histogram, from: IO.STREAM] ~ {
        ENABLE IO.Error => IF ec=SyntaxError THEN Error["Syntax"];
        dimKey: CHAR ~ from.GetChar[];
        IF NOT dimKey IN ['1 .. '2] THEN Error["Syntax"];
        IF h.dimensionality # dimKey - '0 THEN ERROR;
        IF from.GetChar[] # '[ THEN Error["Syntax"];
        SELECT h.dimensionality FROM
            1 => {
                first: BOOL _ TRUE;
                WHILE from.PeekChar[] # '] DO
                    x: REAL;
                    c: CARD;
                    IF first THEN first _ FALSE ELSE IF from.GetChar[] # ', THEN Error["Syntax"];
                    c _ from.GetCard[];
                    IF NOT from.GetTokenRope[IO.IDProc].token.Equal["of"] THEN Error["Syntax"];
                    x _ from.GetReal[];
                    WHILE c>0 DO
                        d: INTEGER ~ MIN[c, INTEGER.LAST];
                        ChangeTransformed[h, x, 0, d];
                        c _ c - d;
                        ENDLOOP;
                    IF from.PeekChar[] = '" THEN [] _ from.GetRopeLiteral[];
                    ENDLOOP;
                h _ h};
            2 => {
                firstI: BOOL _ TRUE;
                WHILE from.PeekChar[] # '] DO
                    firstJ: BOOL _ TRUE;
                    x: REAL;
                    IF firstI THEN firstI _ FALSE ELSE IF from.GetChar[] # ', THEN Error["Syntax"];
                    x _ from.GetReal[];
                    IF from.GetChar[] # ': THEN Error["Syntax"];
                    IF from.GetChar[] # '  THEN Error["Syntax"];
                    IF from.GetChar[] # '[ THEN Error["Syntax"];
                    WHILE from.PeekChar[] # '] DO
                        y: REAL;
                        c: CARD;
                        IF firstJ THEN firstJ _ FALSE ELSE IF from.GetChar[] # ', THEN Error["Syntax"];
                        y _ from.GetReal[];
                        IF from.PeekChar[] = '" THEN [] _ from.GetRopeLiteral[];
                        IF from.GetChar[] # ': THEN Error["Syntax"];
                        c _ from.GetCard[];
                        WHILE c>0 DO
                            d: INTEGER ~ MIN[c, INTEGER.LAST];
                            ChangeTransformed[h, x, y, d];
                            c _ c - d;
                            ENDLOOP;
                        ENDLOOP;
                    IF NOT firstJ THEN IF from.GetChar[] # '] THEN Error["Syntax"];
                    ENDLOOP;
                h _ h};
            ENDCASE => ERROR;
        IF from.GetChar[] # '] THEN Error["Syntax"];
        RETURN};

    Enumerate: PUBLIC PROC [h: Histogram, range: RealRange2 _ fullRealRange2, Consume: PROC [x, y: REAL, n: CARD] RETURNS [BOOL--TRUE 
    to stop--]] RETURNS [found: BOOL _ FALSE, x, y: REAL _ 0.0] ~ {
        d: DataRef ~ h.data;
        min: IntPair ~ Invert2[h, MAX[range[X].min, ExvertI[h, 0]], MAX[range[Y].min, ExvertJ[h, 0]]];
        max: IntPair ~ Invert2[h, MIN[range[X].max, ExvertI[h, h.nI-1]], MIN[range[Y].max, ExvertJ[h, h.nJ-1]]];
        SELECT h.dimensionality FROM
        1 => {
            FOR i: NAT IN [min.i .. max.i] DO
                n: CARDINAL ~ d[i];
                IF n#0 THEN {
                    x _ ExvertI[h, i];
                    IF Consume[x, 0.0, n] THEN RETURN [TRUE, x, 0.0]};
                ENDLOOP;
            RETURN [FALSE]};
        2 => {
            FOR i: INT IN [min.i .. max.i] DO
                index: INT _ i*h.nJ + min.j;
                x _ ExvertI[h, i];
                FOR j: INT IN [min.j .. max.j] DO
                    n: CARDINAL ~ d[index];
                    IF n#0 THEN {
                        y _ ExvertJ[h, j];
                        IF Consume[x, y, n] THEN RETURN [TRUE, x, y]};
                    index _ index + 1;
                    ENDLOOP;
                ENDLOOP;
            RETURN [FALSE]};
        ENDCASE => ERROR};

    Project: PUBLIC PROC [h: Histogram, cx, cy: REAL] RETURNS [p: Histogram] ~ {
        Consume: PROC [x, y: REAL, n: CARD] RETURNS [BOOL] ~ {
            Change1DTransformed[p, cx*x + cy*y, n];
            RETURN [FALSE]};
        xo: REAL ~ ExvertI[h, 0]*cx;
        yo: REAL ~ ExvertJ[h, 0]*cy;
        IF h.dimensionality#2 THEN ERROR;
        IF h.log[X] # h.log[Y] THEN ERROR;
        p _ Create1D[MaybeEx[(h.iFactor*cx+h.jFactor*cy)/(cx+cy), h.log[X]], IF cx=0 THEN yo ELSE IF cy=0 THEN xo ELSE MIN[xo, yo], 
        h.log[X], NIL, NIL];
        IF Enumerate[h, fullRealRange2, Consume].found THEN ERROR;
        RETURN};

    END.