1% This file not only provides access to the module, it also provides several
2% other functions that are easier to implement in slang.  These include:
3%
4%   whist1d     : weighted 1-d histogram
5%   whist2d     : weighted 2-d histogram
6%   hist2d_rebin: Rebins a 2-d histogram
7%
8#ifeval _slang_version < 20000
9if (current_namespace () != "")
10  import ("histogram", current_namespace ());
11else
12#endif
13 import ("histogram");
14
15define whist1d ()
16{
17   variable func = NULL;
18   variable rev = NULL;
19
20   switch (_NARGS)
21     {
22      case 4:
23	rev = ();
24     }
25     {
26      case 5:
27	(rev, func) = ();
28     }
29     {
30	if (_NARGS != 3)
31	  usage ("h = %s (pnts, weights, bin-edges [,rev-indices [,func]])",
32		 _function_name);
33     }
34
35   variable pnts, w, edges;
36   (pnts, w, edges) = ();
37
38   if (func == NULL)
39     func = &sum;
40
41   variable r;
42   () = hist1d (pnts, edges, &r);
43
44   variable n = length (edges);
45   variable h = Double_Type [n];
46
47   _for (0, n-1, 1)
48     {
49	variable i = ();
50	h[i] = (@func) (w[r[i]]);
51     }
52
53   if (rev != NULL)
54     @rev = r;
55
56   return h;
57}
58
59define whist2d ()
60{
61   variable rev = NULL, func = NULL;
62   switch (_NARGS)
63     {
64      case 6:
65	rev = ();
66     }
67     {
68      case 7:
69	(rev, func) = ();
70     }
71     {
72	if (_NARGS != 5)
73	  usage ("img=%s(xpnts, ypnts, weights, xgrid, ygrid [,rev [,func]])",
74		 _function_name ());
75     }
76
77   variable x, y, w, xgrid, ygrid;
78   (x, y, w, xgrid, ygrid) = ();
79
80   if (func == NULL)
81     func = &sum;
82
83   variable r;
84   () = hist2d (x, y, xgrid, ygrid, &r);
85
86   variable nx = length (xgrid);
87   variable ny = length (ygrid);
88   variable img = Double_Type[nx, ny];
89   _for (0, nx-1, 1)
90     {
91	variable i = ();
92	_for (0, ny-1, 1)
93	  {
94	     variable j = ();
95	     img[i,j] = (@func)(w[r[i,j]]);
96	  }
97     }
98
99   if (rev != NULL)
100     @rev = r;
101
102   return img;
103}
104
105define hist2d_rebin ()
106{
107   if (_NARGS != 5)
108     usage ("new_ijhist=%s(new_igrid,new_jgrid,old_igrid,old_jgrid,old_ijhist)",
109	    _function_name);
110
111   variable new_igrid, new_jgrid, old_igrid, old_jgrid, old_hist;
112   (new_igrid, new_jgrid, old_igrid, old_jgrid, old_hist) = ();
113
114   variable new_ilen = length (new_igrid);
115   variable new_jlen = length (new_jgrid);
116   variable old_ilen = length (old_igrid);
117   variable old_jlen = length (old_jgrid);
118
119   variable h_i = hist1d_rebin (new_jgrid, old_jgrid, old_hist[0,*]);
120   variable type = _typeof (h_i);
121
122   variable new_hist = @Array_Type (type, [old_ilen, new_jlen]);
123
124   new_hist[0,*] = h_i;
125   _for (1, old_ilen-1, 1)
126     {
127	variable i = ();
128	new_hist[i,*] = hist1d_rebin (new_jgrid, old_jgrid, old_hist[i,*]);
129     }
130   old_hist = new_hist;
131   new_hist = @Array_Type (type, [new_ilen, new_jlen]);
132
133   _for (0, new_jlen-1, 1)
134     {
135	variable j = ();
136	new_hist[*,j] = hist1d_rebin (new_igrid, old_igrid, old_hist[*,j]);
137     }
138
139   return new_hist;
140}
141
142
143#ifexists add_doc_file
144$1 = path_concat (path_concat (path_dirname (__FILE__), "help"),
145		  "histogram.hlp");
146if (NULL != stat_file ($1))
147  add_doc_file ($1);
148#endif
149
150provide ("histogram");
151