diff --git a/pyHealpix_perftest.py b/pyHealpix_perftest.py
new file mode 100644
index 0000000000000000000000000000000000000000..6a6851c8781e5f5129f1a4009d4b75c54d5d9c5f
--- /dev/null
+++ b/pyHealpix_perftest.py
@@ -0,0 +1,81 @@
+from __future__ import print_function
+import time
+import math
+import numpy as np
+import pyHealpix as ph
+
+def report (name,vlen,ntry,nside,isnest,perf):
+  print (name,": ",perf*1e-6,"MOps/s",sep="")
+
+def random_ptg(vlen):
+  res=np.empty((vlen,2),dtype=np.float64)
+  res[:,0]=np.arccos((np.random.random_sample(vlen)-0.5)*2)
+  res[:,1]=np.random.random_sample(vlen)*2*math.pi
+  return res
+
+def random_pix(nside,vlen):
+  return np.random.randint(low=0,high=12*nside*nside-1,size=vlen,dtype=np.int64)
+
+def dummy(vlen):
+  inp=np.zeros(vlen,dtype=np.int64)
+
+def genperf(func,fname,inp,vlen,ntry,nside,isnest):
+  cnt=0
+  t=time.time()
+  while (cnt<ntry):
+    func(inp)
+    cnt+=1
+  t=time.time()-t
+  p=(vlen*ntry)/t
+  report (fname,vlen,ntry,nside,isnest,p)
+
+def perf_pix2ang(vlen,ntry,nside,isnest):
+  inp=random_pix(nside,vlen)
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  genperf(base.pix2ang,"pix2ang",inp,vlen,ntry,nside,isnest)
+def perf_ang2pix(vlen,ntry,nside,isnest):
+  inp=random_ptg(vlen)
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  genperf(base.ang2pix,"ang2pix",inp,vlen,ntry,nside,isnest)
+
+def perf_pix2vec(vlen,ntry,nside,isnest):
+  inp=random_pix(nside,vlen)
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  genperf(base.pix2vec,"pix2vec",inp,vlen,ntry,nside,isnest)
+def perf_vec2pix(vlen,ntry,nside,isnest):
+  inp=ph.ang2vec(random_ptg(vlen))
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  genperf(base.vec2pix,"vec2pix",inp,vlen,ntry,nside,isnest)
+
+def perf_ring2nest(vlen,ntry,nside,isnest):
+  inp=random_pix(nside,vlen)
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  genperf(base.ring2nest,"ring2nest",inp,vlen,ntry,nside,isnest)
+def perf_nest2ring(vlen,ntry,nside,isnest):
+  inp=random_pix(nside,vlen)
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  genperf(base.nest2ring,"nest2ring",inp,vlen,ntry,nside,isnest)
+
+def perf_neighbors(vlen,ntry,nside,isnest):
+  inp=random_pix(nside,vlen)
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  genperf(base.neighbors,"neighbors",inp,vlen,ntry,nside,isnest)
+
+def suite (vlen,ntry,nside,isnest):
+  print ("vlen=",vlen,", ","NEST" if isnest else "RING",sep="")
+  dummy(vlen)
+  perf_pix2ang(vlen,ntry,nside,isnest)
+  perf_ang2pix(vlen,ntry,nside,isnest)
+  perf_pix2vec(vlen,ntry,nside,isnest)
+  perf_vec2pix(vlen,ntry,nside,isnest)
+  perf_neighbors(vlen,ntry,nside,isnest)
+
+nside=512
+ntry=1000
+print ("nside=",nside,sep="")
+for vlen in (1,10,100,1000,10000):
+  for isnest in (True, False):
+    suite(vlen,ntry,nside,isnest)
+    perf_ring2nest(vlen,ntry,nside,isnest)
+    perf_nest2ring(vlen,ntry,nside,isnest)
+    print()
diff --git a/pyHealpix_test.py b/pyHealpix_test.py
new file mode 100644
index 0000000000000000000000000000000000000000..25099a91128ad29950e156d2bd60c95d80ff23e4
--- /dev/null
+++ b/pyHealpix_test.py
@@ -0,0 +1,109 @@
+import pyHealpix as ph
+import numpy as np
+import math
+
+def random_ptg(vlen):
+  res=np.empty((vlen,2),dtype=np.float64)
+  res[:,0]=np.arccos((np.random.random_sample(vlen)-0.5)*2)
+  res[:,1]=np.random.random_sample(vlen)*2*math.pi
+  return res
+
+def check_pixangpix(vlen,ntry,nside,isnest):
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  cnt=0
+  while (cnt<ntry):
+    cnt+=1
+    inp=np.random.randint(low=0,high=12*nside*nside-1,size=vlen)
+    out=base.ang2pix(base.pix2ang(inp))
+    if (np.array_equal(inp,out)==False):
+      raise ValueError("Test failed")
+
+def check_vecpixvec(vlen,ntry,nside,isnest):
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  cnt=0
+  while (cnt<ntry):
+    cnt+=1
+    inp=ph.ang2vec(random_ptg(vlen))
+    out=base.pix2vec(base.vec2pix(inp))
+    if (np.any(np.greater(ph.v_angle(inp,out),base.max_pixrad()))):
+      raise ValueError("Test failed")
+
+def check_pixangvecpix(vlen,ntry,nside,isnest):
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  cnt=0
+  while (cnt<ntry):
+    cnt+=1
+    inp=np.random.randint(low=0,high=12*nside*nside-1,size=vlen)
+    out=base.vec2pix(ph.ang2vec(base.pix2ang(inp)))
+    if (np.array_equal(inp,out)==False):
+      raise ValueError("Test failed")
+
+def check_pixvecangpix(vlen,ntry,nside,isnest):
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  cnt=0
+  while (cnt<ntry):
+    cnt+=1
+    inp=np.random.randint(low=0,high=12*nside*nside-1,size=vlen)
+    out=base.ang2pix(ph.vec2ang(base.pix2vec(inp)))
+    if (np.array_equal(inp,out)==False):
+      raise ValueError("Test failed")
+
+def check_pixvecpix(vlen,ntry,nside,isnest):
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  cnt=0
+  while (cnt<ntry):
+    cnt+=1
+    inp=np.random.randint(low=0,high=12*nside*nside-1,size=vlen)
+    out=base.vec2pix(base.pix2vec(inp))
+    if (np.array_equal(inp,out)==False):
+      raise ValueError("Test failed")
+
+def check_ringnestring(vlen,ntry,nside):
+  base=ph.Healpix_Base (nside,"NEST")
+  cnt=0
+  while (cnt<ntry):
+    cnt+=1
+    inp=np.random.randint(low=0,high=12*nside*nside-1,size=vlen)
+    out=base.nest2ring(base.ring2nest(inp))
+    if (np.array_equal(inp,out)==False):
+      raise ValueError("Test failed")
+
+def check_pixxyfpix(vlen,ntry,nside,isnest):
+  base=ph.Healpix_Base (nside,"NEST" if isnest else "RING")
+  cnt=0
+  while (cnt<ntry):
+    cnt+=1
+    inp=np.random.randint(low=0,high=12*nside*nside-1,size=vlen)
+    out=base.xyf2pix(base.pix2xyf(inp))
+    if (np.array_equal(inp,out)==False):
+      raise ValueError("Test failed")
+
+def check_vecangvec(vlen,ntry):
+  cnt=0
+  while (cnt<ntry):
+    cnt+=1
+    inp=random_ptg(vlen)
+    out=ph.vec2ang(ph.ang2vec(inp))
+    if (np.any(np.greater(np.abs(inp-out),1e-10))):
+      raise ValueError("Test failed")
+
+check_vecangvec(1000,1000)
+
+for nside in (1,32,512,8192,32768*8):
+  check_ringnestring(1000,1000,nside)
+  for isnest in (False,True):
+    check_vecpixvec(1000,1000,nside,isnest)
+    check_pixangpix(1000,1000,nside,isnest)
+    check_pixvecpix(1000,1000,nside,isnest)
+    check_pixxyfpix(1000,1000,nside,isnest)
+    check_pixangvecpix(1000,1000,nside,isnest)
+    check_pixvecangpix(1000,1000,nside,isnest)
+
+isnest=False
+for nside in (3,7,514,8167,32768*8+7):
+  check_vecpixvec(1000,1000,nside,isnest)
+  check_pixangpix(1000,1000,nside,isnest)
+  check_pixvecpix(1000,1000,nside,isnest)
+  check_pixxyfpix(1000,1000,nside,isnest)
+  check_pixangvecpix(1000,1000,nside,isnest)
+  check_pixvecangpix(1000,1000,nside,isnest)