# GPS L1Cd code construction
#
# Copyright 2018 Peter Monta

import numpy as np
from sympy.ntheory import legendre_symbol

chip_rate = 1023000
code_length = 10230

l1cd_params = {
   1: (5097,181),     2: (5110,359),    3: (5079,72),      4: (4403,1110),
   5: (4121,1480),    6: (5043,5034),   7: (5042,4622),    8: (5104,1),
   9: (4940,4547),   10: (5035,826),   11: (4372,6284),   12: (5064,4195),
  13: (5084,368),    14: (5048,1),     15: (4950,4796),   16: (5019,523),
  17: (5076,151),    18: (3736,713),   19: (4993,9850),   20: (5060,5734),
  21: (5061,34),     22: (5096,6142),  23: (4983,190),    24: (4783,644),
  25: (4991,467),    26: (4815,5384),  27: (4443,801),    28: (4769,594),
  29: (4879,4450),   30: (4894,9437),  31: (4985,4307),   32: (5056,5906),
  33: (4921,378),    34: (5036,9448),  35: (4812,9432),   36: (4838,5849),
  37: (4855,5547),   38: (4904,9546),  39: (4753,9132),   40: (4483,403),
  41: (4942,3766),   42: (4813,3),     43: (4957,684),    44: (4618,9711),
  45: (4669,333),    46: (4969,6124),  47: (5031,10216),  48: (5038,4251),
  49: (4740,9893),   50: (4073,9884),  51: (4843,4627),   52: (4979,4449),
  53: (4867,9798),   54: (4964,985),   55: (5025,4272),   56: (4579,126),
  57: (4390,10024),  58: (4763,434),   59: (4612,1029),   60: (4784,561),
  61: (3716,289),    62: (4703,638),   63: (4851,4353),
  64: (4955,9899),   65: (5018,4629),  66: (4642,669),    67: (4840,4378),
  68: (4961,4528),   69: (4263,9718),  70: (5011,5485),   71: (4922,6222),
  72: (4317,672),    73: (3636,1275),  74: (4884,6083),   75: (5041,5264),
  76: (4912,10167),  77: (4504,1085),  78: (4617,194),    79: (4633,5012),
  80: (4566,4938),   81: (4702,9356),  82: (4758,5057),   83: (4860,866),
  84: (3962,2),      85: (4882,204),   86: (4467,9808),   87: (4730,4365),
  88: (4910,162),    89: (4684,367),   90: (4908,201),    91: (4759,18),
  92: (4880,251),    93: (4095,10167), 94: (4971,21),     95: (4873,685),
  96: (4561,92),     97: (4588,1057),  98: (4773,3),      99: (4997,5756),
 100: (4583,14),    101: (4900,9979), 102: (4574,9569),  103: (4629,515),
 104: (4676,753),   105: (4181,1181), 106: (5057,9442),  107: (4944,669),
 108: (4401,4834),  109: (4586,541),  110: (4699,9933),  111: (3676,6683),
 112: (4387,4828),  113: (4866,9710), 114: (4926,10170), 115: (4657,9629),
 116: (4477,260),   117: (4359,86),   118: (4673,5544),  119: (4258,923),
 120: (4447,257),   121: (4570,507),  122: (4486,4572),  123: (4362,4491),
 124: (4481,341),   125: (4322,130),  126: (4668,79),    127: (3967,1142),
 128: (4374,448),   129: (4553,875),  130: (4641,555),   131: (4215,1272),
 132: (3853,5198),  133: (4787,9529), 134: (4266,4459),  135: (4199,10019),
 136: (4545,9353),  137: (4208,9780), 138: (4485,375),   139: (3714,503),
 140: (4407,4507),  141: (4182,875),  142: (4203,1246),  143: (3788,1),
 144: (4471,4534),  145: (4691,8),    146: (4281,9549),  147: (4410,6240),
 148: (3953,22),    149: (3465,5652), 150: (4801,10069), 151: (4278,4796),
 152: (4546,4980),  153: (3779,27),   154: (4115,90),    155: (4193,9788),
 156: (3372,715),   157: (3786,9720), 158: (3491,301),   159: (3812,5450),
 160: (3594,5215),  161: (4028,13),   162: (3652,1147),  163: (4224,4855),
 164: (4334,1190),  165: (3245,1267), 166: (3921,1302),  167: (3840,1),
 168: (3514,5007),  169: (2922,549),  170: (4227,368),   171: (3376,6300),
 172: (3560,5658),  173: (4989,4302), 174: (4756,851),   175: (4624,4353),
 176: (4446,9618),  177: (4174,9652), 178: (4551,1232),  179: (3972,109),
 180: (4399,10174), 181: (4562,6178), 182: (3133,1851),  183: (4157,1299),
 184: (5053,325),   185: (4536,10206),186: (5067,9968),  187: (3905,10191),
 188: (3721,5438),  189: (3787,10080),190: (4674,219),   191: (3436,758),
 192: (2673,2140),  193: (4834,9753), 194: (4456,4799),  195: (4056,10126),
 196: (3804,241),   197: (3672,1245), 198: (4205,1274),  199: (3348,1456),
 200: (4152,9967),  201: (3883,235),  202: (3473,512),   203: (3669,1078),
 204: (3455,1078),  205: (2318,953),  206: (2945,5647),  207: (2947,669),
 208: (3220,1311),  209: (4052,5827), 210: (2953,15),
}

N = 10223
L = np.array([legendre_symbol(i,N) for i in range(N)])
L[L==-1] = 0
L[0] = 0

def l1cd(prn):
  w,p = l1cd_params[prn]
  W = np.array([L[k]^L[(k+w)%N] for k in range(N)])
  expansion = np.array([0,1,1,0,1,0,0])
  c = np.concatenate((W[0:p-1],expansion,W[p-1:N]))
  return c

codes = {}

def l1cd_code(prn):
  if prn not in codes:
    codes[prn] = l1cd(prn)
  return codes[prn]

def code(prn,chips,frac,incr,n):
  c = l1cd_code(prn)
  idx = (chips%code_length) + frac + incr*np.arange(n)
  idx = np.floor(idx).astype('int')
  idx = np.mod(idx,code_length)
  x = c[idx]
  return 1.0 - 2.0*x

boc11 = np.array([1.0,-1.0])

try:
  from numba import jit
except:
  def jit(**kwargs):
    return lambda x: x

@jit(nopython=True)
def correlate(x,prn,chips,frac,incr,c,boc11):
  n = len(x)
  p = 0.0j
  cp = (chips+frac)%code_length
  bp = (2*(chips+frac))%2
  for i in range(n):
    boc = boc11[int(bp)]
    p += x[i]*(1.0-2.0*c[int(cp)])*boc
    cp = (cp+incr)%code_length
    bp = (bp+2*incr)%2
  return p

# test

def chips2octal(x):
  s = ''
  for i in range(len(x)//3):
    d = 4*x[3*i] + 2*x[3*i+1] + x[3*i+2]
    s = s + '%o'%int(d)
  return s

if __name__=='__main__':
  for prn in range(1,64):
    c = l1cd_code(prn)
    s1 = chips2octal(c[0:24])
    s2 = chips2octal(c[-24:])
    print("%d %s %s"%(prn,s1,s2))
