function cmb_angle_dot,a,b,absolute=absolute ;compute dot product angle between the unit vectors of 'a' and 'b' ;a = a[nc,nt] or a[nc] & b = b[nc,nt] or b=b[nc] -nc number of components in vector, nt-vectors at different time steps if cmb_var_type(a) eq 'DOUBLE' or cmb_var_type(b) eq 'DOUBLE' then idouble =1 si = size(a) if si[0] eq 1 then begin if norm(a) eq 0 or norm(b) eq 0 then return,!pi a0 = a/sqrt(total(a^2, double = idouble)) b0 = b/sqrt(total(b^2, double = idouble)) dot = total(a0*b0, double = idouble) if abs(dot) gt 1 then dot = cmb_sgn(dot) ang = acos(dot) if keyword_set(absolute) then dot = abs(dot) if finite(ang) ne 1 then stop return,acos(dot) endif ;a[3,n],b[3,n] a0 = a b0 = b am = sqrt(total(a^2,1, double = idouble)) bm = sqrt(total(b^2,1, double = idouble)) ;help,am,bm for i=0,si[1]-1 do begin a0[i,*] = a0[i,*]/am b0[i,*] = b0[i,*]/bm endfor dot = total( a0*b0,1, double = idouble) if keyword_set(absolute) then dot = abs(dot) ii=where( abs(dot) gt 1) if ii[0] ne -1 then begin dot[ii] = cmb_sgn(dot[ii]) endif return,acos(dot) end