;+ ; NAME: ; MULTISORT ; ; AUTHOR: ; Craig B. Markwardt, NASA/GSFC Code 662, Greenbelt, MD 20770 ; craigm@lheamail.gsfc.nasa.gov ; ; PURPOSE: ; Perform a sorting operation with multiple sort keys ; ; CALLING SEQUENCE: ; INDICES = MULTISORT(KEY1, KEY2, ..., [/L64, ], [ORDER=order]) ; ; DESCRIPTION: ; ; The function MULTISORT performs a sorting operation with multiple ; sort keys. Unlike the IDL built-in SORT() function, which can ; only sort a single key, MULTISORT can accept multiple keys. In ; cases where the primary key is equal, the sort order is based on ; any secondary keys provided. The return value is an array of ; indices which will place the key arrays into sorted order. ; ; MULTISORT works by building an internal sort key string which can ; be sorted in a single pass. Because MULTISORT is not a built-in ; function, and because it must build these auxiliary strings, it ; cannot be as fast or memory-efficient as the built-in function. ; Users will need several times more memory than the memory used ; to store just the input keys. ; ; MULTISORT() allows the user to choose the sort order for each key ; separately. The ORDER keyword is an N-vector, one order for each ; input key. ORDER[i] is +1 to sort KEYi ascending, and ORDER[i] is ; -1 to sort KEYi descending. ; ; INPUTS: ; ; KEY1, KEY2, ... KEY9 - input sort keys. Any integer, floating ; point or string value is allowed. The ; number of values must be the same for each ; key. ; ; ; ; KEYWORDS: ; ; ORDER - an N-vector, giving the sort order for each key (see ; documentation above). ; Default: REPLICATE(+1,N_PARAMS()) (all keys ascending) ; ; L64 - if set, then return a LONG64 index instead of LONG. ; ; RETURNS: ; ; An array of indices which will place the keys into sorted order. ; I.e., KEYS1[INDICES], KEYS2[INDICES] ... will be in order. ; ; SEE ALSO: ; ; SORT ; ; MODIFICATION HISTORY: ; Written, CM, Jun 2007 ; Document the encoding format, and make some floating point ; operations more efficient, CM, Jan 2008 ; ; $Id: multisort.pro,v 1.2 2008/01/13 01:07:01 craigm Exp $ ; ;- ; Copyright (C) 2007, 2008, Craig Markwardt ; This software is provided as is without any warranty whatsoever. ; Permission to use, copy, modify, and distribute modified or ; unmodified copies is granted, provided this copyright and disclaimer ; are included unchanged. ;- ; ======================================================== ; Utility function to transform an integer into a string key ; Integers are coded like this: ; 123 = 'B0000123' ; -4567 = 'A0004567' ; where the number of digits is given by LEN. The A and B ; ensure proper sort order for signed integers. function multisort_intkey, x, len, unsigned=u, order=order COMPILE_OPT strictarr n = n_elements(x) if order LT 0 then x1 = -temporary(x) $ ;; Reverse order else x1 = temporary(x) slen = strtrim(len,2) fmt = '(I'+slen+'.'+slen+')' ;; (In.n) - zero-padded key = string(abs(x1),format=fmt) if NOT keyword_set(u) then begin key = 'B'+temporary(key) ;; Prefix to indicate positive values wh = where(x1 LT 0, ct) if ct GT 0 then begin ;; change the digits 0 -> 9, 9 -> 0 bb = byte(temporary(key)) bb[1:*,wh] = 105b - bb[1:*,wh] bb[0,wh] = 65b ;; 'A' - Prefix to indicate negative values key = string(temporary(bb)) endif endif return, key end ; ======================================================== ; Utility function to transform a float into a string key ; Floating point numbers are coded as: ; +X.XXXXXXE+YYY = 'BBYYYX.XXXXXX' ; -X.XXXXXXE+YYY = 'ABYYYX.XXXXXX' ; +X.XXXXXXE-YYY = 'BAYYYX.XXXXXX' ; -X.XXXXXXE-YYY = 'AAYYYX.XXXXXX' ; +Infinity = 'BB9991.111111' (similar for -Infinity) ; NaN = 'BBaaa1.111111' ; In other words: ; * Character 0 = A or B, depending on sign of X (B = +) ; * Character 1 = A or B, depending on sign of exponent ; * Character 2-5 = exponent ; * Character 6-* = mantissa ; function multisort_fltkey, x, type, order=order COMPILE_OPT strictarr n = n_elements(x) if type EQ 4 then len = 17 else len = 8 if order LT 0 then x1 = -temporary(x) $ ;; Reverse order else x1 = temporary(x) mlen = strtrim(len,2) ;; Mantissa length tlen = strtrim(len+8,2) ;; Total length including exponent fmt = '(E'+tlen+'.'+mlen+')' ;; (Et.m) xkey = string(abs(x1),format=fmt) mant0= string(bytarr(mlen)+49b) ;; Default mantissa: all 1's key = strarr(n) nan = x1 NE x1 ;; Is this a NaN? (cannot equal itself) wh = where(NOT finite(x1) AND nan EQ 0, ct) ;; Infinity (but not NaN) if ct GT 0 then key[wh] = 'BB999'+mant0 wh = where(nan,ct) ;; NaN (doesn't equal itself) if ct GT 0 then key[wh] = 'BBaaa'+mant0 nan = 0 iexp = strpos(xkey,'E') ;; Normal scientific notation wh = where(iexp GE 0, ct) if ct GT 0 then begin iexp = reform(iexp[wh],1,ct) xkeyi = (temporary(xkey))[wh] wh = 0 ;; Mantissa key[wh] = 'B' + $ ;; Default positive multisort_intkey(fix(strmid(xkeyi,iexp+1)),3,order=+1) + $;; Exponent strtrim(strmid(temporary(xkeyi),0,temporary(iexp)),2) endif xkey = 0 xkeyi = 0 mant1 = 0 expo1 = 0 iexp = 0 wh = where(x1 LT 0, ct) if ct GT 0 then begin ;; change the digits 0 -> 9, 9 -> 0 bb = byte(temporary(key)) bb[1:*,wh] = 105b - bb[1:*,wh] bb[0,wh] = 65b ;; 'A' - Prefix to indicate negative values key = string(temporary(bb)) endif return, key end ; ======================================================== ; Utility function to transform a string into a string key function multisort_strkey, x, order=order COMPILE_OPT strictarr len = strlen(x) maxlen = max(len, min=minlen) if maxlen GT minlen then begin ;; Pad out to the maximum string length (i.e. left-align the strings) pad = string(bytarr(maxlen-minlen)+32b) key = strmid(x+pad,0,maxlen) endif else begin key = x endelse ;; Reverse order if requested if order LT 0 then begin key = string( (255b - byte(temporary(key))) > 1b ) endif return, key end ; ======================================================== ; MAIN ROUTINE ; ======================================================== function multisort, x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, L64=L64, $ keys=keys0, order=order0 COMPILE_OPT strictarr nkeys = n_params() if nkeys EQ 0 then begin message, 'USAGE: INDICES = MULTISORT(KEY1[,KEY2,KEY3,...])', /info return, -1L endif order = intarr(nkeys) + 1 if n_elements(order0) GT 0 then order[0] = round(order0) ;; Special case: only one term, no need to do complicate sort key ;; manipulations. if nkeys EQ 1 AND order[0] EQ +1 then begin return, sort(x0, L64=L64) endif ;; Master key mkey = '' for i = 0, nkeys-1 do begin xi = 0 & dummy = temporary(xi) case i of 0: xi = x0 1: xi = x1 2: xi = x2 3: xi = x3 4: xi = x4 5: xi = x5 6: xi = x6 7: xi = x7 8: xi = x8 9: xi = x9 endcase if n_elements(xi) EQ 0 then begin message, string(i,format='("ERROR: no data was in parameter X",I0)') return, -1L endif sz = size(xi) tp = sz[sz[0]+1] o = order[i] case tp of 1: mkey = temporary(mkey) + multisort_intkey(temporary(xi),3,/u,o=o) ;; BYTE 2: mkey = temporary(mkey) + multisort_intkey(temporary(xi),6,o=o) ;; INT 3: mkey = temporary(mkey) + multisort_intkey(temporary(xi),10,o=o) ;; LONG 4: mkey = temporary(mkey) + multisort_fltkey(temporary(xi),4,o=o) ;; FLOAT 5: mkey = temporary(mkey) + multisort_fltkey(temporary(xi),5,o=o) ;; DOUBLE 7: mkey = temporary(mkey) + multisort_strkey(temporary(xi),o=o) ;; STRING 12: mkey = temporary(mkey) + multisort_intkey(temporary(xi),5,/u,o=o) ;; UINT 13: mkey = temporary(mkey) + multisort_intkey(temporary(xi),10,/u,o=o) ;; ULONG 14: mkey = temporary(mkey) + multisort_intkey(temporary(xi),20,o=o) ;; LONG64 15: mkey = temporary(mkey) + multisort_intkey(temporary(xi),20,/u,o=o) ;; ULONG64 else: begin message, string(tp, i, $ format='("ERROR: data type ",I0," for parameter X,",I0," is not sortable")') return, -1L end endcase xi = 0 endfor return, sort(mkey, L64=L64) end