code: def get_nnsoap(obj, first_shell, alphas, betas, rcut=6, nmax=10, lmax=9, all_atomtypes=[]):
      
  
      soap_vector = []
      nnn = len(first_shell)
      for tbh in range(0,3):
          try:
              atom_idx = first_shell[tbh]
          except:
              soap_vector.append(soap_zero)
          else:
              Hpos = []
              print(atom_idx)
              pos = obj.get_positions()[atom_idx]
              Hpos.append(pos)
              x = soapPy.get_soap_locals(obj, Hpos, myAlphas, myBetas, rCut=rcut, NradBas=nmax, Lmax=lmax,crossOver=False, all_atomtypes=all_atomtypes)
              soap_zero = np.zeros(x.shape)
              soap_vector.append(x)
      print(len(soap_vector), soap_vector[0].shape, soap_vector[1].shape, soap_vector[2].shape)
      print("exemplary soapvalues",soap_vector[0][0,1], soap_vector[1][0,1], soap_vector[2][0,1])
      soap_array = np.hstack(soap_vector)
      return soap_array
docs: Takes cluster structure and nearest neighbour information of a datapoint,
      Returns concatenated soap vectors for each nearest
      neighbour (up to 3). Top, bridge, hollow fill the initial
      zero soap vector from left to right.