Line data    Source code

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